全球变暖背景下内蒙古地区沙尘暴频次变化的预估

杨诗妤 闻新宇

1. 北京大学物理学院大气与海洋科学系, 北京 100871; 2. 全球变化与中国绿色发展协同创新中心,北京 100875; † 通信作者, E-mail: xwen@pku.edu.cn

春季沙尘暴是对我国北方自然环境和农业生产, 甚至人体健康带来严重影响的灾害性天气, 其危害表现在以下方面: 1) 造成短时的较强的跨越省际规模的区域性空气污染事件, 或能见度大幅度下降[1]; 2) 降低土壤肥力, 影响农业生产力[2]; 3) 沙尘携带过敏原和细菌微生物, 加速疾病的大面积传播[3]。从全球范围看, 沙尘暴主要发生在干旱、半干旱地区, 北非撒哈拉沙漠、中国西北部、澳大利亚西部、北美西南部是沙尘暴的 4 个主要源地[4],其中塔克拉玛干沙漠和内蒙古戈壁是中国两个主要的沙尘源[5]。塔克拉玛干沙漠地处南疆无人区, 南侧和北侧被昆仑山脉和天山山脉包围, 尽管沙尘活动频繁, 但对社会经济和人类生产、生活的影响很小。虽然内蒙古戈壁地区(包括甘肃河西走廊、宁夏、内蒙古中-西部)的沙尘暴次数少于塔克拉玛干沙漠地区, 但该区域人口众多, 工农业活动强度大,因此内蒙古戈壁沙尘暴造成的影响远大于塔克拉玛干沙漠。不仅如此, 内蒙古地区的沙尘暴除造成本地损失外, 还波及其下游人口更稠密的华北(陕西、山西、京津冀)和东北地区(内蒙古东部), 甚至随着西风急流输送到东北亚和太平洋地区[6]。1995—2008 年, 内蒙古地区共发生沙尘暴 998 次, 其中特强沙尘暴就有 59 次[7]。2006 年 4 月 16 日这一天中,来自内蒙古的沙尘在北京的沉降量高达 33 万吨左右[8]。

全球变暖是中国沙尘暴长期变化的大背景。人类自工业活动以来排放的过量温室气体导致全球温度上升, 进而造成冰冻圈融化、海平面上升和极端天气。未来, 随着全球变暖进一步加剧, 中国沙尘暴会出现怎样的长期变化？对这一基本问题的研究尚不充分。关于降水的长期趋势, Li 等[9]的研究表明, 到 21 世纪末, 如果全球增温 2℃, 中国出现百年一遇极端降水的概率将增加 2.4 倍; Zhang 等[10]指出, 2000—2100 年, 中亚地区帕默尔干旱指数呈下降趋势, 中亚干旱越来越严重。关于沙尘暴的长期趋势, 相关研究十分有限。Wang 等[11]、Shao 等[12]和 Guo 等[13]都指出过去 50 年中国沙尘暴频次在减少, 但未讨论未来的长期变化; Zhu 等[14]根据 1955—2005 年 258 个台站的沙尘暴观测数据和 NCEP/NCAR 再分析资料等数据, 认为贝加尔湖温度升高使蒙古气旋频次减少, 未来中国沙尘暴频次可能减少。总的来说, 针对未来中国沙尘暴长期变化趋势这一基本问题, 现有研究不仅未从沙尘暴频次增减方面进行定性的研究, 更没有降尺度到沙尘暴日数据(daily)的层面进行定量的研究, 而 daily 是沙尘暴最敏感、最关键的时间尺度。因此, 本文旨在通过重建内蒙古地区沙尘暴未来 100 年的 daily 指数, 研究其长期变化趋势。

