一种地表水质对旱涝事件的响应分析方法与流程

本发明涉及水污染预防与应对领域，尤其涉及一种地表水质对旱涝事件的响应分析方法。

背景技术：

根据环境保护部统计，2006年至2012年间941起突发环境事件中因自然灾害引发的比例为9％，暴雨、洪涝和干旱等极端气候事件已经成为水污染事故的重要原因，并且多数污染事故是由于工业点源污染和农业面源污染因暴雨冲刷将污染物带入水体而引发的，旱涝等极端气候事件已经严重威胁到地表水环境安全。极端天气如暴雨、洪涝、干旱等发生的强度增大和频率升高，将会进一步增加污染事故发生的可能性和严重性，因此，研究地表水质与旱涝事件之间的关系具有重要的意义。

目前旱涝事件等极端气候事件呈现广发、频发的态势，这些极端气候事件的发生会对地表水质造成突发性的影响。干旱发生后，不仅会导致水量短缺，同时会伴随水质恶化的危害。干旱对地表水环境的影响不仅通过影响其水文要素过程，还会因为干旱导致的缺水对地面植被土壤、水生动植物及微生物等各个方面产生带来不利的影响。1992年，Margarida等研究发现发生60年一遇的干旱时，葡萄牙地区50个主要水库由于河道径流量减少导致水库水质下降。1997年，Mulholland发现干旱发生后，因为河道径流量的减少导致水体中营养盐、有机物和污染物等浓度上升而污染水质。2009年，Delpla等分析土壤中C、N、P等元素的矿化及分解过程，提出干旱时期这些元素会加速分解，在干旱期结束后随着降水及土壤侵蚀进入水体中。

在地质复杂、降水量大的山区，强降水甚至会导致泥石流和暴洪灾害的发生，不仅会导致大量的人畜伤亡和财产损失，还会严重的污染地表水环境。2007年，Polyakov等提出突发性强降水不仅造成土壤侵蚀，还会将土壤中大量的有机物、营养盐以及病菌等物质带入水体，造成水质污染严重。

技术实现要素：

本发明的目的就在于针对现有旱涝事件对地表水质的影响分析方法的不足，提出一种地表水质对旱涝事件的响应分析方法。该方法基于气候模式预估结果、SWAT模型和统计分析方法，结合标准化降水指数、24小时降水量和连续无雨日等指标，筛选出干旱事件、洪涝事件和旱涝急转事件的典型代表月份，对比分析旱涝事件对典型月水质的影响情况，进而分析地表水质变化与气候变化的响应关系。

为了实现上述目的，本发明采用的技术方案是：一种地表水质对旱涝事件的响应分析方法，方法步骤如下，

步骤一旱涝事件历史演变及发展趋势分析

1)旱涝事件历史演变规律

根据研究区历史月降水量数据计算标准化降水指数，从月尺度与年尺度两个方面分析当地旱涝事件的历史演变规律。

2)气候模式适用性评价与修正

选取国际耦合模式比较计划第五阶段(CMIP5)提供的五种气候模式：GFDL-ESM2M、HadGEM2-ES、IPSL-CM5A-LR、MIROC-ESM-CHEM和NorESM1-M，分别模拟RCP2.6、RCP4.5、RCP8.5三种代表性排放情景下的降水量，经过校准和降尺度后，与实测值进行对比，根据计算出的相对误差、均方根误差和相关系数特征值结果，以及每个模式模拟变化趋势，综合分析和评估不同气候模式的降水模拟能力。基于上述分析结果，选取模拟结果相对较好的气候模式做集合平均，进行气候模式修正。

3)旱涝事件发展趋势预估

利用修正后的气候模式模拟出研究区未来月降水量，计算标准化降水指数，从月尺度与年尺度两个方面分析当地未来旱涝事件发展趋势。

步骤二SWAT模型构建与校准

1)模型输入数据准备及格式化处理

