SWMM模型参数全局敏感性分析

段明印，李传奇，韩典乘，杨幸子

(山东大学土建与水利学院，山东 济南 250061)

0 引 言

城市化的进程使得我国城市雨洪灾害问题日益严峻，因此我国在防治雨洪灾害的同时正在大力推行海绵城市的建设，暴雨管理模型(SWMM)可以模拟城市降雨地表径流的产生和在排水管网中的输送，能够用来指导海绵城市的建设，因此在国内得到了较为广泛的运用[1]。

准确率定SWMM的参数是精确构建SWMM模型、指导海绵城市建设的前提。但是SWMM参数较多，在率定时为提高效率可仅率定对输出变量影响较为显著的参数，那些对输出变量影响较小的参数可直接取经验值。因此，需要对SWMM参数进行敏感性分析，识别出影响较大的参数[2]。参数敏感性分析方法分为局部敏感性分析方法和全局敏感性分析方法，局部分析方法只是分析单个参数对输出变量的影响，因其计算速度快，得到了较广泛的运用，但是忽略了参数之间相互作用对输出变量的影响，在“异参等效”的情况下会使一些敏感参数无法识别；全局分析方法可以分析参数之间相互作用对输出变量的影响，能够在“异参等效”的情况下识别出敏感参数，适用于参数众多的模型[3]。

常见的全局敏感性分析方法有多元逐步回归法[4]、Sobol法[5]、FAST法[6]、偏相关法、偏秩相关法等[3]。这些方法都是通过分析参数和输出变量的线性相关关系来识别敏感参数，忽略了参数和输出变量的非线性相关关系，可能导致某些敏感参数无法识别。因此，本文以某小区为例，选取能够分析变量之间非线性关系的互信息法分析参数全局敏感性，将以互信息法分析的结果与能够分析变量之间线性关系的偏秩相关法分析的结果相结合，精确识别所有敏感参数，为下一步高效率定打下基础。

1 材料与方法

1.1 SWMM模型参数

在利用SWMM模型对兴隆山校区进行降雨径流模拟时，需要输入的参数共计14个。这14个参数中的管渠长度和不透水面积可由设计资料读取或者直接测量，所以仅对余下的12个参数的敏感性进行分析。其中子汇水区的面积、宽度、坡度这3个输入参数具有显著空间特性，输入时存在测量误差，所以定义3个修正系数：K-Area(面积修正系数)、K-With(宽度修正系数)和K-Slope(坡度修正系数)，确定这3个修正系数的分布范围，则SWMM模型的输入值为修正系数和对应参数的乘积[7]。依据SWMM使用手册和国内相关文献[7-11]确定参数的范围，依据研究区域内的实测值确定参数初始值，参数分布及取值如表1所列。

表1 SWMM参数取值表Tab.1 The parameters table of the storm water management model

1.2 SWMM模型构建

以山东大学兴隆山校区为研究区域，该区总面积为63.82 万m2，其中绿地面积较大，绿地面积高达30.82 万m2，建筑占地面积为9.2 万m2，道路占地面积8.9 万m2，块石铺装占地面积1.9 万m2，硬质铺装面积8.7 万m2，水体面积2.8 万m2，足球场面积1.5 万m2，研究区域概况如图1所示。

图1 研究区域概况Fig.1 The overview of study area

结合兴隆山校区的地势条件、地下管网资料、下垫面类型以及校区内的河道等资料，构建研究区的SWMM模型。SWMM模型包括45个子汇水区域，55个管道铰点，55条雨水管道，3个排水口，兴隆山SWMM模型如图2所示。

图2 研究区域SWMM模型Fig.2 Storm water management model of the study area

1.3 抽样方法

常见的抽样方法有简单随机抽样、分层抽样和整群抽样，分层抽样方法的原理是将样本空间分成互不交叉的若干层，然后从各层中抽取样本，相对于简单随机抽样和整群抽样，分层抽样能够使得样本点均匀地分布在样本空间内，避免样本集中的问题。因此本文选用拉丁超立方抽样法在SWMM模型参数的分布空间内进行抽样，拉丁超立方抽样法属于多维分层抽样方法，能够高效的抽取分布较为均匀的样本[8]。

首先通过拉丁超立方抽样法对有关降雨径流的12个SWMM模型参数随机抽样1 000次，得到1 000组12个SWMM模型参数，然后通过MATLAB编程将这1 000组参数导入到SWMM模型的输入文件中，生成1 000组新的SWMM输入文件，利用SWMM模型将这1 000组输入文件依次模拟，得到这1 000组参数对应的SWMM输出文件，提取这1 000组输出文件中某排放口的3个重要水文变量：流量峰值、流量峰值发生时间、总产流量，因为研究区域主要通过1号排放口排放总产流，因此提取1号排放口的3个水文变量作为3个因变量，以12个输入参数作为自变量，分别以偏秩相关法和互信息法分析自变量和因变量之间的关系。

1.4 偏秩相关分析方法