数值模式和统计模型是获得沙尘暴频次的两大类方法。尽管现有的气候模式对气候变化的模拟有很好的代表性, 但是对沙尘暴这一天气尺度系统的长时间模拟存在很大的不确定性。Pu 等[15]指出,CMIP5 模式模拟的多年平均沙尘光学厚度和观测数据存在 1～2 个数量级的差异, 很难直接使用其输出结果对沙尘暴进行预估。与数值模式相比, 统计模型对沙尘暴这类高频问题的模拟表现更佳。在有限的研究中, 使用过的方法主要有广义线性回归和人工神经网络两种。Huang 等[16]和 Kaboodvandpour等[17]利用风温压湿数据, 获得台站的沙尘暴频次神经网络模型, 并指出神经网络模型的结果优于线性回归方法。在回归因子选择方面, 前人对与沙尘暴相配合的气象要素和环流形势进行了较广泛的研究,局地天气要素包括地面风速、温度和降水等[11,18],高空天气系统包括气旋槽脊[19-20]和极涡[21-22]等,气候要素包括北极涛动[23]和西太平洋遥相关[24]等。本文在前人工作的基础上筛选出 9 个关键因子, 用神经网络方法重建 1860—2100 年内蒙古地区沙尘暴频次日变化趋势, 并尝试对内蒙古地区沙尘暴频次未来长期变化进行预估, 为探究全球变暖背景下中国地区沙尘暴的频次变化、沙尘暴的防治以及相关环境政策的修订提供科学参考。

1 数据与方法

1.1 数据

本文所用内蒙古地区(36°—44°N, 97°—114°E)沙尘暴频次数据来源于中国气象局最新编制的 1954—2007 年沙尘暴台站观测日数据集[25]。该数据集包含全国 490 个台站的数据, 本研究使用内蒙古地区受沙尘暴影响最严重的 115 个台站的数据(图 1 中黑实线框)。本文中沙尘暴日频次为该区域 115 个台站记录次数的总和, 如果同一台站一天中有多个记录, 则只计算一次。内蒙古沙尘暴主要发生在春季[11,26], 沙尘暴年频次为春季 3 个月(3—5 月)日频次的加和。

图1 台站观测的1954—2007 年平均春季沙尘暴频次Fig. 1 Observed annual mean number of springtime dust storms in 1954-2007

沙尘暴频次统计模型中所用 9 个关键因子的数据来源于 1954—2007 年 NCEP/NCAR 再分析资料的日数据集[27], 包括 500 hPa 位势高度、海平面气压、700 hPa 纬向风速、700 hPa 经向风速、700 hPa温度、700 hPa 比湿、250 hPa 纬向风速和 250 hPa经向风速, 数据的水平分辨率为 2.5°×2.5°。

耦合模式比较计划(CMIP5)提供 1860—2100 年历史及未来 9 个关键因子的数据。本文出于对 1860—2100 年 9 个因子 daily 数据可用性的考量, 选用GFDL-CM3, IPSL-CM5A-LR 和 HadGESM-ES 这 3个模式的数据(在所有 CMIP5 模式中, 只有 5 个模式数据完整, 其中 MRI-CGCM3 和 MIROC5 两个模式的数据存在较大的系统性误差, 只有上述 3 个模式可用)。数据的水平分辨率如下: GFDL-CM3 为144×90 个网格点, IPSL-CM5A-LR 为 96×96 个网格点, HadGESM-ES 为 144×92 个网格点。特别地, 因HadGESM-ES 模式的历史试验数据不完整, 1860—2005 年历史试验数据(Historical Runs)只来源于GFDL-CM3 和 IPSL-CM5A-LR 这两个模式。2006—2100 年的未来数据则来源于 3 个模式的低辐射强迫 RCP2.6 情景和高辐射强迫 RCP8.5 情景的试验数据。

1.2 方法

建立统计模型的本质是建立影响沙尘暴的关键因子与沙尘暴频次之间的关系。地表风速大, 湿度小, 易起沙; 高空风速大, 有利于沙尘的远距离传输。蒙古气旋和高空槽脊影响各高度层的风速, 极涡强弱影响蒙古气旋频次。本研究在前人工作的基础上, 凝练与沙尘暴存在物理联系的天气学特征,利用相关系数法, 选出影响内蒙古地区沙尘暴的 9个关键因子: 500 hPa 位势高度(Z500)、海平面气压(SLP)、700 hPa 纬向风速(U700)、700 hPa 经向风速(V700)、700 hPa 温度(T700)、700 hPa 比湿(q700)、250 hPa 纬向风速(U250)、250 hPa 经向风速(V250)和极涡指数(AO)。通过历史上 390 个沙尘暴日距平值的合成分析, 选定前 8 个因子的关键区, 对每个因子 2～3 个关键区的指纹(fingerprint)信息进行线性叠加, 得到前 8 个因子的日数据。关键回归因子的定义和计算方法参见表 1。第 9 个因子是大尺度极涡指数, 在区域平均的 30°—60°N 与 60°—90°N 海平面气压差值基础上, 计算春季平均值。最后, 利用上述 9 个回归因子各自的 1971—2010 年气候态平均值和标准差, 进行标准化(图2)。