SWAT模型为了使模拟结果更加准确，需要输入的数据有土地利用类型、土壤类型、温度降雨等气象条件、河流水质等，格式需要经过标准化处理后才可以输入到模型中去，其中数字高程模型图、土壤类型图、土地利用类型图需要使用ArcGIS采用GRID或shp格式生成地形图，而降雨、气温等气象数据和水质数据等一般以dbf格式生成输入文件，具体格式见下表。

表1 SWAT模型输入数据及格式

2)模型校准和验证

①确定模型评价指标

选取相对误差(Re)、决定系数(R²)和Nash-Suttcliffe效率系数(Ens)三个指标用于评价模型的适用性。根据模拟标准，模拟结果与实测值误差应小于20％，决定系数R²大于0.6且Ens大于0.5，当模拟结果满足以上条件时，模拟结果才可以进行使用。

②参数敏感性分析

敏感性分析是评估各项参数对模拟结果的影响程度，优先选择影响程度大的参数进行调整。SWAT模型敏感性分析是模型调参和校验的前提和基础工作。根据模拟结果进行敏感性分析后，从SWAT模型众多参数中筛选敏感性较高参数，优先进行参数率定。

③模型校准

基于上述步骤，对模型进行校准和验证。

步骤三总氮、总磷计算

根据SWAT模型对非点源污染的迁移转化过程和原理，可以对非点源污染物进行模拟，模拟的对象主要是氮元素和磷元素负荷量，由于SWAT模型在对河道总氮和总磷进行模拟时，将氮和磷分成不同形态进行循环转化，模拟得到的是各种形态的氮和磷的负荷，总氮和总磷负荷量采用下公式进行计算：

TN＝ORGN_OUT+NO₃_OUT+NH₄_OUT+NO₂_OUT (1)

式中：TN为总氮负荷；ORGN_OUT为有机氮负荷；NO₃_OUT为硝态氮负荷；NH₄_OUT为氨态氮负荷；NO₂_OUT为亚硝态氮负荷。

TP＝ORGP_OUT+MINP_OUT (2)

式中：TP为总磷负荷；ORGP_OUT为有机磷负荷；MINP_OUT为矿物磷负荷。

(4)地表水质对旱涝事件的响应

1)旱涝事件对氮磷的影响分析

基于修正后的气候模式预估结果和校准过的SWAT模型，结合标准化降水指数(SPI)，分析旱涝事件对氮磷的影响。根据模拟结果，首先分别挑选出总氮、总磷负荷量最多的20个月份，根据SPI值进行旱涝频率分析。继续筛选出总氮、总磷含量大于多年平均负荷的月份，对比其中发生旱涝事件的比例。

2)典型干旱事件对氮磷的影响

通过对比研究标准化降水指数、平均无雨日和降水量指标，筛选出10个典型干旱月份，模拟典型干旱月总氮、总磷污染负荷变化情况。

3)典型洪涝事件对氮磷的影响

根据标准化降水指数和24小时降水量等指标，筛选出10个典型洪涝月，模拟暴雨情景下典型月份总氮、总磷污染负荷变化情况。

4)典型旱涝急转事件对氮磷的影响

根据时间的连续性、标准化降水指数、24小时降水量、连续无雨日等指标，筛选出四个典型旱涝急转月，模拟旱涝急转情景下典型月份总氮、总磷污染负荷变化情况。

作为优选，步骤一中，标准化降水指数SPI采用Γ分布概率来描述径流量的变化，将偏态概率分布的径流量进行正态标准化处理，最终用标准化降水累积频率分布来划分干旱等级。假设某一时间尺度下的径流量为x，则满足伽马分布的概率密度函数f(x)为：

$f\left(x\right)=\frac{1}{{\mathrm{\β}}^{\mathrm{\α}}\mathrm{\Γ}\left(\mathrm{\α}\right)}{x}^{\mathrm{\α}-1}{e}^{-x/\mathrm{\β}},(x\>0,\mathrm{\α}\>0,\mathrm{\β}\>0)---\left(3\right)$

对伽马分布概率进行整体标准化处理求解得到SPI：

$SPI=S\frac{t-(c_{2}t+c_1)t+c_0}{\[(d_{3}t+d_2)t+d_1\]+1}---\left(4\right)$