偏秩相关分析是一种对变量分布不做要求的非参数统计分析方法，能够对等级数据进行线性相关分析，在分析两变量的相关关系时能够控制其他变量的影响，因此能够精确识别变量间的线性相关关系，可用做分析参数的全局敏感性[9]。

令输出结果为Y，输入参数为X1，X2，…，Xn，设这n+1个变量之间的Spearman秩相关系数r组成的矩阵为T，求矩阵T的逆矩阵得矩阵C。

(2)

则得某输入参数Xi与输出结果Y之间的偏秩相关系数PXi为：

(3)

式中：r为两变量间Spearman秩相关系数；R2Xi为X1，X2,…,Xi-1,Xi+1,…,Xn线性回归Xi的可决系数；RY2为它可决系数；Bi为X1，X2，…，Xn线性回归Y的标准回归系数；PXi为偏秩相关系数，它的绝对值大小为参数敏感性大小，它的正负表示呈正相关还是负相关。

1.5 互信息分析

互信息是指两个变量之间或者多个变量之间共享的信息量，能够度量变量之间的相关性，相关性越强，互信息越大。基于互信息的R统计量可以用来分析参数的全局敏感性，R统计量能够精确识别参数和输出结果之间的非线性相关关系，R统计量的大小代表参数的敏感性大小，但R统计量只能识别参数敏感性的大小，无法识别参数和输出结果之间是正相关还是负相关[10]。

两离散型随机变量X与Y之间的互信息定义为：

(5)

R统计量表示为：

(6)

式中：p(x)为离散型随机变量X取值x的概率；p(y)为离散型随机变量Y取值y的概率；p(x,y)为X取值x、Y取值y同时发生的概率；I(x,y)为变量X、Y的互信息；b为互信息的量纲。本文互信息的量纲采用自然常数e，将SWMM模型的输入参数视为离散型随机变量X，将SWMM模型输出结果视为离散型随机变量Y，求解变量X与Y之间的互信息以及R统计量，分析参数对输出结果影响的大小。在计算时为提高计算效率，将变量X分为i等份，将变量Y分为i等份，用pi表示变量X第i份中的数量占变量X总数量的概率，同时用pj表示变量Y第j份中的数量占变量Y总数量的概率，将分成i等份的变量X与分成j等份的Y组成联立表，求得变量X和变量Y的联合概率，用pij表示变量X在第i份取值、自变量Y在第j份取值同时发生的概率，则互信息的公式变为：

(7)

式中：i、j为输入参数和输出结果分成的等份数，本文都分为10等份；

1.6 降雨资料

降雨资料采用兴隆山校区内雨量站(测站编码：41822450)测得的降雨数据，选取20120708场次的降雨，该场次降雨历时为13个小时，降雨量为47.5 mm，雨峰系数约为0.4，降雨强度过程如图3所示。

2 结果与分析

2.1 初步分析

为初步分析SWMM模型的12参数和3个水文输出结果的关系，分别画出这1000组12个参数和3个输出结果的残差散点图，其中散点图的横坐标为某一参数与控制参数回归分析得到的残差，纵坐标为输出结果与控制参数回归分析得到的残差。各个参数和峰值流量的关系如图4所示，观察发现Manning-N与峰值流量的残差散点图几乎在一条一直线上，具有高度的线性关系，呈负相关。各个参数和峰值时间的关系如图5所示，观察发现N-perv和峰值时间具有较强的线性关系，呈负 相关。各个参数和总产流量的关系如图6所示，观察发现Manning-N、Min-Rate和总径流量具有较强的负相关关系，其中Manning-N和总产流量的负相关性强于Min-Rate。

图3 研究区域20120708场次降雨强度过程Fig.3 Rainfall intensity processes of label 20120708 in the study area

图4 参数和峰值流量残差散点图Fig.4 Residual scatter diagram of parameters and peak flow

图5 参数和峰值时间残差散点图Fig.5 Residual scatter diagram of parameters and peak time

图6 参数和总产流量残差散点图Fig.6 Residual scatter diagram of parameters and total runoff volume

2.2 偏秩相关分析

通过MATLAB编程利用偏秩相关分析方法分析参数和输出结果的相关性，分析结果如表2所列。由表2可以得出参数和输出结果之间线性相关关系的强弱，即参数的敏感性大小。对峰值流量来说，最敏感的参数是Manning-N，它与峰值流量的偏秩相关系数P高达-0.996，具有高度的线性关系；第二敏感参数是N-perv，P值的绝对值位0.670，排名前两位的敏感参数和峰值流量都是负相关关系，这是因为管渠和透水地表的糙率系数增大会使滞蓄雨水量增大，从而使得峰值流量减小；其他参数的敏感性相对较小，说明Manning-N和N-perv对峰值流量起着决定性作用。对峰值时间来说，最敏感的参数是N-perv，P值为-0.417，呈负相关，其他参数的敏感性都较小。对总产流量来说，最敏感的参数是Manning-N，P值高达-0.936，具有高度的线性关系，第二敏感参数是Min-Rate，P值为-0.838，其与总产流量也具有较强的线性关系；Manning-N、Min-Rate与总产流量都是负相关关系，Manning-N和Min-Rate的增大将会使总产流量减小，这是因为最小渗透率增加将会使雨水下渗量增大，管渠糙率增大会使管渠滞蓄雨水量增大，这些都会使总产流量减小。