图2 3 个CMIP5 模式平均的9 个关键因子1860—2100 年标准化值Fig. 2 Time series of annual mean standardized value of 9 key factors in 1860-2100 derived from three CMIP5 models

表1 关键回归因子的定义与计算方法Table 1 Definitions and computing methods of nine key factors

关键因子与沙尘暴频次之间的统计关系并非简单的线性关系。由简单的神经元与神经元互联组成的人工神经网络模型(BP 神经网络)能较好地拟合出大数据之间的非线性关系[28]。本文使用的 BP 神经网络由输入层(第一层)、隐藏层(第二层)、输出层(第三层)构成(图 3)。其中, 隐藏层可以是一层, 也可以是多层, 每层的节点数按需设置。BP 神经网络的基础是迭代算法, 包括激励传输和权重更新。隐藏层神经元将输入数据加权求和, 再带入激活函数, 将结果传递到输出层。根据误差最小原则, 按梯度下降和误差反向传播算法进行权重更新[29]。当隐藏层只用一层且神经元数趋于无穷时, 可替代任何非线性映射关系, 适用于训练样本不大的情形。本文测试了几种不同的神经网络, 经权衡收敛速度和训练集表现等多种指标后, 使用的神经网络输入层为 9 个关键因子标准化数据; 隐藏层为一层,共有 500 个节点, 激活函数为 ReLU; 输出层为沙尘暴日频次, 对 1954—2007 年的数据做 90%随机训练和 10%随机验证, 共训练 100 轮。

图3 BP 神经网络结构示意图Fig. 3 Topology of BP neural network structure

2 对 1954—2007 年沙尘暴频次的重建与检验

由 9 个关键因子构建的神经网络能很好地再现观测的沙尘暴频次。图 4 对比 1954—2007 年台站观测的沙尘暴频次、BP 模型基于再分析资料计算得到的沙尘暴年频次以及 BP 模型基于 CMIP5 模式历史试验数据计算得到的沙尘暴年频次。

图 4(a)显示, 观测的沙尘暴年频次 1954—1997年为下降趋势, 1998—2001 年为上升趋势, 高值出现在 1958, 1966 和 1976 年, 这种先减少后增加的变化趋势与 Guo 等[13]和 Guan 等[30]的研究结果一致。利用再分析资料重建的 1962—1980 年沙尘暴年频次可以很好地再现台站观测结果, 两者的量级、峰值和谷值年份均能较好地对应, 相关系数达到 0.67,证明了 BP 神经网络模型对重建沙尘暴频次日尺度数据的有效性; 通过 CMIP5 多模式平均而重建的沙尘暴频次与观测数据及再分析重建数据有较大的差异, 且平均值的量级高于观测值。这是因为耦合模式历史试验中缺少观测资料中的年际信号强迫, 所以不能与观测数据进行简单的年际变率数值对比。尽管如此, CMIP5 模式计算结果仍能再现 1954—2007 年沙尘暴频次下降的趋势。

图4 1954—2007 年沙尘暴年频次时间序列Fig. 4 Time series of annual mean dust storm frequency in 1954-2007

图 4(b)显示, 两个模式的结果中沙尘暴年频次都存在下降趋势, 但 GFDL-CM3 的沙尘暴年频次量级高于 IPSL-CM5A-LR 和观测值。对于沙尘暴这一天气尺度的系统, 气候模式对 daily 尺度信息的模拟能力欠佳, 计算得出的沙尘暴频次与实际情况差异较大, 但能利用影响沙尘暴关键因子的长期变化趋势, 计算得到沙尘暴频次的趋势性信息。另外, 尽管模式之间存在差异, 但模式集合平均的结果能减少单一模式的高估或低估, 使结果更加可信。