式中：c₀＝0.802853，c₁＝0.0100328，c₂＝0.010328，d₁＝1.432788，d₂＝0.189269，d₃＝0.001308。当F＞0.5时，S＝1；当F≤0.5时，S＝-1。

表2 SPI旱涝等级划分

作为优选，步骤一中，进行气候模式适用性评价时，代表性温室气体浓度路径(RCPs)是IPCC在2011年发布的第五次评估报告中提出新的排放情景，将气候、大气、碳循环预估与排放和社会经济情景相结合，为气候模式提供了更加可靠的背景依据。根据辐射由高到低分为4种情景：RCP8.5、RCP6.0、RCP4.5、RCP2.6。具体概况如下表：

表3 RCPs情景概况

作为优选，步骤一中，进行气候模式适用性评价时，由于全球气候模式对区域的气候要素变化描述分辨率不足，为此在区域气候变化研究中需要进行降尺度处理，采用统计降尺度方法，结合双线性插值，将输出数据统一插值到0.5°×0.5°的网格上，并基于概率分布进行的统计偏差进行修正。

作为优选，步骤一中，进行气候模式适用性评价时，将全球气候模式数据通过ArcGIS转化为空间栅格，从大尺度气候信息转化提取为区域尺度，将提取出来的气候模式结果根据IDW插值计算月平均降水。经过校准和降尺度的GFDL-ESM2M、HadGEM2-ES、IPSL-CM5A-LR、MIROC-ESM-CHEM、NorESM1-M五个模式分别在不同的RCPs情景下统计模拟降水量、多模式集合平均降水量结果，将其与实测值进行对比分析。

作为优选，步骤二中，模型校准使用SWAT模型官方网站提供的专业SWAT模型参数校准和调节工具SWAT-CUP(Calibration and Uncertainty Programs)，首先设置SWAT-CUP运行的相关参数，包括SWAT模型初步模拟、主要校准的参数、参数变化的范围、满足条件的评价标准、最大模拟次数等等，然后将研究区域内实测径流数据输入SWAT-CUP中作为验证，SWAT-CUP的原理是调节参数在允许变化的范围内任意取值进行参数的优化，每次优化后的参数自动调用SWAT模型后进行模拟。将每一次模拟值与实测值进行对比，计算各项评价指标，即完成一次调参，在设定的模拟次数结束后评价各项指标是否满足研究需要，当指标结果不能满足时，在新得到各项参数继续使用SWAT-CUP循环的验证，直至参数达到合理值。

与现有技术相比，本发明的优点在于：本发明方法基于气候模式预估结果、SWAT模型和统计分析方法，结合标准化降水指数、24小时降水量和连续无雨日等指标，筛选出干旱事件、洪涝事件和旱涝急转事件的典型代表月份，对比分析旱涝事件对典型月水质的影响情况，可有效地分析地表水质变化与旱涝事件之间的响应关系，可广泛应用于研究地表水质对旱涝事件的响应，特别适用于河套地区地表水质对旱涝事件的响应分析。

附图说明

图1为地表水质对旱涝事件的响应分析方法流程图；

图2为河套地区典型干旱情景下氮磷负荷模拟结果；

图3为河套地区典型洪涝情景下氮磷负荷模拟结果；

图4为河套地区典型旱涝急转情景下氮磷负荷模拟结果。

具体实施方式

下面将对本发明作进一步说明。

实施例：参见图1，一种地表水质对旱涝事件的响应分析方法，选取灌溉系统发达的河套地区进行实证研究，方法步骤如下：

(1)旱涝事件历史演变及发展趋势分析

1)旱涝事件历史演变规律

根据河套地区1971-2010年月降水量数据计算标准化降水指数，从月尺度与年尺度2个方面分析当地旱涝事件的历史演变规律。

