基于REOF的珠江三角洲河网区水位时空变异特征分析

杭宜铖，范敏韬，刘智勇

(1.中山大学土木工程学院，广东 广州 510275； 2. 中山大学水资源与环境研究中心,广东 广州 510275；3.南方海洋科学与工程广东省实验室(珠海), 广东 珠海 519082)

国内外一些学者从不同角度对全球不同区域的水位进行了研究。Männikus等[1]研究了波罗的海东部的拉脱维亚海岸的平均水位、季节性水位和极端水位的变化，计算了不同时间段内不同水位出现的经验概率分布并使用标准的Kolmogorov-Smirnov检验作为显著性检验，根据每周和每月的水质研究了平均和极端水位的季节性变化并对趋势斜率进行了评估。Mudersbach等[2]研究了海平面的升高是否显著影响德国库克斯港的水位变化，分析了不同百分位海平面时间序列的线性趋势并将其与全球平均海平面趋势进行比较，研究了第99个海平面百分位序列与长期记录的库克斯港的全球平均海平面之间的非线性趋势差异，最后分析了长期记录的库克斯港的潮汐成分和非潮汐残差(波动)的变化。

珠江三角洲位于广东省中南部、珠江下游，是由珠江水系的东江、西江、北江及其支流组成的复合型三角洲，是世界上水系最复杂的三角洲之一[3]。珠江三角洲有着“三面环山、一面临海，三江汇合、八口分流”的独特地形地貌，河网区河道密集，河网密度达0.68～1.07 km/km2[4]。Zhang等[5]研究了珠江三角洲地区年最高水位和最低水位的时空变化规律，采用Mann-Kendall趋势检验和Pettitt法检验趋势和突变点，采用去趋势预置白方法(trend-free pre-whitening)消除数据的自相关性的影响,认为水位极值的变化趋势在不同地区表现出不同的特征，总体上珠江三角洲上游水位呈下降趋势，而中下游水位呈上升趋势。陈晓宏等[3]认为珠江三角洲河网区水文特征变异主要是由人类活动造成的。陈晓宏等[6]基于区域化变量理论和克里格方法分析珠江三角洲水位的空间特征，发现水位间的相关性与水系走向有关，丰水期存在13.02 km的相关半径，超过这个距离，水位变化的相关性较低。唐亦汉等[7]分析了珠江三角洲的高水位集中期特征，得出了珠江三角洲入口及口门站点的高水位集中期比河网区站点的时间更长的结论。杨清书等[8]通过计算得出珠江三角洲水位总体呈上升趋势，平均变化率为3.1 mm/a。蒋陈娟等[9]应用功率谱分析和正、反傅氏变换方法对河网区月均水位序列的周期变化和趋势变化进行了研究，发现水位序列中包含0.2～11.7 a多尺度的显著周期波动，水位变化率从-3.8～3.7 mm/a不等，珠江三角洲河网区水位总体呈上升趋势。时翠等[4]运用累积距平法、趋势分析法等对1月份和7月份的水位时间序列进行突变点分析，结果表明在1974—1982年和1989—1995年月平均高水位发生变异点集中现象，且小部分潮位站存在2个变异点，大多数仅存在1个。谢凌峰等[10]研究了人类采沙对珠江三角洲河网区河性演变的影响，结果表明在人类大规模采沙影响下，河网区河床下切严重，来流量变化不大但来沙量大幅下降，同流量下水位下降明显。孔兰等[11]利用Spearman相关分析对珠江三角洲5个代表站的水文年内变化规律进行研究，结果表明珠江口水位年内变化不均匀系数和集中度呈现由珠江口深处向近河口处减少的趋势。

以往学者对珠江三角洲的水位研究，大多是分别分析水位的时间尺度和空间尺度的特征，但是从空间分型的角度来分析珠江三角洲水位演变趋势的研究还鲜有报道。本文通过空间分型识别出珠江三角洲河网区存在的各水位要素场，并采用旋转经验正交函数法(rotated empirical orthogonal function, REOF)分析其时间变异特征，以期为珠江三角洲不同地区的水情监测、防洪防涝措施的制定、水资源的分配与调度提供参考。

1 数据来源与研究方法

1.1 数据来源