3 对未来沙尘暴频次的预估

对 1954—2007 年沙尘暴频次日数据的重建结果表明, 将 BP 神经网络应用于 CMIP5 数据, 可以反映沙尘暴频次的长期变化趋势。下面将 BP 神经网络应用于 3 个 CMIP5 模式(GFDL-CM3, IPSLCM5A-LR 和 HadGESM-ES)的未来情景试验(RCP2.6和 RCP8.5)结果中, 分析到 21 世纪末, 内蒙古地区沙尘暴频次对全球变暖的响应, 作为对其预估结果的一阶近似值。

3.1 GFDL-CM3 模式计算结果

GFDL-CM3 模式计算结果显示, 在 RCP2.6 和RCP8.5 两种未来情景中, 沙尘暴年频次的高值出现次数减少, 即未来沙尘暴事件的总频次减少, 影响范围较大的沙尘暴事件频次也减少。从 1860—2100 年沙尘暴年频次的时间序列(图 5(a))可知, 据GFDL-CM3 模式数据计算得出的沙尘暴年频次过高(高于 900 的较多)。对比 RCP2.6 情景、RCP8.5情景和历史试验可知, 未来沙尘暴年频次超过 900的年份减少, 且 2004—2026 年、2060—2075 年和2080—2100 年 RCP8.5 情景的沙尘暴年频次低于RCP2.6 情景。表 2 显示, 1860—2005 年沙尘暴年频次多年平均值为 433, 2006—2100 年 RCP2.6 情景为395, RCP8.5 情景为 336, 表明未来沙尘暴频次减少。比较历史试验(图 5(b))、RCP2.6 情景(图 5(c))和 RCP8.5 情景(图 5(d))的沙尘暴年频次概率分布可知, 1860—2005 年沙尘暴年频次高于 600 的年份多于 2006—2100 年, RCP8.5 情景沙尘暴年频次高于700 的年份比 RCP2.6 情景少。图 5(e)～(g)中, 从历史试验到 RCP2.6 情景, 再到 RCP8.5 情景, 红点(沙尘暴日频次多于 15 次)逐渐减少, 表明未来出现影响范围较大的沙尘暴频次逐渐减少。

图5 GFDL-CM3 模式中1860—2100 年沙尘暴频次Fig. 5 Frequency of dust storms from GFDL-CM3 model results in 1860-2100

3.2 IPSL-CM5A-LR 模式计算结果

IPSL-CM5A-LR 模式计算结果显示, RCP8.5 情景中, 2006—2100 年出现沙尘暴年频次和日频次高值的年份少于 1860—2005 年, 但 RCP2.6 情景中沙尘暴年频次和日频次出现高值的年份多于 1860—2005 年。图 6(a)显示, RCP8.5 情景中 2006—2100年沙尘暴年频次高于 600 的年份比 1860—2005 年少, 且 2040—2100 年 RCP8.5 情景中沙尘暴年频次比RCP2.6 情景少。从图 6(b)～(d)可以看出, RCP2.6情景中沙尘暴年频次高于 600 的年份比 1860—2005年多。总体而言, IPSL-CM5A-LR 模式结果与 GFDLCM3 模式结果反映的趋势相同。

图6 IPSL-CM5A-LR 模式中 1860—2100 年沙尘暴频次Fig. 6 Frequency of dust storms from IPSL-CM5A-LR model results in 1860-2100

3.3 HadGESM-ES 模式计算结果