从时间尺度上看，SPI值对时间的敏感性不同，月尺度下旱涝变化比年尺度更为显著，但总体规律表现为河套地区的旱涝事件呈交替出现。按照对季节划分，以3-5月为春季、以6-8月为夏季、以9-11月为秋季、以12-翌年2月为冬季，根据不同季节统计河套地区的旱涝情况。1971-2010年各季节旱涝频次分布如表4，在过去的40年里，春季为旱涝事件的高发时段，发生频率高达38.3％，其他季节发生旱涝事件的频率基本相同，极端干旱发生率较高，全年各季节均有可能出现，极端洪涝事件发生相对较少。从年尺度上看，旱涝事件发生年份高达40％，体现出河套地区历史旱涝事件的频发。

表4 1971-2010年旱涝频次分布

2)气候模式适用性评价与修正

选取国际耦合模式比较计划第五阶段(CMIP5)提供的五种气候模式：GFDL-ESM2M、HadGEM2-ES、IPSL-CM5A-LR、MIROC-ESM-CHEM和NorESM1-M，分别模拟RCP2.6、RCP4.5、RCP8.5三种代表性排放情景下的降水量，经过校准和降尺度后，与实测值进行对比，根据计算出的相对误差、均方根误差和相关系数特征值结果，以及每个模式模拟变化趋势，综合分析和评估不同气候模式的降水模拟能力。基于上述分析结果，选取模拟结果相对较好的气候模式做集合平均，进行气候模式修正。

结果显示，5个气候模式对降水均具有一定的模拟能力，但是整体上相关性不强，相关系数主要集中在0.6-0.8之间，绝大多数模拟结果比实测值高，并且差异比较明显，因此，需对模式进行修正。

选择RCP8.5情景下HadGEM2-ES和NorESM1-M这两个模式做集合平均，将修正后的模拟结果继续与实测值进行对比分析，可以发现修正后的模拟结果降水量的波动程度变低，变化趋势基本保持一致，相关系数达到0.806，属于高度线性相关，模拟结果与实测值十分接近，相对误差只有3.61％，基本符合实际降水情况，修正后的气候模式可以用于河套地区未来旱涝事件预估。

3)旱涝事件发展趋势预估

利用修正后的气候模式模拟出河套地区2021-2050年月降水量，计算标准化降水指数，从月尺度与年尺度2个方面分析当地未来旱涝事件发展趋势。

河套地区的旱涝事件依然呈交替出现，从年尺度上来看，发生旱涝事件的年份占总年份的近1/3，旱涝事件发生的频率较高，并且发生时间主要集中在2025-2035时间段。从月尺度上来看，发生洪涝事件的月份增多，干旱事件发生的月份更加集中在春冬季，这可能与未来降水时空分布改变有关。

2021-2050年河套地区各季节旱涝频次分布见表5，可以看出春冬两季发生旱涝事件最多，其发生频率高达39％，其他季节发生旱涝事件频率也接近34％。预计未来发生严重干旱以上的频率近7％，发生严重洪涝以上的频率为5.2％。未来旱涝事件发生的频率很高，会进一步引发水资源短缺、水质恶化等一系列生态环境效应，加剧水资源的供需矛盾。

表5 2021-2050年旱涝频次分布

(2)SWAT模型构建与校准

1)模型输入数据准备及格式化处理

SWAT模型为了使模拟结果更加准确，需要输入的数据有土地利用类型、土壤类型、温度降雨等气象条件、河流水质等，格式需要经过标准化处理后才可以输入到模型中去，其中数字高程模型图、土壤类型图、土地利用类型图需要使用ArcGIS采用GRID或shp格式生成地形图，而降雨、气温等气象数据和水质数据等一般以dbf格式生成输入文件，具体格式见下表。

表1 SWAT模型输入数据及格式

地形数据来自STRM30米分辨率的DEM数据高程数据，根据DEM进行子流域划分、河网水系的生成；土地利用数据根据模型需要划分为6种类型：耕地、林地、草地、水域、城乡工矿居民用地和未利用土地；土壤数据包括土壤类型和土壤理化性质，土壤类型根据查询当地土壤数据资料得到，土壤理化性质根据计算和查阅资料得到包括土壤分层、土壤水文分组、饱和导水率、土壤粒径组成等一系列数据；气象数据根据河套地区及周边13个气象站点逐日实测数据，包括降水、温度、辐射、湿度和风速。

