磨刀门水道枯水期咸潮变化趋势及成因

2019-03-29王艺霖苏程佳陈晓宏王佳雯

水资源保护 2019年2期

王艺霖，苏程佳，陈晓宏，王佳雯

(1.中山大学水资源与环境研究中心，广东广州 510275； 2.华南地区水循环与水安全广东普通高校重点实验室，广东广州 510275)

咸潮上溯作为河口地区的一大全球性关键水文问题，严重威胁到了当地的供水安全[1-3]。我国珠江河口区河道纵横交错、洪潮作用复杂，是咸潮现象的多发地段[4-6]。以磨刀门水道为例，1992年咸潮上溯至大涌口水闸，持续时间约60 d；1995年上溯至神湾港，持续时间约30 d；1998年上溯至南镇水厂；1999年上溯至全禄水厂[7]；2004—2005年广昌泵站水体含氯度曾连续超标达38 d；2005—2007年枯水期连续遭遇较为严重的咸潮灾害，珠海、澳门等城市最高长达48 d无法正常取水[8]。而最近几年，例如2009—2010年及2011—2012年枯水期期间，磨刀门水道也是爆发了严重的咸潮灾害，其中前者咸潮导致平岗泵站含氯度曾累积超标长达110 d[9]，后者造成马角水闸关闸时间高达73 d[10]。严重的咸潮灾害极大地影响了河口地区正常的供水秩序，因此，研究咸潮的变化规律及其成因，对咸潮的应对与预报有重要的现实意义。

关于磨刀门水道咸潮上溯的变化趋势，目前已有了一定的研究成果。例如，周文浩[11]基于珠江河口区主要代表站的含盐度分析，发现20世纪90年代以前的盐度总体表现为下降趋势；闻平等[12]通过对平岗泵站1998—2006年枯水期氯化物的超标情况进行统计发现，以2002—2003年枯水期为界，前4年咸潮逐渐减弱，后3年则逐步增强。而就磨刀门水道咸潮变化趋势的影响因素而言，前人也作了一定的研究。大量研究表明，径流及潮汐是控制磨刀门水道含氯度变化的关键要素[8,13-14]；部分研究基于实测资料分析，认为河口地形变化可能会影响枯水期咸潮的入侵规律[15-16]。然而，已有的研究大多针对单场次咸潮事件[17-18]，且研究时段较短，未能对磨刀门水道的咸潮演变规律及其影响要素的量化有一个全面深刻的认识，特别是在全球气候变暖的大背景下，极端气候及海平面上升也极大程度地影响着咸潮上溯的强度，因而需要对咸潮上溯的变化规律进行更为全面的掌握。

基于此，本文以磨刀门水道为研究对象，在统计分析1998—2015年磨刀门水道咸潮上溯变化趋势特征的基础上，选取流量、潮位、潮差、降雨以及海平面上升等要素，探究其与咸度超标时数之间的相关关系，并利用主成分分析法量化了各要素对咸度超标时数的贡献率，旨在识别出造成1998—2015年磨刀门水道咸潮活动变化趋势的关键要素，为区域水资源管理及咸潮预测提供参考。

1 研究区域概况及数据来源

磨刀门水道属于珠江八大口门之一，流经江门、中山、珠海等市，是西江的主要入海口之一，同时也是珠海、澳门等城市和地区的主要用水水源地，其地理位置如图1所示。

图1 磨刀门水道位置

