工业革命以来, 化石能源得到了大规模的开发和使用, 随着经济的飞速增长, 全球温室气体含量迅速增加^[1]。据世界气象组织2022年10月26日发布的《2021年全球温室气体公报》, 2021年CO₂气体含量突破历史记录, 达到(415.7±0.2)×10^–6, 增幅约为2.5×10^–6, 略高于近10 a的平均增幅2.46×10^–6。与此同时, 全球主要温室气体浓度相继突破历史记录, CH₄和N₂O气体浓度达新高: CH₄质量分数达到(1908±2)×10^–9, 增幅为1.8×10^–8, 是1983年系统记录以来最大增幅; N₂O质量分数达(334.5±0.1)×10^–9, 增幅达1.3×10^–9[2]。温室气体浓度的迅速增长造成全球气候的改变, 在全球气候变暖的影响下, 全球海洋发生热膨胀效应以及冰川、冰盖减少的现象, 最终造成全球海平面整体呈现上升趋势^[3-4]。

海平面变化高度作为气候变化的指标之一, 其对全球气候变化响应能够反映出气候变化的整个过程^[5], 积累了海内外学者以及研究人员的大量成果。通过全球各海域的潮位站提供的1990—2000年海平面高度数据, Miller等^[6]对全球海平面变化趋势进行了研究, 分析表明海平面高度的上升速率为1.7±0.4 mm/a^[6]。其他学者的研究结果也表明了1900—1979年全球平均海平面上升速率为1.75±0.5 mm/a^[7], 1880—1980年间为1.8 mm/a^[8]。联合国政府间气候变化专门委员会(IPCC)在2007年报告指出, 20世纪全球海平面高度上升速率达到1.7 mm/a, 1993—2010年的平均上升速率接近1901—2010年的2倍^[9]。利用多源卫星数据, 李大炜等^[10]对全球海平面的变化趋势进行了分析和研究, 研究结果显示1993—2011年间全球平均海平面高度上升速率达到了3.12±0.4 mm/10 a^[10]。

海平面高度的上升将对人类社会以及生态系统等会造成巨大的损失。海平面上升将造成海水盐度的增加, 逐渐淹没湿地, 预计到2080年, 因海平面上升将造成全球约五分之一的湿地退化。同时海平面高度的上升会造成风暴潮灾害次数增加以及强度增大。海平面高度的变化与人类生活密不可分, 近海的生物、生态环境及人类自身的活动等将会受到海平面高度变化的影响^[11]。此外, 海平面升高还会加剧其他海洋灾害发生次数及灾害强度, 如海岸侵蚀、海水入侵威胁, 破坏沿海地区的经济发展、人类的生存环境^[12-13]。

中国的海域广阔, 海岸带地势较低, 因此海平面高度的上升将带来十分深刻的影响^[14-15], 如海水入侵沿海低地, 大面积淹没潮滩、湿地, 损害土地以及旅游、生物等资源^[16]。因此研究海平面历史未来的变化规律, 有利于正确地评估灾害以及及时对人类活动发出预警, 以此减少海平面升降带来的影响和损失。基于上述回顾可知, 通过研究历史及未来高排放情景下海平面高度的变化, 掌握海平面变化的规律迫在眉睫。

根据世界气候研究计划(World Climate Research Programme, WCRP), 海平面变化存在全球尺度和区域性多尺度的变化特征, 且区域海平面变化幅度高于全球平均海平面变化幅度。因此我们在关注全球尺度海平面变化的同时, 注重对区域海平面变化的研究^[17-18]。本文基于欧洲中期天气预报中心(European Centre for Medium-Range Weather Forecasts, ECMWF) ERA5数据集提供的1993—2023年历史数据以及2023—2050年SSP5-8.5高排放情景下的海平面日均高度, 对中国海域前后近30 a的海平面高度变化进行研究。以往研究中, 对海平面变化的研究常采用线性拟合方法。线性拟合方法便于实现和理解, 但在非线性数据上的表现较差, 同时非常敏感于异常值。本文采用Mann-Kenddall检验法(MK方法), 部分数据缺失不会对结果造成影响, 不受少数异常值的干扰, 适用性强。不但可以检验时间序列的变化趋势, 还可以检验时间序列是否发生了突变。