收集1962—2016年珠江三角洲河网区12个水文站的月平均水位数据，水文站分布如图1所示。所有站点的数据完整，通过三性审查。资料来源于广东省水文年鉴《珠江流域水文资料》和广东省水文局。

图1 研究区域及水文站分布

1.2 研究方法

1.2.1REOF分析

经验正交函数(empirical orthogonal function, EOF)分解方法最初由Pearson[12]于1901年提出，该方法可将空间上随时间变化的要素场分解为正交的空间函数和时间函数两部分，用方差最大的几个时间函数与其对应的空间函数乘积，以此结果描述要素场的主要特征。Horel[13]指出如果第一个EOF在其域上具有恒定的符号，则第二个EOF通常具有两个符号且零线穿过第一个EOF的最大值。这导致了域依赖性和非局部性问题[12-15]，由此产生的一系列问题使模型产生的结果变得难以解释[16-18]。

为了解决这些问题，研究人员基于旋转的EOF的线性变化，提出了REOF的概念[13,19-20]。REOF分析是基于EOF分析，选取多个解释方差达到一定标准的特征向量进行正交旋转转换，旋转后的空间模态结构更加清晰，可以更好地体现地域差异，便于将水文要素一致的地区划分为同一区域。相比于EOF方法，REOF在局部模态提取方面的准确性和有效性都有大幅提升，提取出来的空间模态具有更好的物理意义[21]。REOF取样误差小，分离出的空间模态可以有效地体现要素场在空间分布上的相关性[22]，在反映局部变量之间关系和分布方面有着更明显的优势[23]，能更好地反映水位的时空分布特征。

珠江三角洲地区河网密布交错，同时又受到涨潮的影响，虽然每个水文站受到上游来水影响不一，但来自海相的潮汐影响比较一致，因而河网内的水文站具有一定的空间联系和连续性，比如都受到珠江三角洲河口的涨潮影响。基于这样的前提和假设引用REOF分析来探讨其空间模块。本文使用反距离插值法对该区域河道水位进行均匀插值，不考虑陆地区域，仅对河道内空间模态进行特征分析。

1.2.2Mann-Kendall检验

Mann-Kendall非参数秩次检验法通常用于评估水文气象时间序列(如流量、温度、降水等)的趋势和突变检验[24]。与参数统计检验相比，非参数检验被认为更适合于水文气象时间序列中经常遇到的非正态分布的数据和缺失数据，且不受少数异常值干扰。关于Mann-Kendall检验中统计值S与标准化检验统计值ZMK的计算方法详见文献[24]。

ZMK是一个服从正态分布的统计量，反映了时间序列数据的变化趋势。本文中，若ZMK>0，表明平均水位数据随着时间的推移呈上升趋势；若ZMK<0，则表明平均水位数据随着时间的推移呈下降趋势。当|ZMK|≥1.96时，表示通过了α=0.05的显著性水平检验，即通过了置信度为95%的显著性检验。

1.2.3多项式曲线拟合

数据拟合是一种重要的数据处理办法，被广泛应用于水文资料的整编中。数据拟合中最常用的方法是多项式曲线拟合。多项式曲线拟合有许多方法，数据量较少时适合使用最小二乘法[25]，其主要原理是保证拟合值和原数据的差的平方和为最小。对方程式求偏导后，使用高斯消元法求得拟合系数，确定拟合系数后得出拟合方程。本文对REOF分析结果所得的时间序列数据进行二项式拟合，可以看出整个研究时段中参数变化趋势的突变。

2 结果与分析

2.1 平均水位空间分布特征

从多年平均水位的空间分布来看，珠江流域下游站点的平均水位由西到东、由北向南递减(表1)。整体递减方向与水系的方向相同，其中东西方向的递减更明显。水位高值区位于研究区域上游，包括紫洞、南华和澜石3个站点，平均水位都在0.5 m以上，其中紫洞站的多年平均水位最高，达到了0.8 m。水位低值区位于东江下游接近入海口处，包括三沙口和黄埔两个站点，平均水位均小于零，其中黄埔的平均水位值最低，达到-0.03 m。旱季和雨季水文站的平均水位与多年平均水位的空间分布及空间变化大致相同。同一个站点旱季和雨季的平均水位差值与该站点的平均水位呈正相关，其中位于高值区的紫洞和南华2个站点旱季雨季平均水位差距明显，差值达0.9 m以上，说明其受上游来水和降雨影响较大。