表2 偏秩相关分析结果Tab.2 Result of partial rank correlation analysis

2.3 互信息分析

通过MATLAB编程利用互信息相关原理分析参数和输出结果之间的非线性相关关系，结果如表3所列。根据表3中的R统计值可以判断参数和输出结果之间非线性关系的强弱，及参数敏感性的大小。对峰值流量来说，最敏感的参数是Manning-N，R值高达0.886，这说明Manning-N和峰值流量之间具有较强的非线性相关关系；其他参数的R值较小，说明其他参数和峰值流量的非线性关系较弱。对峰值时间来说，N-perv是最敏感的参数，R值为0.457，说明N-perv和峰值时间的非线性关系较弱；第二敏感参数是K-Width，R值仅为0.337，对峰值时间影响较小，其他参数对峰值时间的影响也都较小。对总产流量来说，Manning-N是最敏感的参数，R值为0.574；第二敏感参数是Min-Rate，R值仅为0.372，说明其对总产流量的影响较小。

将偏秩相关分析的结果和互信息分析的结果进行对比发现：两种方法所求的对峰值流量最敏感的参数都是Manning-N，但是用偏秩相关法算出的敏感性大小为0.997，大于用互信息法算出的敏感性大小，说明Manning-N和峰值流量的线性关系强于非线性关系；对峰值时间而言，两种方法求得的最敏感参数都是N-perv，但是互信息求得的R值大于偏秩相关法求得的P值绝对值，说明N-perv与峰值时间的非线性关系强于线性关系；对总产流量来说，两种方法求得的最敏感参数都是Manning-N，且偏秩相关法求得的P值绝对值大于互信息法求得的R值，说明Manning-N和总产流量的线性关系强于非线性关系；两种方法识别出的较敏感参数存在差异，这是因为偏秩相关法是通过分析线性相关关系来识别敏感参数，互信息法是通过分析非线性相关关系来识别敏感参数，因此二者识别的敏感参数存在差异，但是这两种方法可以共同使用，从而全面识别敏感参数，避免只用一种方法而导致的疏漏。

3 结 语

(1)偏秩相关分析方法能够通过分析参数和输出结果之间的线性关系来识别敏感参数，互信息法能够通过分析参数和输出结果之间的非线性关系来识别敏感参数，这两种方法都是全局敏感新分析方法，将两种方法结合使用能够精确识别所有的敏感参数。

(2)通过两种分析方法得知，Manning-N是对峰值流量影响最大的参数，Manning-N的增大将会使峰值流量减小；N-perv是对峰值时间影响最大的参数，N-perv的增大将会使峰值时间提前；Manning-N是对总产流量影响最大的参数，Manning-N的增大将会使总产流量减小。

(3)通过敏感性分析，精确识别所有敏感参数，能够在率定过程中提高效率。

□

[1] 官奕宏,吕 谋,王 灿,等. 低影响开发技术的雨洪控制效果及水质影响分析——基于SWMM模型[J]. 中国农村水利水电,2017,411(1):84-87.

[2] 朱嘉祺,徐向阳,何 爽. 基于LH-OAT的SWMM模型参数敏感性分析[J]. 中国农村水利水电,2014,377(3):21-26.

[3] 宋晓猛,张建云,占车生,等. 水文模型参数敏感性分析方法评述[J]. 水利水电科技进展,2015,35(6):105-112.

[4] ZHOU N, PIERRE J W, TRUDNOWSKI D. A stepwise regression method for estimating dominant electromechanical modes[J]. IEEE Transactions on power systems. 2012,27(2):1 051-1 059.

[5] 王海霞,李 昱,周惠成,等. 基于Sobol’敏感性分析法的水库多目标调度模型求解质量与效率[J]. 水电能源科学,2016,34(6):43-47.

[6] 任启伟,陈洋波,舒晓娟. 基于Extend FAST方法的新安江模型参数全局敏感性分析[J]. 中山大学学报(自然科学版),2010,49(3):127-134.

[7] 王浩昌,杜鹏飞,赵冬泉,等. 城市降雨径流模型参数全局灵敏度分析[J]. 中国环境科学,2008,28(8):725-729.

[8] Chen Y, Wen J, Cheng S. Probabilistic load flow method based on nataf transformation and latin hypercube sampling[J]. IEEE Transactions on sustainable energy, 2013,4(2):294-301.

[9] 曹 渊,韩 峰,王铁良,等. 爆炸容器力学安全的敏感度分析方法[J]. 科技导报,2013,31(3):28-32.

[10] 熊剑智. 城市雨洪模型参数敏感性分析与率定[D]. 济南:山东大学,2016.

[11] 侯改娟. 绿色建筑与小区低影响开发雨水系统模型研究[D]. 重庆:重庆大学,2014.