1 数据来源及方法 1.1 研究区域及数据来源

本文研究区域为中国近海海域(3°~45°N, 108°~131°E), 研究时间范围1993—2050年。此次研究的数据来源于ERA5数据, 由哥白尼气候变化服务中心生成。该数据集基于水动力模型, 提供了全球海平面相关变量的时间序列, 包括1950年至2050年的潮汐、风暴潮和海平面上升。历史时期(1993—2023年)的模拟数据空间分辨率为0.25° × 0.25°, 时间分辨率为1 d; 未来时期(2023—2050年)基于使用高排放情景SSP5-8.5的全球气候预测, 仅潮汐模拟(MSL)的平均水位, 数据空间分辨率为0.1° × 0.1°, 时间分辨率为1 a。在ERA5海平面数据集中, 海平面数据并非实际海面高度, 而是以1993—2012年共20a的全球海平面平均高度为参考, 计算得到的海平面相对高度。ERA5提供的数据为评估全球海平面变化、沿海洪水和海岸侵蚀的研究提供了基础, 对目前缺乏详细建模研究的数据稀缺区域具有重大意义。

1.2 研究方法 1.2.1 MK趋势检验

本文采用MK变化趋势分析法研究平均海平面高度变化趋势, 得到历史海平面高度变化趋势以及海平面高度的突变年份, 具体方法简述如下:

(1) 构建1993—2050年海平面高度的时间序列(x₁, x₂, ..., x_n), 假设H₀表示该时间序列中的数据是独立的, 随机变量同分布的样本, 没有线性变化趋势。

(2) 计算检验的统计变量S:

$ S = \sum\limits_{i = 1}^{n - 1} {\sum\limits_{j = i{\text{ + }}1}^n {{\rm sgn} ({x_j} - {x_i})} } , $

(1)

其中sgn(x)为阶跃函数,