表1 珠江三角洲水文站平均水位

2.2 REOF分区

对珠江三角洲河网区12个水文站1962—2016年月平均水位进行REOF分析，结果显示前4个空间模态的方差贡献率依次为23.59%、18.23%、10.44%和8.89%。由此可知，第一空间模态集中了珠江三角洲河网区水位要素最主要的信息，前4个空间模态累计方差贡献率为61.17%，基本可以代表珠江三角洲河网区水位要素的时空变化特征。将4个空间模态所对应的时间序列进行月统计及旱季雨季划分，并在旱季和雨季2种时间序列下讨论各个空间模态的时空特征(图2)。

由图2(a)可知，第一空间模态正值中心位于研究区域南部，以白蕉、竹银2个水文站为主，其中白蕉站的载荷向量正值最高，达到0.91，定义为第一空间异常型，结合图3、图4对应的时间序列图，第一空间异常型的旱季平均水位时空分布特征为：1962—2016年站点平均水位有上升趋势。在1962—1964年、1980—1984年2个时间段，水位的平均年际变化率最大。平均水位最低值出现在2005年，最高值出现在1983年。第一空间异常型的雨季平均水位时空分布特征为：1962—2016年站点的平均水位有下降趋势。1991年之前，水位的平均年际变化率较小，1991年之后年际变化率较大。在1962—1964年、1971—1974年、1991—2003年3个时间段，水位的平均年际变化率最大。雨季平均水位最低值出现在1963年，最高值出现在1994年。

(a) 第一模态

由图2(b)可知，第二空间模态的正值中心位于研究区域中部，以南华、竹银两个水文站为主，其中南华站的载荷向量正值最高，达到0.74，定义为第二空间异常型。结合图3、图4对应的时间序列图，第二空间异常型的旱季平均水位时空变化特征为：1962—2016年正值中心区站点的水位变化趋势为先上升再下降，转折点位于2003—2005年之间。1984年之前，水位的平均年际变化率较大，1984年之后平均年际变化率较小。1980—1984年，水位的平均年际变化率最大。平均水位最低值出现在1963年，最高值出现在1983年。第二空间异常型的雨季平均水位时空分布特征为：1962—2016年正值中心区站点的水位变化趋势为先上升再下降，在1962—1964年、1967—1974年、1991—2003年3个时间段，水位的平均年际变化率较大，1974—1991年水位的平均年际变化率较小。平均水位最低值出现在1963年，最高值出现在1973年。

由图2(c)可知，第三空间模态的正值中心位于研究区域的北部，以广州浮标厂、三善滘站为主，其中广州浮标厂的载荷向量正值最高，达到0.92，定义为第三空间异常型。结合图3、图4对应的时间序列图，第三空间异常型的旱季平均水位时空变化特征为：1962—2016年正值中心区站点的水位变化趋势为先上升再下降，转折点位于1991—1993年之间。在1969—1973年、1982—1984年2个时间段，水位的平均年际变化率较大。平均水位最低值出现在1963年，最高值出现在1983年。第三空间异常型的雨季平均水位时空分布特征为：1962—2016年站点的水位变化趋势为先上升再下降，在1962—1964年、1991—2003年这2个时间段，水位的平均年际变化率较大，1974—1991年水位的平均年际变化率较小。平均水位最低值出现在1963年，最高值出现在1994年。

(a) 第一模态

由图2(d)可知，第四空间模态的正值中心位于研究区域东南部，北江下游接近入海口处，以灯笼山站为主，载荷向量正值最高，达到0.91，定义为第四空间异常型。结合图3对应的时间序列图，第四空间异常型的旱季平均水位时空变化特征为：在1962—2016年时间区间正值中心区站点的水位变化趋势为先上升再下降，转折点位于1991—1993年之间。1980—1984年，水位的平均年际变化率较大，1985—1995年平均水位的年际变化率较小。平均水位最低值出现在1963年，最高值出现在1983年。第四空间异常型的雨季平均水位时空变化特征为：1962—2016年时间区间正值中心区站点的水位变化趋势为先上升再下降，转折点位于1986—1988年之间。在1962—1964年、1991—2003年这2个时间段，水位的平均年际变化率较大，1974—1991年水位的平均年际变化率较小。平均水位最低值出现在1963年，最高值出现在1994年。