2)模型校准和验证

①确定模型评价指标

选取相对误差(Re)、决定系数(R²)和Nash-Suttcliffe效率系数(Ens)三个指标用于评价模型的适用性。根据模拟标准，模拟结果与实测值误差应小于20％，决定系数R²大于0.6且Ens大于0.5，当模拟结果满足以上条件时，模拟结果才可以进行使用。

②参数敏感性分析

敏感性分析是评估各项参数对模拟结果的影响程度，优先选择影响程度大的参数进行调整。SWAT模型敏感性分析是模型调参和校验的前提和基础工作。根据模拟结果进行敏感性分析后，从SWAT模型众多参数中筛选敏感性较高参数，优先进行参数率定。

根据敏感性分析后得到了河套地区SWAT模型模拟过程中最敏感的10个参数，在接下来的调参过程中以优点调节这个10个参数为主，利用SWAT-CUP进行调参过程，选择以这10个参数为主要调节参数以及每个参数变化的范围进行校准，调节参数在允许变化的范围内任意取值进行参数的优化，每次优化后的参数自动调用SWAT模型后进行模拟，直至参数达到合理值。参数调节过程结束后，计算各项评价指标，分析SWAT模型是否适用于河套地区。

表6敏感性分析参数排序

③模型校准

基于上述步骤，对模型进行校准和验证。结果显示，河套地区模拟月径流量和实测月径流量的变化趋势基本一致，模拟结果达到了模型的评价标准，SWAT模型可以应用于河套地区径流、氨氮的长时间序列模拟。河套地区月径流量率定期和验证期的模拟结果和实测值各评价指标值见表7。

表7河套地区月径流量模拟评价指标

(3)总氮、总磷计算

根据SWAT模型对非点源污染的迁移转化过程和原理，可以对非点源污染物进行模拟，模拟的对象主要是氮元素和磷元素负荷量，由于SWAT模型在对河道总氮和总磷进行模拟时，将氮和磷分成不同形态进行循环转化，模拟得到的是各种形态的氮和磷的负荷，总氮和总磷负荷量采用下公式进行计算：

TN＝ORGN_OUT+NO₃_OUT+NH₄_OUT+NO₂_OUT (1)

式中：TN为总氮负荷；ORGN_OUT为有机氮负荷；NO₃_OUT为硝态氮负荷；NH₄_OUT为氨态氮负荷；NO₂_OUT为亚硝态氮负荷。

TP＝ORGP_OUT+MINP_OUT (2)

式中：TP为总磷负荷；ORGP_OUT为有机磷负荷；MINP_OUT为矿物磷负荷。

(4)地表水质对旱涝事件的响应

1)旱涝事件对氮磷的影响分析

基于修正后的气候模式预估结果和校准过的SWAT模型，结合标准化降水指数(SPI)，分析旱涝事件对氮磷的影响。

根据模拟结果，首先分别挑选出总氮、总磷负荷量最多的20个月份，参见表8，根据SPI值进行旱涝频率分析，发现其中旱涝事件发生的月份分别有10个和12个，氮磷负荷污染最严重的月份中有超过一半发生了旱涝事件，说明发生旱涝事件后有可能会造成严重的水体污染。

继续筛选出总氮、总磷含量大于多年平均负荷的月份，对比其中发生旱涝事件的比例，可以得出氮负荷大于均值的月份有137个，发生旱涝事件的月份有50个；磷负荷大于均值的月份有54个，发生旱涝事件的月份有25个；氮、磷负荷同时大于均值的月份46个，发生旱涝事件的月份有20个。当氮、磷负荷高于均值时旱涝事件发生比例分别为36.5％、46.3％和43.5％。

由此可见，地表水质对干旱事件和洪涝事件均呈现出较为显著的响应，旱涝事件发生后与水体中氮磷负荷量增多有着明显的联系。

表8总氮、总磷负荷量排序

2)典型干旱事件对氮磷的影响

