|
【水文泥沙】 基于互信息的SWMM模型参数全局敏感性分析王汉明1，李传奇1，熊剑智1，宋苏林2 （1.山东大学 土建与水利学院，山东 济南250061；2.济南市水文局，山东 济南 250000） 摘 要：采用互信息方法对SWMM模型输入参数与输出结果之间的非线性关系进行了全局敏感性分析。以山东大学千佛山校区降雨径流模拟为例，采用拉丁超立方体抽样，研究了SWMM模型模型参数对输出结果的影响程度。结果表明：峰值流量最敏感的参数为透水区曼宁糙率系数；峰现时间最敏感的参数为管道曼宁糙率系数，其次为最小渗透率；总产流最敏感的参数为最小渗透率。 关键词：信息熵；互信息；全局敏感性分析；参数选择；SWMM模型 Abstract： Mutual information method was used to conduct the overall sensitivity analysis on the complex non⁃linear relationship between the input parameters and outputs of Storm Water Management Model（SWMM） .Taking rainfall⁃runoff simulation in an area of Shandong Univer⁃sity as an example， the Latin Hypercube Sampling （LHS） method was used to sample values for input parameters of SWMM model， and the influence of input parameters to the model results was studied.The results show that the most sensitive parameter for peak discharge is perme⁃able zone Manning value and the most sensitive parameter for peak time is pipeline Manning value， and then the minimum permeability； the most sensitive parameter of total output is the minimum permeability. Key words： information entropy； mutual information； global sensitivity analysis； SWMM Model 随着城市化进程的不断加快，城区不透水区比例增大，城市下垫面条件的改变使得暴雨时产生洪涝灾害的概率显著增大。SWMM（Storm Water Management Model）模型是美国环保局开发的城市暴雨管理模型［1］，该模型可以对区域内降雨及产生的径流输送过程进行模拟，目前已被广泛地应用到城市暴雨洪水管理中。 参数率定是降雨径流模型研究中的重要步骤，通过率定后的参数值可以使模型的结果更加接近实测值，从而使模拟结果更加可靠。与众多的水文模型类似，SWMM模型参数较多，难以率定。若盲目地对全部参数进行调节，则不仅费时费力，而且会出现不同的参数组合产生相同结果的情况，降低模型的可靠性。如何快速、准确地调整参数，使模拟结果更加接近实测值，是水文模型研究的难题。 参数敏感性分析是模型参数率定的辅助方法。通过研究模型参数对模型输出结果的影响，识别关键参数，可为模型的参数率定提供重要的参考。参数敏感性分析方法可分为局部敏感性分析和全局敏感性分析。局部敏感性分析通过控制其他参数取值不变，反映单个参数对模型输出结果的影响；全局敏感性分析则可以反映所有参数对输出结果的影响，适合存在“异参同效”现象的水文模型［2］。常用的敏感性分析方法 有 多 元 回 归 法、 Sobol 法［2］、 Morris 法［3］、 Glue法［4］、互信息法［5］等，其中互信息法可以表示两个变量或多个变量之间共享的信息量，能够很好地刻画变量间的非线性相关关系，适合处理非线性动力学问题。Srikanta Mishra等用互信息法分析了不同参数对地下水模型输入输出的影响，结果表明互信息法可快速、有效地识别模型中的敏感参数。 本研究以山东大学千佛山校区为例，利用拉丁超立方体抽样（LHS）方法对SWMM模型的输入参数进行随机抽样，并采用互信息法对输入参数进行敏感性分析，量化模型参数对模型结果产生的影响，找出对模型结果不确定性影响较大的模型参数，以期为有效提高模型模拟精度提供可靠的依据。 1 模型与方法1.1 SWMM模型及参数 SWMM模型是美国环保局为解决日益严重的城市排水问题所开发的城市暴雨管理模型，此模型可有效分析城市降雨径流问题［6］。根据排水系统的水流特性，SWMM模型分为地表径流和管道输水两部分，分别利用非线性水库模型和圣维南方程来模拟，并通过一定的数值分析规则求解［4］。 SWMM模型中与产汇流计算有关的参数共有14个，其中汇水区面积、汇水区坡度、不透水率、汇水区宽度、管道长度这5个参数可以通过直接测量求得，另外9个参数可分为与汇水区糙率有关的参数（2个）、与洼地蓄积量有关的参数（3个）、与Horton渗透有关的参数（3个）、与管道糙率有关的参数（1个）四大类。本文选取不易直接测量的9个参数作为研究对象进行参数敏感性分析。由于各个参数的概率分布未知，因此假定所有参数都服从均匀分布［7］，参数的取值参考了SWMM模型用户手册，并分析了某些参数在国内的适应情况，对参数的取值范围进行了适当调整，具体数值见表1。 表1 SWMM主要模型参数及取值范围  参数编号 参数类别 参数符号 物理意义 取值范围1汇水区曼宁糙率参数 N－Imperv 不透水区曼宁系数 0.005～0.050 2 N－Perv 透水区曼宁系数 0.05～0.50 3 S－Imperv 不透水区洼蓄量／mm 0.2～10.0 4洼地蓄积量参数S－Perv 透水区洼蓄量／mm 2～10 5 Pct－Zero 不透水区中无洼地面积百分比／% 5～85 6 Max.Infil.Rate 最大渗入速率／（mm·h－1） 26～80 7 Horton渗透参数Min.Infil.Rate 最小渗入速率／（mm·h－1） 0～10 8 Decay 渗透衰减系数／（h－1） 2～7 9管道曼宁糙率参数 Manning－N 管道曼宁系数 0.011～0.024 1.2 互信息法 1.2.1 信息熵 信息熵用来描述随机变量所包含的信息量。当某随机变量值混乱时，人们需要更多的信息描述清楚该随机变量所包含的信息量，此时信息熵就大；当某随机变量值有序时，人们只需要较少的信息就可以描述清楚该随机变量所包含的信息量，此时信息熵就小［4］。信息熵的公式为  式中：H（x）为随机变量X的信息熵；p（x）为随机变量 X为x时发生的概率；b为信息熵的单位。 1.2.2 互信息 互信息表示两个变量或多个变量之间共享的信息量，是信息论中度量随机变量间相互依赖程度的指标，也是一种分析两个变量间相关程度的敏感性分析方法。互信息可以看成是一个随机变量中包含另一个随机变量的信息量，进而引申为两个变量的非线性相关程度［8］。两个离散型随机变量X、Y的互信息可以定义为  式中：I（x，y）为两个变量的互信息；p（y）为随机变量Y为y时发生的概率；p（x，y）为随机变量X和Y的联合概率。 互信息、联合熵和信息熵的关系见图1。  图1 互信息、联合熵和信息熵的关系 假如变量 X、Y 互相独立，则 I（x，y） ＝ 0。 互信息还具有对称性，I（x，y） ＝ I（y，x）。 当变量 X、Y 完全相关时，I（x，y） ＝ 0.5［H（x） ＋ H（y）］。 U 和 R 是基于互信息量度相关性的两个指标［9－10］，定义为  U和R的取值范围为0～1，U和R的值越大说明变量之间的非线性相关性越强。 2 算例与分析2.1 研究区概况 山东大学千佛山校区南院位于济南市经十路与经十一路之间，千佛山路以西，舜耕路以东。地势南高北低，东高西低，地面最高处海拔为59.95 m，最低处海拔为53.95 m，最大高差为6 m。受地形影响，该校区降雨时无其他地区径流汇入，是一片独立的区域。本项目研究区域总面积为79 887 m2（约8 hm2），根据研究内容可将下垫面分为透水区和不透水区。透水区为草坪、花坛等，总面积12 890 m2，约占总面积的16%；不透水区包括沥青路面、建筑屋面以及操场等混凝土地面，总面积66 997 m2，约占总面的84%，其中沥青路面总面积12 656 m2，建筑屋面总面积11 300 m2，操场等其他混凝土地面面积43 041 m2。 整个模拟区域共划分子汇水区域24个、节点23个、管道23个（其中9条为地下管线、14条为道路，共有3个排放口，分别为 PFK1、PFK2、PFK3，由于地形以及地下排水管道布置较少，因此该校区70%以上降雨形成的径流主要沿道路由PFK3排出，见图2（SUB为子汇水区，JUN为节点，GD为管道，PFK为排放口）。  图2 SWMM模型研究区域 2.2 设计降雨 本次模拟降雨强度采用济南市1 a一遇设计降雨的雨强，降雨过程使用芝加哥降雨过程线［11］。根据2014年最新的全国各城市暴雨强度公式，求得济南市1 a一遇降雨强度，利用芝加哥降雨过程线模型求得时间间隔为6 min的2 h降雨过程线，见图3。  图3 济南市1 a一遇设计降雨2 h过程线 我国暴雨强度公式一般采用式（5）计算：  式中：q为平均暴雨强度，mm／min 或 L／（s·hm2）；A1为重现期为1 a的设计降雨量；C为雨量变动参数，是反映设计降雨各历时不同重现期强度变化程度的参数之一；P 为设计降雨重现期，a；td为暴雨历时，min；b、c为常数，共同反映重现期的设计降雨随历时延长其强度递减变化的情况。 式（5） 中的经验参数A1、C、b、c可直接通过地方降雨强度公式确定，于是式（5）可改写为 q＝ 2.3 拉丁超立方体抽样 参数敏感性分析依赖于统计学，需要对参数进行抽样。拉丁超立方体抽样方法是一种分层随机抽样方法，通过划分概率相等的间隔，从每个间隔中选取一个样本点。由于样本点间具有独立性，因此不需要更多维度的样本就可以得到相同的效果。采用LHS法进行随机抽样，生成1 000组SWMM模型参数，运行1 000次SWMM模型，得到1 000组输出结果，量化每个参数对模拟结果不确定性的贡献率。 2.4 互信息的Matlab实现 进行互信息计算前，首先用抽样得到的1 000组参数替换SWMM模型中.inp文件的参数，生成1 000个.inp文件，接着分别利用SWMM模型运算得到1 000个.rpt输出文件，利用Matlab提取输出文件所需的输出结果。用互信息法分析模型参数敏感性。 SWMM模型运算最重要的3个输出结果分别为峰值流量、峰现时间和总产流。本次模拟区域内的径流主要经过PFK3排出，因此利用Matlab从1 000个SWMM输出文件中提取PFK3处的峰值流量、峰现时间以及总产流信息。 本文采用自然常数e为对数底数，把9个模型参数作为自变量X，把模型的输出结果作为因变量Y。自变量X可以分成i等份，同时因变量Y可以分成j等份。Pi为X中第i份的数量占总数量的概率，Pi＝Ni／N，Pj为Y中第j份的数量占总数量的概率，Pj＝Nj／N。变量X、Y可以组成联立表，联立表是由一个事件中每个组合的行和列形成的整数表，继而得到两个变量的联合概率，Pij为在变量X第i等份取值时Y中第j份因变量的个数占总数量的概率，Pij＝Nij／N。根据以上定义，信息熵和互信息的公式分别为  利用Matlab求解互信息的程序：首先利用9个抽样参数与3个提取结果（峰值流量、峰现时间和总产流）分别生成1 000×10的矩阵，利用sortrows函数对9个抽样参数分别进行秩排序，并用hist函数把参数分成10等份，把3个输出结果也相应分成10等份，求得各个参数的信息熵H（x）和输出变量的信息熵H（y）；利用zeros函数生成10×10的零阶矩阵，并把各个参数第i份取值时第j份输出结果中的个数占总数量的概率赋值到这个零阶矩阵，计算联合熵H（x，y）；根据求得的信息熵 H（x） 、 H（y） 以及 H（x，y） ，求得互信息I（x，y） 的大小，并根据式（3）、式（4）得到互信息的两个相关性指标U和R。 2.5 参数敏感性分析 利用Matlab对选取的模型参数和输出结果之间的关系进行互信息分析，得到统计结果。由结果可知输出结果对各个参数的敏感性不尽相同。利用气泡图展现选取的模型参数与输出结果的关系，可使互信息的结果更加直观。首先对选取的9个模型参数分别进行秩排序，并且等分成10份，相应的输出结果也等分成10份，每个点的大小表示参数在第i等份取值时，第j份输出结果中的个数占总数量的概率。气泡越大概率越大，形成一个10×10的气泡图。 由表2可知，参数N－Perv（透水区曼宁系数）的U值为 0.228 2、R 值为0.605 3，其他 8 个参数的 U 值为0.019 7～0.042 8，R 值为 0.196 7～0.286 5。 由图 4 可知，参数N－Perv与输出结果（峰值流量）的相关性较强，其他参数概率分布比较离散。因此，峰值流量最敏感的参数是N－Perv，其他参数则比较平衡，差异较小。这说明参数N－Perv对峰值流量起决定性作用。 表2 峰值流量互信息分析结果  参数编号 参数名 H（x，y） I（x，y） U（x，y） R（x，y）1 N－Imperv 4.109 7 0.086 9 0.041 4 0.282 0 2 N－Perv 3.717 8 0.478 8 0.228 2 0.605 3 3 S－Imperv 4.142 5 0.054 1 0.025 8 0.224 3 4 S－Perv 4.145 8 0.050 8 0.024 2 0.217 3 5 Pct－Zero 4.155 2 0.041 4 0.019 7 0.196 7 6 Max.Infil.Rate 4.115 5 0.081 1 0.038 7 0.272 7 7 Min.Infil.Rate 4.148 8 0.017 8 0.022 8 0.211 1 8 Decay 4.122 3 0.074 3 0.035 4 0.261 5 9 Manning－N 4.106 6 0.089 9 0.042 8 0.286 5 图4 峰值流量气泡图 由表3可知，峰现时间最敏感的参数为Manning－N（管道曼宁糙率系数），其 R值为 0.580 4、U值为0.205 4；其次为 Min.Infil.Rate（最小渗入速率），其 R值为0.341 5、U 值为0.062 0；其他参数的 R 值均在0.3以下、U值均小于0.05。由图5可知，参数Manning－N取值较小时，峰现时间的计算结果普遍偏小，随着Manning－N取值增大，峰现时间的计算结果增大，而其余参数的气泡分布比较离散，参数取值变化对峰现时间计算结果影响较小。 表3 峰现时间互信息分析结果 编号 参数名 H（x，y） I（x，y） U（x，y） R（x，y）1 N－Imperv 4.237 5 0.101 1 0.046 6 0.298 3 2 N－Perv 4.301 3 0.037 4 0.017 2 0.184 0 3 S－Imperv 4.280 6 0.058 1 0.026 8 0.228 4 4 S－Perv 4.284 7 0.054 0 0.024 9 0.220 3 5 Pct－Zero 4.284 0 0.054 6 0.025 2 0.221 6 6 Max.Infil.Rate 4.270 0 0.068 7 0.031 7 0.247 7 7 Min.Infil.Rate 4.204 2 0.134 5 0.062 0 0.341 5 8 Decay 4.277 4 0.061 2 0.028 2 0.234 3 9 Manning－N 3.893 1 0.445 5 0.205 4 0.580 4 图5 峰现时间气泡图 由表4可知，总产流最敏感的参数为Min.Infil.Rate（最小渗入速率），其 R 值为 0.506 3，U 值为0.148 1。 3个参数 N－Imperv（不透水区曼宁糙率系数）、N－perv（透水区曼宁糙率系数）、Manning－N（管道曼宁糙率系数）的R值均在0.2以下、U值均在0.02以下，尤其是参数Manning－N的 R值为仅0.089 0，U值仅为0.010 2。这说明3个曼宁粗糙率参数对总产流影响都较小。由图6可知，参数Min.Infil.Rate取值较小时，总产量的计算结果普遍偏小，随着参数Min.Infil.Rate取值增大而增大，其余参数的气泡分布比较离散，取值变化对总产流计算结果的影响较小。 表4 总产流互信息分析结果 参数编号 参数名 H（x，y） I（x，y） U（x，y） R（x，y）1 N－Imperv 4.214 3 0.038 7 0.018 2 0.189 0 2 N－Perv 4.209 7 0.040 3 0.019 9 0.194 2 3 S－Imperv 4.196 6 0.056 4 0.026 5 0.227 3 4 S－Perv 4.200 4 0.052 6 0.024 7 0.219 7 5 Pct－Zero 4.195 1 0.058 0 0.027 3 0.230 3 6 Max.Infil.Rate 4.185 1 0.068 0 0.032 0 0.248 9 7 Min.Infil.Rate 3.938 1 0.314 9 0.148 1 0.506 3 8 Decay 4.181 3 0.071 7 0.033 7 0.255 4 9 Manning－N 4.228 6 0.023 5 0.010 2 0.089 0 图6 总产流气泡图 3 结 论（1）采用信息熵和互信息法对SWMM模型参数进行了敏感性分析，结果显示峰值流量最敏感的参数为透水区曼宁糙率系数；峰现时间最敏感的参数为管道曼宁糙率系数，其次为最小渗透率；总产流最敏感的参数为最小渗透率。 （2）利用气泡图展现选取的模型参数与输出结果的关系，可使敏感性分析结果更加直观清晰。 （3）互信息法能够较好地识别出敏感参数，尤其适用于复杂的非线性水文模型的参数敏感性分析，可减少模型率定过程中参数调整的盲目性。 参考文献： ［1］ HUBER W C， DICKINSON R E.Storm Water Management Model User′s Manual， （Version 4） ［M］.Georgia： Environ⁃mental Protection Agency， 1992：1－10. ［2］ 任启伟，陈洋波，周浩澜，等.基于Sobol法的TOPMODEL模型全局敏感性分析［J］.人民长江，2010，41（19）：91－94. ［3］ 史蓉，庞博，赵刚，等.SWMM模型在城市暴雨洪水模拟中的参数敏感性分析［J］.北京师范大学学报，2014，50（5）：456－460. ［4］ 赵冬泉，王浩正，陈吉宁，等.城市暴雨径流模拟的参数不确定性研究［J］.水科学进展，2009，20（1）：45－51. ［5］ WESTRA S， BROWN C， LALL U， et al.Modeling Multiva⁃riable Hydrological Series：Principal Component Analysis or Independent Component Analysis［J］.Water Resources Re⁃search， 2007， 43（6）：W064296. ［6］ 董楠.基于GIS的城市排水管网系统模拟研究［D］.天津：天津大学，2009：18－30. ［7］ 王浩昌，杜鹏飞，赵冬泉，等.城市降雨径流模型参数全局灵敏度分析［J］.中国环境科学，2008，28（8）：725－729. ［8］ 龚伟.基于信息熵和互信息的流域水文模型不确定性分析［D］.北京：清华大学，2012：19－28. ［9］ LI Chuanqi， WANG Wei， XIONG Jianzhi.Sensitivity Anal⁃ysis for Urban Drainage Modeling Using Mutual Information［J］ .Entropy， 2014，16（11）： 5738－5752. ［10］ SRIKANTA Mishra， NEIL Deeds， GREG Ruskauff.Global Sensitivity Analysis Techniques for Probabilistic Ground Water Modeling ［J］.Ground Water， 2009， 47（5）：727－744. ［11］ 张大伟，赵冬泉，陈吉宁，等，芝加哥降雨过程线模型在排水系统模拟中的应用［J］.给水排水，2008，34（增刊）：354－358. 【责任编辑 翟戌亮】 Analysis of SWMM Parameters Global Sensitivity Based On Mutual Information WANG Hanming1， LI Chuanqi1， XIONG Jianzhi1， SONG Sulin2 中图分类号：TV122 文献标志码：A doi：10.3969／j.issn.1000－1379.2017.10.005 收稿日期：2016－09－13 基金项目：水利部公益性行业科研专项（201401003，201301057）；山东省重大水利科研与技术推广基金资助项目（SDSLKY201403）。 作者简介：王汉明（1990—），男，山东聊城人，硕士研究生，研究方向为城市雨洪。 通信作者：李传奇（1963—），男，河南沈丘人，教授，博士生导师，主要研究方向为水环境与水生态、城市雨洪模拟及预报。 E⁃mail：76052514＠ qq.com |
|
|