HadGESM-ES 模式计算结果显示, RCP2.6 情景和 RCP8.5 情景的沙尘暴年频次和日频次出现高值的次数少于 1860—2005 年这一历史时期。图 7(a)显示, RCP2.6 和 RCP8.5 两个情景中, 未来沙尘暴年频次数值区间为 200～600, 小于 GFDL-CM3 和 IPSLCM5A-LR 两个模式的计算结果, 且 2040～2100 年RCP8.5 情景中沙尘暴年频次少于 RCP2.6 情景, 即随着温度升高, 沙尘暴频次减少。图 7(b)和(c)显示,RCP2.6 和 RCP8.5 两个情景中, 沙尘暴年频次大于600 的年份数不足 10 个, 少于 1860—2005 这一历史时期。图 7(d)和(e)显示, RCP2.6 和 RCP8.5 两个未来情景中, 沙尘暴日频次高值(大于 15)出现的次数较少。

图7 HadGESM-ES 模式中1860—2100 年沙尘暴频次Fig. 7 Frequency of dust storms from HadGESM-ES model results in 1860-2100

3.4 模式集合平均值

3 个模式的集合平均值可以在一定程度上修正单个模式的不确定性。总体而言, 沙尘暴年频次和日频次高值出现的次数 1860—2005 年的历史时期多于 RCP2.6 情景, RCP2.6 情景又多于 RCP8.5 情景,即全球变暖会导致沙尘暴频次减少。图 8(b)显示,1860—2100 年沙尘暴年频次高于 900 的年份较少,且 1860—2005 年沙尘暴年频次高于 600 的年份多于2006—2100 年; 2055—2100 年, RCP8.5 情景的沙尘暴年频次均小于 RCP2.6 情景。对比图 8(b)、(c)和(d)可知, 沙尘暴年频次高于 600 的年份数 1860—2005 年的历史时期多于 RCP2.6 情景, RCP2.6 情景又多于 RCP8.5 情景, 且 RCP2.6 情景和 RCP8.5 情景中几乎没有沙尘暴年频次高于 700 的年份。未来,影响范围较大的沙尘暴出现频次减少, 且增温越多,越不利于该类沙尘暴事件发生。表 2 列出沙尘暴日频次模式集合多年平均值, 1860—2005 年为 4.7,2006—2100 年 RCP2.6 情景和 RCP8.5 情景中分别为4.0, 为 3.5, 表明未来沙尘暴会减少。

表2 历史(1860—2005 年)和未来(2006—2100 年)各模式多年平均沙尘暴年频次和日频次Table 2 Multi-year average yearly frequency and daily frequency of dust storms for three models in history (1860-2005) and future (2006-2100)

图8 CMIP5 多模式平均的700 hPa 温度和沙尘暴频次Fig. 8 Multi-model mean 700 hPa temperature and dust storm frequency

从沙尘暴出现的月份看, 台站观测资料显示,1954—2007 年这一历史时期内蒙古地区沙尘暴主要发生于 4 月[30-31]; 1954 年之前的历史时期和 2007年之后的未来, 沙尘暴也多发于 4 月, 这一季节性特点未发生变化。从图 9(a)可知, 不论是历史时期还是未来, 4 月的沙尘暴频次高于其他月份, 仍表明内蒙古地区沙尘暴多发于 4 月; 年平均沙尘暴日频次历史时期高于 RCP2.6 情景, RCP2.6 情景又高于RCP8.5 情景, 未来沙尘暴频次减少, 与前面沙尘暴年频次结果一致。从旬这一时间尺度看, 图 9(b)～(d)分别显示历史时期的 1860—2005 年、RCP2.6 情景的 2006—2100 年以及 RCP8.5 情景的 2006—2100年的沙尘暴旬频次, 其中灰度条为每年旬频次最高值, 多在 4 月出现。同时, 图 9(b)～(d)中沙尘暴旬频次高值出现次数 5 月多于 3 月。对比图 9(b)、(c)和(d)可知, 1860—2005 年沙尘暴旬频次高于 100 的年份多于 2006—2100 年。

图9 沙尘暴频次的季节分布Fig. 9 Seasonal distribution of dust storm frequency

4 结论与讨论

本文使用 BP 神经网络方法, 结合 NCEP/NCAR再分析资料、气象台站沙尘暴观测数据, 构建内蒙古地区沙尘暴频次统计模型, 并利用 GFDL-CM3,IPSL-CM5A-LR 和 HadGESM-ES 这 3 个 CMIP5 模式的历史试验、RCP2.6 情景和 RCP8.5 情景数据, 获得 1860—2100 年内蒙古地区沙尘暴日频次的长期序列, 初步探讨全球变暖背景下内蒙古地区沙尘暴频次的长期变化趋势, 得到以下结论。