通过对比研究标准化降水指数、平均无雨日和降水量指标，筛选出10个典型干旱月份，参见表9，模拟典型干旱月总氮、总磷污染负荷变化情况。

通过10个典型月的模拟结果来看，参见图2，根据干旱发生时间不同分为2种情形，一种是干旱发生在上半年，另一种是干旱发生在下半年，可以看出干旱发生的时间不同，对水质造成的影响是不一样的。干旱发生在上半年时，干旱发生的典型月污染物的负荷量很小，虽然在此阶段河流水体中污染物较少，但是随着农业耕作的开始，氮磷负荷已经在土壤中大量积累，随着干旱过后夏季的到来，降水量增多，径流增加，使得接下来的月份水体中氮磷负荷量急剧增加。而干旱发生在下半年时，虽然氮磷等有机物同样在土壤中积累，但是由于即将进入冬季缺水时期，降雨量逐渐减少，温度降低，河流开始结冰进入冰封期，土壤中的污染物不再随着径流进入水体中，导致河流中的氮磷元素反而开始减少。

表9典型干旱月

3)典型洪涝事件对氮磷的影响

根据标准化降水指数和24小时降水量等指标，筛选出10个典型洪涝月，参见表10，模拟暴雨情景下典型月份总氮、总磷污染负荷变化情况。

通过10个典型月的模拟结果显示，参见图3，同样根据发生暴雨的时间分为2种情形。一种是暴雨发生在夏季，夏季的典型月氮磷负荷量不仅是一年中含量最高的月份，而且相比上一个月的污染负荷量有一个飞速的提升，速度使得水体中污染物负荷量达到峰值；另一种情形是暴雨发生在其他季节，该类型的典型月氮磷负荷量没有达到一年中最大值，但是会导致该典型月氮磷负荷量相比正常情况下相同月份的污染负荷量高出很多，由此可以看出暴雨对河流中氮磷含量的累积有着重大影响。

表10典型洪涝月

4)典型旱涝急转事件对氮磷的影响

根据时间的连续性、标准化降水指数、24小时降水量、连续无雨日等指标，筛选出4个典型旱涝急转月，参见表11，模拟旱涝急转情景下典型月份总氮、总磷污染负荷变化情况。

通过4个典型月的模拟，参见图4，其中红点表示发生干旱，黑点表示发生洪涝，根据干旱和洪涝的发生时间分为2种情形：先旱后涝和先涝后旱。发生先旱后涝的模拟结果与典型干旱的模拟结果类似，根据干旱发生时间不同导致的结果不同，在春旱结束后，随着农业耕作的开始和整个冬季污染物的积累，再遭遇暴雨洪涝后会迅速加重水体污染；在夏末秋初时发生干旱，由于夏季的降水较大已经将土壤中的污染物冲刷进入水体，干旱后即使发生了洪涝不会使河流湖泊中的污染物有显著增加，但是由图可以发现，干旱结束由于发生了洪涝事件，总氮负荷量的下降速度减缓，而总磷负荷量有少量增加。先涝后旱的模拟结果表明，洪涝过后河流水体中氮磷负荷量已经达到峰值，随后突发干旱灾害虽然不会增加污染负荷量，但是由于干旱会导致河流水量减少，而氮磷负荷量的降解速度相对较慢，使得水体中污染物浓度有一个较大的提升。

表11典型旱涝急转月

本发明该方法基于气候模式预估结果、SWAT模型和统计分析方法，结合标准化降水指数、24小时降水量和连续无雨日等指标，筛选出干旱事件、洪涝事件和旱涝急转事件的典型代表月份，对比分析旱涝事件对典型月水质的影响情况。上述结论说明利用本发明提出的方法能够较准确地分析地表水质变化与旱涝事件之间的响应关系，可为下一步的管理措施制定提供辅助支撑，可应用于实际的项目分析中。

以上对本发明所提供的一种地表水质对旱涝事件的响应分析方法进行了详尽介绍，本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，对本发明的变更和改进将是可能的，而不会超出附加权利要求所规定的构思和范围，综上所述，本说明书内容不应理解为对本发明的限制。