$\operatorname{sgn}(x)=\left\{\begin{array}{lll} ~~~-1 ~~~~, ~~~~ \text { 若 } x>0 \\ ~~~~0~~~~~~~, ~~~~ \text { 若 } x=0 . \\ ~~~~1~~~~~~, ~~~~\text { 若 } x<0 \end{array}\right. $

(2)

(3) 对于有10个或者更多数据点的时间序列, 使用正态近似。统计值S的方差计算如下:

$ {\rm var} (S) = \frac{1}{{18}}[n(n - 1)(2n - 5) - \sum\limits_i^n {{e_i}({e_i} - 1)(2{e_i} - 5)} ], $

(3)

其中海平面高度相同的样本的个数记为i, i值相同的样本中降水量不同的样本的个数记为e_i。

(4) 当n≥10, 则统计量S近似服从正态分布, 其正态分布的检验统计量Z:

$ Z_S=\left\{\begin{array}{lc} \frac{S-1}{\sqrt{\operatorname{vars}(S)}} & , \text { 若 } S>0 \\ 0 & , ~~\text { 若 } S=0 . \\ \frac{S+1}{\sqrt{\operatorname{vars}(S)}} & , \quad \text { 若 } S<0 \end{array}\right. $

(4)

(5) 在一定的α置信水平上, 如果Z≥|Z₁-α/2|则拒绝原假设, 即在α置信水平下海平面高度时间序列存在明显的趋势。

(6) 海平面高度变化趋势的大小用β表示, 计算如下:

$ \beta=\operatorname{Median}\left(\frac{x_k-x_j}{k-j}\right), \forall j<k . $

(5)

其中Median函数能够返回给定数值的中值; 若β>0, 表示呈海平面高度呈上升趋势; 若β<0, 表示海平面高度呈下降趋势^[19-20]。

1.2.2 MK突变检验

根据海平面高度的时间序列(x₁, x₂, $ \cdots $, x_n), 然后构造一个秩序列:

$ {s_k} = \sum\limits_{i = 1}^k {{r_i}} ,k = 2,3, \cdots ,n, $

(6)

其中,

$ {r}_{i}=\{\begin{array}{l}1,{x}_{i}>{x}_{j},\\ 0,{x}_{i}\le {x}_{j}\text{，}\end{array}\text{j}=1,2,\cdots ,i . $

(7)

定义统计量:

$ U{F_k} = \frac{{[{s_k} - E({s_k})]}}{{\sqrt {{\rm var} ({s_k})} }},k = 1,2, \cdots ,n, $

(8)

其中, s_k的均值和方差分别为E(s_k), var(s_k), 时间序列独立分布的情况下, 其表达式如下:

$ \left\{\begin{array}{l} E\left(s_k\right)=\frac{k(k-1)}{4}, \\ \operatorname{var}\left(s_k\right)=\frac{k(k-1)(2 k+5)}{72}, k=2,3, \cdots, n . \end{array}\right. $

(9)

将构建的时间序列逆序排列, 再通过式(8)求出UF_k, 同时令UF_k=UB_k(k=n, n–1, $ \cdots $1), UB₁=0。得到统计序列UF_k和UB_k, UF_k的正负表明x_k(海平面平均高度)的变化趋势, 正值表明序列呈上升趋势, 即海平面高度呈上升趋势, 负值表明序列, 即海平面高度呈下降趋势; 当UF_k和UB_k超过显著性临界直线时, 说明序列的变化趋势显著。若出现在临界线之间, 则交点对应的时刻为海平面高度发生突变的时间^[21]。

2 数据验证

为保证研究数据的准确性, 本文选取了Aviso网站的卫星数据对使用数据进行验证。法国Aviso(Archiving, Validation and Interpretation of Satellite Oceanographic)提供的格网化数据(MSLA)。该数据融合了T/P, Jason-1, ERS和ENVISAT等多颗卫星的测高资料. 数据时间分辨率7 d, 空间分辨率为1°×1°, 时间范围为1993年1月至2023年7月, 数据经过了必要的标准改正, 如电离层延迟改正、对流层干湿分量改正、固体潮和海潮改正、海洋负荷潮汐改正、极潮改正、电磁偏差改正、仪器改正和反变压计改正等。国内外大量有关海平面变化的研究采用了Aviso网站提供的卫星高度计数据: Soumya等^[22]利潮汐计站获得的每月海平面异常数据, 合并网格化了Aviso网站中的SLA海平面异常数据, 研究了中国南海海域的海平面年际变化。此外, Feng等^[23]利用1900—2013年潮汐计和1993—2013年网格卫星高度计产品Aviso对西北太平洋观测到的季节性海平面变化进行了研究。

为验证ERA5数据的准确度^[24], 采用技术分数(skill score, SS)、均方根误差(root mean square error, rmse)和相关系数(correlation coefficient, CC)共三个参数进行统计分析^[25], 分别用K、C、R表示其计算公式如下所示。

$K=1-\frac{\Sigma\left(X_{\text {ERA5 }}-X_{\text {Aviso }}\right)^2}{\sum\left(X_{\text {ERA5 }}-\bar{X}_{\text {Aviso }}\right)^2}, $

(10)

$ C = \frac{{\sum \left( {{X_{{\text{ERA5}}}} - {{\overline X }_{{\text{ERA5}}}}} \right)\left( {{X_{{\text{Aviso}}}} - {{\overline X }_{{\text{Aviso}}}}} \right)}}{{\sqrt {\left[ {\sum {{\left( {{X_{{\text{ERA5}}}} - {{\overline X }_{{\text{ERA5}}}}} \right)}^2}\sum {{\left( {{X_{{\text{Aviso}}}} - {{\overline X }_{{\text{Aviso}}}}} \right)}^2}} \right]} }} , $

(11)

$ R = \sqrt {\frac{{\sum {{\left( {{X_{{\text{ERA5}}}} - {X_{{\text{Aviso}}}}} \right)}^2}}}{N}} . $

(12)

其中, X_ERA5和X_Aviso分别为ERA5和Aviso 2011—2020年海平面高度的月平均数据; $ {\overline X _{{\text{ERA5}}}} $和$ {\overline X _{{\text{Aviso}}}} $分别为ERA5数据和Aviso数据的算术平均; N为统计变量的个数。K>0.65表示数据精度非常好, 0.5~0.65表示很好, 0.2~0.5表示好, K<0.2则表示差。

本文选取了经纬度范围38°~39°N, 120°~121°E, 2011年1月至2020年12月连续120个月的月平均海平面卫星高度计数据进行验证, 验证结果见图 1, ERA5数据与卫星数据基本一致。经计算验证, ERA5数据技术分数大于0.98, 均方根误差为0.017, 两种数据的相关系数达到了0.99, 说明ERA5数据集较可靠, 与前人结论一致^[26], 可以用于后续的计算。

3 研究结果 3.1 突变分析

通过前述方法得到海平面高度的年均高度变化趋势, 统计参数见表 1, 海平面年均高度随年份的变化趋势见图 2。

利用MK检验法研究分析海平面高度的1993—2050年的时间序列, 在0.05的显著性水平下, Z≥|Z_1-α/2|拒绝原假设具有统计学意义。α=0.05时, 临界值为1.96, 拒绝原假设说明海平面年均高度存在显著的变化趋势; α=0.01时, 临界值等于2.58, 拒绝原假设说明海平面年均高度存在非常显著的变化趋势。在1993—2023年海平面高度的MK趋势检验统计参数中, β=2.33>0, 说明该时间段平均海平面高度呈现上升趋势, 1993—2023近30年海平面高度的MK趋势检验值Z=12.23>2.58, 通过置信度99% 的显著性趋势检验, 因此, 1993—2023年海平面高度呈显著性上升趋势。对2023—2050年未来SSP5-8.5高排放情景下平均海平面进行估计, β=5.76>0, 海平面高度变化呈上升趋势, 2023—2050年SSP5-8.5高排放情景下海平面高度的MK趋势检验值Z=8.01>2.58, 通过了置信度99%的显著性趋势检验, 因此2023—2050年未来SSP5-8.5高排放情景下海平面高度呈现十分显著的上升趋势。

1993年至2050年平均海平面高度见图 2。计算得到1993年至2050年的平均海平面的增长速率, 1993年至2023年增长速率为2.33 mm/a, 海平面变化整体上升; 2023年至2050年中国海域海平面增长速率为4.98~6.92 mm/a, SSP5-8.5高排放情景下平均海平面高度上升趋势相对较快。

通过前文叙述的MK突变检验方法, 绘出的UF_k和UB_k曲线图如图 3所示。

图 3中黑色虚线之间的区域表示0.05显著水平的区间, 若交点位于该区间内, 则表示该年份的海平面上升具有突变性, 若交点位于区间外, 则表示该年份海平面上升不具备突变性。从1993—2050年海平面高度MK检验的UF趋势线中可以清晰看出海平面高度变化趋势: 1993—1996年, UF值小于零, 海平面高度呈下降趋势; 1996—2023年, UF值大于零, 海平面高度持续增长, 且在1998年之后, UF的值位于0.05显著性水平之上, 通过了0.05显著性检验。由图 3(a)可知, 2009年UF和UB曲线产生交点, 但该交点没有位于检验范围内, 说明该交点没有通过0.05的检验, 所以该年份的海平面高度的上升不具有突变性。2024—2050年, UF值大于零, 海平面高度持续增长, 且在2028年之后, UF的值位于0.05显著性水平之上, 通过了0.05显著性检验。图 3(b)可以看出, 2038年UF和UB曲线产生交点, 但该交点没有在置信水平区间[-1.96, 1.96]内, 说明该交点没有通过0.05的检验, 所以该年份的海平面高度的同样上升不具有突变性。

3.2 区域分析

1993—2023年及2023—2050年中国海域平均海平面高度如图 4所示。根据先前研究可知, 黄渤海区域在过去的19 a中整体呈现上升趋势, 且有季节性和地域性的特点。在过去的30 a中, 渤海区域辽东湾与渤海湾附近平均海平面高度较高, 黄海海域及渤海其他区域较低, 在未来SSP5-8.5排放情景下, 黄海区域的海平面高度普遍高于渤海区域, 其中海平面高度的峰值出现在黄海南部区域, 黄渤海区域在东南向西北方向呈现出两侧海平面高度高, 中间较低的趋势; 东海海域的海平面高度整体呈现从北向南逐渐升高的分布特点, 峰值出现在台湾岛东北部区域。根据相关研究可知, 东海海平面高度同样呈现明显的季节变化, 呈现冬、春季为海平面低值季节, 夏、秋季是高值季节的特点。SSP5-8.5高排放情景下海平面平均高度整体有所上升。南海海域的海平面高度整体同样呈现从北向南逐渐升高的分布特点, 峰值出现在南沙群岛区域。SSP5-8.5高排放情景下南海海平面平均高度整体明显升高, 峰值区域由原来南沙群岛区域北移至东沙群岛附近。

南海海平面高度明显高于黄渤海区域及东海区域。海平面高度呈现自西北向东南递增的分布特征, 在未来SSP5-8.5高排放情景下, 海平面上升幅度总体上呈现东海最高, 南海次之, 黄渤海最小的分布特点。

3.3 季节分析

中国海域区域各季节平均如图 5所示。由图 5(a)可以看出, 渤海春季辽东湾区域海平面高度最高, 渤海湾区域次之, 莱州湾区域平均海平面高度最低; 黄海区域南部及北部区域平均海平面高度较低, 中间区域较高, 山东半岛南部区域平均海平面高度最高; 东海区域由北向南平均海平面高度由北向南逐渐升高, 浙江、江苏附近海域平均海平面高度较低, 台湾岛东部海域平均海平面较高, 在其东北方向出现峰值; 春季南海区域平均海平面最高值出现在南海中心区域, 由该区域向四周逐渐降低, 最低值出现在海南及福建附近海域。

夏季渤海海域海平面高度由西北向东南方向逐渐降低; 夏、秋两季黄海平均海平面高度由北向南逐渐降低, 山东半岛南部海域同样出现峰值, 相同范围内, 夏秋两季平均海平面高度高于春季; 春季夏季平均海平面高度由东海中心区域向四周逐渐降低, 福建附近海域平均海平面高度最低, 夏季南海最高值出现在南海东部海域同时呈现出由东向西海平面高度逐渐降低的趋势。

秋季沿岸地区海域海平面高度较高, 离岸较远的区域海平面高度相对较低。秋冬两季除辽东湾部分区域, 随着纬度的上升, 渤海海域海平面高度变低; 秋季在南部沿海城市附近海域出现多处峰值, 除以上靠近陆地海域, 东海由北向南海平面高度逐渐降低, 其中台湾岛附近海域最低; 秋季平均海平面高度峰值同样出现在海南及广东附近海域, 整体呈现由西北向东南方向逐渐降低的趋势。

夏、秋、冬三个季节中, 渤海北部平均海平面高度低于渤海南部。渤海同一区域内, 夏秋两季平均海平面高度高, 冬春两季平均海平面高度低; 其中, 相同范围内夏季海平面高度最高, 冬季最低。冬季黄海区域由北向南平均海平面高度逐步上升, 同时在山东半岛南部区域及靠近朝鲜半岛海域出现峰值。相同季节内, 黄海平均海平面高度高于渤海海域。冬季福建附近海域平均海平面高度较高, 中国浙江、江苏以及日本九州岛附近海域平均海平面高度较低, 整体呈现由北向南逐渐升高的趋势。相同区域内, 东海海域夏秋两季平均海平面高度较高, 冬季次之, 春季最低。与黄渤海区域相比, 春、夏、秋三个季节平均海平面高度相差不大, 但冬季东海区域平均海平面高度高于黄渤海区域。冬季南海平均海平面高度整体由东向西逐渐降低。相同区域内, 南海海域整体平均海平面高度春季最高, 夏季次之, 秋冬两季最低。与黄渤海及东海海域相比, 春季南海平均海平面高度高于以上海域, 秋季略低于其他海域, 夏季冬季与以上海域相差不大。

4 结论

置信度99%的显著性趋势检验, 海平面高度呈现十分显著的上升趋势, 但是上升不具备突变性。1993—2023年海平面上升速率约为2.33 mm/a; 2023—2050年平均海平面高度上升趋势相对较大, 增速在4.98~6.92 mm/a。

(2) 除个别年份外, 1993—2050年海平面整体呈上升趋势, 中国海域海平面高度呈现自西北向东南递增的分布特征。对比历史数据, 未来SSP5-8.5高排放情景下, 海平面上升幅度总体上呈现东海最高, 南海次之, 黄渤海最小的分布特点。

(3) 中国海域海平面变化具有显著的季节性变化, 黄渤海以及东海海域夏秋两季上升幅度较大、冬春两季幅度小, 南海海域春夏两季上升幅度大、秋冬两季幅度。