1) 用 9 个关键因子的日数据作为输入、500 个节点作为隐藏层的 BP 神经网络模型是一种重建沙尘暴频次的有效方法。将再分析资料输入该模型,计算得出的 1954—2007 年沙尘暴年频次, 其量级、峰值和谷值年份以及变化趋势均与观测值相近。将CMIP5 模式历史试验数据输入该模型, 也能再现过去内蒙古地区沙尘频次下降的趋势。将该模型应用于 CMIP5 数据, 多模式集合平均值可以在一定程度上修正单一模式高估或低估的问题。

2) 在 RCP2.6 和 RCP8.5 两个未来情景中, 内蒙古地区沙尘暴频次高值的出现次数明显减少, 即未来沙尘暴事件的总频次减少, 影响范围较大的沙尘暴事件频次也减少, 且增温越多, 越不利于该类沙尘暴事件发生。3 个模式集合平均值显示, 沙尘暴年频次高于 600 的年份数历史时期多于 RCP2.6 情景, RCP2.6 情景又多于 RCP8.5 情景。2055—2100年, RCP8.5 情景中沙尘暴年频次均低于 RCP2.6 情景, 且 1860—2005 年平均沙尘暴日频次高于 2006—2100 年。就单个模式来看, GFDL-CM3 模式计算得出的沙尘暴频次值过高。

3) 沙尘暴的季节性规律不随全球变暖而改变。从沙尘暴分布的月份来看, 历史试验数据中 1860—2005 年以及 RCP2.6 情景和 RCP8.5 情景中 2006—2100 年 4 月的沙尘暴日频次和旬频次值高于 3 月和5 月, 历史时期和未来沙尘暴都多发于 4 月。历史时期和未来, 沙尘暴年平均日频次 5 月多于 3 月,沙尘暴旬频次高值的出现次数也是 5 月多于 3 月。

前人有关未来沙尘暴频次变化预估的研究较少。Pu 等[32]用多元线性回归的统计方法, 得到RCP8.5 情境下美国南部沙尘暴增加、北部沙尘暴减少的结果。本研究在前人工作的基础上, 基于CMIP5 数据重建内蒙古地区沙尘暴频次 daily 指数的人工神经网络统计模型, 据此对未来全球变暖背景下内蒙古地区沙尘暴的长期变化趋势进行预估,初步得出未来(2006—2100 年)中国沙尘暴频次低于历史时期(1860—2005 年)的结论, 与 Zhu 等[14]基于1954—2007 年观测数据的研究结果(贝加尔湖温度升高, 大气斜压性减弱, 气旋减少, 过去 50 年中国沙尘暴减少, 推测未来中国沙尘暴频次受同样机制影响而保持下降趋势)一致。

本文工作的不确定性和局限性主要来自以下 3个方面。

1) 近年来人工神经网络有了长足的进步和跨越式发展, 本文只利用了较为基础的 BP 神经网络,未来可以尝试卷积神经网络等更多的算法, 改进统计模型。

2) 限于目前可用数据, 本研究只分析了 3 个CMIP5 模式中内蒙古地区历史时期以及两个未来情景的沙尘暴频次变化。今后, 可将本文模型应用于不同的模式、时段和情景, 更广泛地探讨气候变化与中国沙尘暴频次变化之间的关系。

3) 近地面的植树造林人为活动主要影响起沙过程。裸露地表被植被覆盖后, 沙尘不易被扬起,沙尘源减少。但是, 沙尘远距离输送依赖于 3～5 km的高空风场, 天气形势起主导作用。本文在假定沙尘源无显著变化的情况下研究内蒙古地区沙尘暴频次的长期变化, 仅分析天气过程对沙尘输送的影响,未考虑人为植树造林等活动的影响。后续工作中可以加入植被变化这一影响因素, 进一步从沙尘源以及沙尘输送两个方面, 更全面地评估未来沙尘暴频次的变化趋势。