研究用到的主要数据包括西江上游马口及三水站的逐日流量和灯笼山站的逐时潮位资料，数据均来源于广东省水文局；广昌泵站及平岗泵站的逐时咸度资料来源于水利部珠江水利委员会。根据GB3838—2002《地表水环境质量标准》的有关规定，集中式生活饮用水地表水源地氯化物的质量浓度不得超过250 mg/L，因此本文采用达到或超过250 mg/L 的小时数，即咸度超标时数表征咸潮上溯活动的强弱。马口及三水站以上西、北江流域的气象站逐日降雨量数据来源于中国气象数据网。此外，还采用了国家海平面上升观测站闸坡站逐月的海平面变化数据，数据可从平均海平面永久服务中心下载(http://www.psmsl.org/)。上述数据的年限均为1998—2015年。

2 研究方法

2.1 Spearman秩相关分析

选用了适用范围较为广泛的Spearman秩相关分析对磨刀门水道各影响因素与咸度超标时数的相关性进行研究[19]。设存在时间序列X(x1,x2,…,xn)和Y(y1,y2,…,yn)，Spearman秩相关系数的计算公式为

(1)

其中

di=R1i-R2i

式中：di为时间序列X和Y的同一组样本位次差；R1i为要素X的位次序号；R2i为要素Y的位次序号。

2.2 主成分分析

主成分分析法是把原来多个变量划为少数几个综合指标的一种统计分析方法[19]。设存在n×p阶的数据矩阵，n代表样本数，p为每个样本的变量个数，则

(2)

基于降维的思想，将原来的变量进行线性组合，通过调整系数使新的变量之间相互独立且代表性最好。记原来的变量为x1，x2，…，xp，新的变量为z1，z2，…，zm(m≤p)，则

(3)

式中:lji为原来变量xj(j=1,2,…,p)在主成分zi(i=1,2,…,m)上的载荷，即x1,x2,…,xp的相关系数矩阵的m个较大的特征值所对应的单位特征向量。详细计算步骤见参考文献[19]。

研究中通常取特征值接近于1且累计贡献率超过80%的前m(m≤p)个主成分[20]。根据主成分贡献率以及各主成分中不同影响因素的主成分载荷，即可确定所研究问题中的主要影响因子。

3 结果分析

3.1 珠江河口区咸潮上溯基本特性

3.1.1 咸潮上溯的年际变化

图2(a)为磨刀门水道1998—2015年平岗泵站的咸度超标时数统计结果以及变化趋势。从图中可以发现，平岗泵站咸度超标时数的最高年份出现在2009—2010年枯水期，高达1 864 h，最低年份出现在2002—2003年枯水期，仅为1 h，多年平均超标时数达到了751.4 h，意味着枯水期内平均长达31 d无法正常取水。高超标时数的时间点大致对应于强咸潮事件的年份，例如2005—2006年、2009—2010年以及2011—2012年枯水期的咸度超标持续时间都很长，均在1 500 h(约63 d)以上，是1998—2015年发生的主要几次强咸潮事件[17, 21]。

由图2(b)可见，根据累积距平曲线的拐点将咸度超标时数划分为两个阶段，划分后的结果见图2(a)。可以看出，平岗泵站咸度超标时数的两个阶段有明显的变化特征，具体表现为：第一阶段的咸度超标时数整体表现为“V”型变化，即前5年咸度超标时数逐渐下降，后2年逐步上升，多年平均超标时数为365.8 h。第二阶段咸度超标时数的线性回归结果表明，咸度超标时数总体上存在下降趋势，下降幅度约为63 h/a；多年平均超标时数达961.7 h，较第一阶段的多年平均超标时数大幅增加。

(a)统计结果以及变化趋势

(b)累积距平曲线

3.1.2 咸潮上溯的月际变化

图3为磨刀门水道平岗泵站1998—2015年枯水期不同月份内的咸度超标情况。从图3(a)～(f)可见，高超标时数主要集中在每年12月至次年2月，1月的咸度超标时数最大，为214.8 h。平岗泵站出现咸度超标的时间大约在枯水期开始后的1～2个月，而咸潮完全消退，即咸度低于标准值(250 mg/L)的出现时间通常比枯水期的结束时间提前约1个月左右，表明平岗泵站出现咸度超标现象的时间约为3～5个月。由图3(c)～(e)可以看出，平岗泵站12月至2月的咸度超标时数表现出了与年际变化一致的过程，均以2004—2005年枯水期为时间节点表现出了明显的两阶段特征，其中，2004—2005年枯水期以前的咸度超标时数主要呈“V”型变化；而2004—2005年枯水期以后的咸度超标时数则表现为下降趋势，且1月的下降幅度较大，约为26 h/a。

此外，值得注意的是，对比1998—2015年枯水期3场重大咸潮事件的咸度超标时数统计结果(表1)，发现2009—2010年枯水期咸度超标的出现时间较早，10月的咸度超标时数已达334 h，且持续至次年3月仍有较为严重的咸度超标现象发生，3月的咸度超标时数为267 h；而2005—2006年枯水期以及2011—2012年枯水期出现持续性咸度超标的时间较晚，且高超标时数普遍集中在12月至次年2月，咸潮完全消退的时间相对提前。

(a)10月

(b)11月

(d)1月

(e)2月

(f)3月

表1 1998—2015年枯水期3场重大咸潮事件平岗泵站咸度超标时数统计 h

3.2 关键因素分析识别

根据磨刀门水道的咸潮特性并结合前人的研究结果[22-24]，选取与逐月咸度超标时数对应的枯水期上游来水(马口+三水)、马口及三水以上西北江流域的降雨、海平面、灯笼山的最高和最低潮位以及平均潮差等6个影响要素作为研究对象，采用Spearman秩相关系数对各影响要素与咸度超标时数之间的相关性进行分析，同时采用主成分分析方法量化各影响要素的贡献程度(以百分比表示)，从而确定影响磨刀门水道枯水期咸潮上溯变化趋势的关键因素。

3.2.1 Spearman秩相关分析结果

表2为1998—2015年枯水期咸度超标时数与各影响要素的Spearman秩相关分析结果。

表2 影响要素与咸度超标时数的Spearman相关系数

注：*表示通过95%显著性水平检验；**表示通过99%显著性水平检验。

从表2可以看出，各影响要素与咸度超标时数均表现为负相关关系，其中流量与咸度超标时数的相关性最大，最低潮位、降雨量和海平面次之，都通过了99%置信水平检验，表明这4个要素与咸度超标时数关系密切；而平均潮差和最高潮位与咸度超标时数的相关性较小，未能通过显著性水平检验。

3.2.2 基于主成分分析法的量化分析

为进一步量化不同影响要素的贡献程度，本文对所选取的影响要素进行主成分分析，各主成分的特征值、主成分贡献率及累计贡献率见表3。

表3 特征值、主成分贡献率及累计贡献率

由表3可见，第一、第二主成分特征值大于1，第三主成分特征值接近于1，且前3个主成分累计贡献率大于80%，说明前3个主成分已具备了反映整体信息的能力，故选取特征值为2.539、1.369以及0.964，累计贡献率达81.2%的第一、第二和第三主成分作后续的主成分载荷计算，结果见表4。

表4 主成分载荷计算结果

综合表3、表4可知，第一主成分的贡献率为42.3%，各影响要素中以流量、海平面和最低潮位占有的主成分载荷较大，说明三者与第一主成分的关系密切，在第一主成分中处于同等重要的位置；第二主成分的贡献率为22.8%，与降雨量有较大的正相关，表征降雨因素对咸潮活动的影响；第三主成分的贡献率为16.1%，以平均潮差的主成分载荷最大，反映了平均潮差的动力特征。主成分分析结果表明，咸潮上溯活动受流量、海平面和最低潮位的影响较大，是影响磨刀门水道枯水期咸潮活动变化趋势的关键要素，这与Spearman秩相关分析的结果是一致的。

4 讨论

上节分析结果表明，磨刀门水道近20年来的咸潮上溯活动呈现较为明显的两阶段变化，流量、海平面和最低潮位是影响咸潮活动强弱的主要因素，这与孔兰等[25]的研究结论一致。图4为1998—2015年咸度超标时数、降雨和流量的逐月变化趋势对比。从图4可以看出，高咸度超标时数始终发生在上游来水较枯的月份，表明上游来水对咸潮上溯活动的变化起重要作用。但进一步对比流量和降雨的变化过程可以发现，2003年末至2008年初的降雨变化过程与流量过程并不一致，考虑到降雨是珠江流域径流的主要来源，因此推测引起降雨与径流变化过程不一致的原因是径流受到了人类活动的影响。经调查发现，自2005年起水利部珠江水利委员会每年均会组织实施不同程度的枯季水量调度，其中以2004—2005年以及2005—2006年枯水期的调水压咸活动规模最大。这两个枯水期期间，珠江口磨刀门水道遭受了严重的咸潮灾害，导致珠海、中山和澳门地区的供水受到严重的影响。为缓解地区严峻的供水形势，水利部珠江水利委员会利用上游的骨干水库进行联合调度，加大下泄流量，有效地保障了地区供水安全[26-27]。这一系列的人类活动极大地改变了径流的分配过程，导致2003—2008年月降雨量与流量的变化过程不一致。由此可以看出，磨刀门水道枯水期咸潮的两阶段变化中，引起第一阶段咸度超标时数呈“V”型变化的主要原因是自然因素(降雨)的影响；而引起第二阶段咸度超标时数呈下降趋势的原因除了受到自然要素的影响外，还受到了人为活动的干扰，这是导致两阶段磨刀门水道枯水期咸潮上溯活动变化趋势不同的原因。

图4 1998—2015年咸度超标时数、降雨和流量的逐月对比

孔兰等[25]的研究表明，珠江口咸潮与海平面呈正相关关系。关帅等[28]进一步论证了海平面上升与咸潮强弱的相关关系：海平面上升增强了咸潮上溯的潮汐动力，导致咸潮上溯距离整体向上游移动，两者呈正相关关系。因此，当海平面上升时，咸度超标时数理应有所增加。但本文上节的相关性分析结果表明，海平面上升与咸度超标时数呈负相关，与理论结果相反。图5为2003—2015年枯水期广昌泵站和平岗泵站的咸度超标时数对比。从图5中可以看出，广昌泵站和平岗泵站咸度超标时数的年际变化并不完全一致，如广昌泵站2004—2005以及2010—2011年枯水期期间的咸度超标时数均比上一年增加，但平岗泵站对应枯水期的咸度超标时数却比上一年下降。考虑到平岗泵站作为珠海市的重要取水点，是调水压咸确保的重点对象，而广昌泵站位于平岗泵站的下游且不作为取水点，故上游调水无须保证其咸度低于标准值。因此，平岗泵站的超标时数过程与广昌泵站的不一致是受到了人类活动的干扰。正如前文所提，水利部珠江水利委员会自2005年起实施了一系列的调水压咸活动，保证了平岗泵站的供水安全，人为地降低了平岗泵站的咸潮上溯强度，这造成的直接影响是咸度超标时数与海平面上升要素的正负相关性与理论存在不一致。

图5 2003—2015枯水期广昌泵站和平岗泵站的咸度超标时数对比

综上所述，1998—2015年咸潮上溯活动的变化与枯水期上游来水的增加密切相关，但其中也有人为因素的影响，如大规模的调水压咸活动。但从整体上看，磨刀门水道近10多年来的咸潮灾害仍然十分严重，因此在无法断定自然因素还是人为因素导致其咸潮活动近几年出现下降趋势的前提下，科学的调水压咸仍旧是保障磨刀门河口地区供水安全的重要举措。此外，需要注意的是，海平面上升对咸潮上溯的影响有增强的趋势。受数据资料的限制，因人工采砂、河道疏浚等工程引起的河口地形变化等因素对咸潮上溯变化的影响还有待进一步的研究。