(a) 第一模态

2.3 水位变化趋势的时空分布

图5给出了1962—2016年珠江三角洲河网区12个站点各月的平均水位序列的Mann-Kendall趋势检验结果，图中各站点按照纬度从高到低排列。研究区域内呈现水位降低趋势的月份主要集中在8—9月，其中有5个站点的水位在8—9月有显著的下降趋势。呈现水位升高趋势的月份主要集中在1—3月和12月，超过一半的站点水位在这几个月有显著的上升趋势。在旱季(1—3月和10—12月)，各月份站点的水位变化趋势以上升为主，其中大部分站点的月平均水位上升趋势显著；在雨季(4—9月)，各月份站点水位的趋势变化中上升和下降趋势相对均匀分布，小部分站点的月平均水位呈显著的上升或下降趋势。由图6可得各站点全年、旱季和雨季平均水位趋势变化的空间分布。在旱季，10个站点的水位有上升趋势，仅有2个站点有下降趋势，其中研究区域的中部、东部，西江和东江下游接近入海口处的7个站点水位上升趋势显著，有下降趋势的站点位于研究区域西北部，且下降趋势不显著。在雨季，有4个站点水位有上升趋势，8个站点水位有下降趋势，水位有上升趋势的站点主要位于研究区域的最东部，东江与西江入海口处，水位有显著下降趋势的6个站点位于北江支流下游和研究区域西北部。全年来看，有5个站点水位有上升趋势，且上升趋势显著；7个站点水位有下降趋势，水位有显著上升趋势的站点主要位于研究区域的最东部，东江与西江入海口处。水位有显著下降趋势的2个站点位于研究区域西北部。下降趋势的出现可能是因为近年来人类大规模采砂，造成珠江三角洲河网区河道明显下切，且部分河道疏浚也加剧了河床的下切，水位随之大幅下降[26]。而表现为上升趋势的站点全部靠近河口，受全球海平面上升的影响较大[27]，加上这些地区人类采砂规模相对较小，前者产生的影响比后者大，河道总体呈淤积状态，水位上升。全年来看，水位呈上升趋势的站点分布与旱季的类似，主要位于东江和西江的下游接近入海口处，且大部分站点的上升趋势显著。水位呈下降趋势的站点位于北江下游接近入海口处和研究区域西北部。大部分站点旱季与雨季的水位变化趋势基本一致。但澜石站和竹银站旱季水位有显著上升趋势，雨季水位有显著下降趋势。白蕉站和灯笼山站旱季水位有上升趋势，雨季水位有显著下降趋势。这类异常现象可能是人类活动与全球海平面上升共同作用，且人类活动在空间分布上具有不均匀性导致的[11]。

图5 各水文站水位Mann-Kendall趋势值分布

(a) 旱季

3 结 论

a.珠江三角洲河网区12个水文站的月平均水位高值区位于研究区域上游，低值区位于东江下游接近入海口处。站点的平均水位由西向东、由北向南递减。整体递减方向与水系的方向相同，其中东西方向的递减更明显。

b.对珠江三角洲河网区的12个水文站的月平均水位的REOF分析表明，前4个空间模态方差贡献率达61.17%，可以基本表现出该区域的水位时空变化特征。其中第一模态是珠江三角洲河网区水位变化的主要特征，方差贡献率达到了23.59%。4个空间模态对应4个空间异常型，分别对应研究区域的南部、中部、北部和东南部。分区结果体现了较强的地域特征，涵盖了整个珠江三角洲河网区的大部分区域。各空间异常型在旱季和雨季表现的水位变化趋势有所区别，具体为第一空间异常型表现为旱季有水位上升的趋势，雨季有水位下降的趋势。第二、三、四空间异常型均表现为旱季和雨季水位有先上升再下降的趋势。

c.Mann-Kendall趋势分析结果表明，水位呈显著上升趋势的月份主要集中1—3月以及12月。水位呈显著下降趋势的月份主要集中在8—9月。水位呈上升趋势的月份占主导。从空间上来看，位于西江与东江下游区域的站点平均水位主要呈上升趋势，位于研究区域西北部和北江下游的站点，平均水位主要呈下降趋势。