基于Morris、Sobol和EFAST的LID設(shè)施模型參數(shù)全局敏感性分析

金夢瀟，田 勇，Michele Lancia，于 江，鄭春苗

（1.哈爾濱工業(yè)大學(xué)環(huán)境學(xué)院，哈爾濱 150090；2.南方科技大學(xué)環(huán)境科學(xué)與工程學(xué)院地表水-地下水污染綜合治理國家環(huán)境保護(hù)重點(diǎn)實(shí)驗(yàn)室，廣東 深圳 518055；3.北京大學(xué)工學(xué)院水科學(xué)中心，北京 100871）

0 引 言

不斷加快的城市化進(jìn)程增加了不透水面的覆蓋率，影響了自然水文循環(huán)過程，進(jìn)而威脅到城市水生態(tài)環(huán)境。我國在2014年首次提出“海綿城市”的管理理念來解決城市化導(dǎo)致的一系列水環(huán)境問題［1，2］。海綿城市的建設(shè)重點(diǎn)是結(jié)合低影響開發(fā)（Low Impact Development，LID）的理念［3］，通過分散的、源頭的設(shè)施使雨水就地下滲，控制徑流量等，來達(dá)到調(diào)蓄雨水的目的。LID設(shè)施包括雨水花園、透水鋪裝、生物滯留池、植被草溝、下沉綠地等，不同的LID設(shè)施對增加雨水下滲、控制徑流總量和削減徑流峰值有著不同的效能。為了評估各類LID 設(shè)施的效能，可利用SWMM 構(gòu)建包含LID 設(shè)施的水文模型，模擬LID 設(shè)施的水文過程及控制效果，這是海綿城市規(guī)劃設(shè)計(jì)中不可或缺的重要環(huán)節(jié)。

在海綿城市規(guī)劃設(shè)計(jì)中，LID 設(shè)施模型的參數(shù)眾多。識別參數(shù)的敏感性，對篩選LID 設(shè)施及篩選后進(jìn)行的SWMM 模擬計(jì)算和結(jié)果分析具有重要意義。敏感性分析是從眾多輸入?yún)?shù)中識別對輸出結(jié)果有顯著影響參數(shù)的方法，對提高模型應(yīng)用效率有重要幫助［4］。敏感性分析方法眾多，大致可分為局部敏感性分析和全局敏感性分析兩大類。局部敏感性分析［5］是研究在局部范圍內(nèi)，單個(gè)輸入變量（其他變量保持不變）變化對輸出結(jié)果的影響，常用方法包括OAT 法、DA 法。局部敏感性分析原理簡單、計(jì)算量小，但不能充分描述模型參數(shù)的空間分布形態(tài)，忽略了參數(shù)之間的相互作用。全局敏感性分析［6］是研究多個(gè)輸入變量同時(shí)變化對模型輸入結(jié)果的影響，此種方法把輸入變量的變化范圍擴(kuò)展到多個(gè)參數(shù)的定義域，可得到不同參數(shù)之間的相互作用，常用方法包括Morris 法、Sobol 法、EFAST 法等［7］。由于LID 設(shè)施模擬涉及復(fù)雜的非線性水文過程，因此本研究選用三種全局敏感性分析方法對LID 設(shè)施的各參數(shù)進(jìn)行對比研究，為提高模型優(yōu)化及LID設(shè)施效能提供技術(shù)支撐。

1 研究方法

1.1 SWMM模型

SWMM（Storm Water Management Model）是由美國環(huán)境保護(hù)署（Environmental Protection Agency，EPA）開發(fā)的一款可以進(jìn)行城市雨洪管理的分布式水文模型。此模型一般多應(yīng)用于模擬城市某一單一或長期降雨事件的地表徑流水量和水質(zhì)變化［8］，近年來有很多研究將SWMM 模型應(yīng)用于海綿城市的規(guī)劃和設(shè)計(jì)［9，10］，以及分析和評估LID 設(shè)施對不同頻率城市雨洪的徑流控制效果［11-13］。

1.1.1 SWMM 模型構(gòu)建

研究區(qū)域位于深圳市寶安區(qū)新橋河流域，流域總面積約為34.67 km2，地勢東南高、西北低，雨水從東南向西北匯集，最后通過管道流入新橋河。該區(qū)域土地利用主要以居住用地和工業(yè)用地為主，不透水面面積較大，約占總面積的67%；而以林地為主的透水面面積則占總面積的33%。流域內(nèi)有1個(gè)氣象觀測站。圖1為研究區(qū)域的SWMM 概化圖。依據(jù)實(shí)際的勘察資料，將研究區(qū)域劃分為92 個(gè)子匯水區(qū)，112 個(gè)雨水收集管道。基于收集的90 m分辨率數(shù)字高程地圖和土地利用類型，使用ArcGIS軟件計(jì)算每個(gè)子匯水區(qū)的面積、坡度、不透水比例和排水管道長度等參數(shù)，作為SWMM模型的輸入。

圖1 研究區(qū)域SWMM概化圖Fig.1 SWMM model of the study area

模擬采用長期降雨事件。降雨數(shù)據(jù)為深圳市氣象局提供的萬豐站日降雨量數(shù)據(jù)，監(jiān)測日期為2018.01.16-2018.10.16（9個(gè)月），總降雨量為1 779 mm。

1.1.2 SWMM 模型校正

根據(jù)實(shí)地勘測資料和SWMM 手冊里的參數(shù)參考范圍，確定模型參數(shù)的初始數(shù)值，對不透水面和透水面的曼寧系數(shù)（N-Imperv、N-Perv）和洼地蓄水（Dstore-Imperv、Dstore-Perv），以及Horton 下滲公式的相關(guān)參數(shù)進(jìn)行調(diào)參。分別選取2018.04.04-2018.04.23（20 d）和2018.06.04—2018.06.23（20 d）進(jìn)行模擬，通過監(jiān)測點(diǎn)河水水位的模擬值和監(jiān)測值的對比進(jìn)行調(diào)參。最終模擬結(jié)果如圖2所示，模擬值和監(jiān)測值的河水水位基本相同，波動(dòng)情況基本相似。

圖2 模擬值和監(jiān)測值對比Fig.2 Comparison of simulated and monitored values

1.1.3 LID設(shè)施模型參數(shù)

本研究在新橋河流域現(xiàn)狀條件基礎(chǔ)上開展了海綿城市設(shè)計(jì)（圖3），在校正后的SWMM 上加入相應(yīng)的LID 設(shè)施模型，來模擬海綿城市水文過程并進(jìn)行參數(shù)敏感性分析。LID 設(shè)施模型是SWMM 里的一個(gè)模塊［10］，此模塊可以模擬不同類型LID 設(shè)施在暴雨情景下的蓄滯、蒸發(fā)、下滲及排水過程。在SWMM 里，將LID 設(shè)施垂直方向概化成地表層、土壤層、蓄水層及相關(guān)排水系統(tǒng)，不同的LID 設(shè)施進(jìn)行不同層的組合［14］。本研究主要選擇了4 種LID 設(shè)施：植被草溝（下沉綠地）、綠色屋頂、生物滯留池（雨水花園）和透水鋪裝。根據(jù)各匯水區(qū)的土地功能及LID 設(shè)施的地形需求，通過海綿城市規(guī)劃設(shè)計(jì)的方法，將4 種LID 設(shè)施分配設(shè)置在每個(gè)子匯水區(qū)。規(guī)劃完成后每種LID設(shè)施的面積及在該區(qū)域海綿城市設(shè)計(jì)面積中的占比分別為，植被草溝2.472 km2（37.9%）；綠色屋頂2.436 km2（37.2%）；生物滯留池1.178 km2（18.0%）；透水鋪裝0.451 km2（6.9%）。

圖3 研究區(qū)域的海綿城市設(shè)計(jì)圖Fig.3 Design of sponge city in the study area

4 種LID 設(shè)施的部分模型參數(shù)采用設(shè)計(jì)值，例如土壤地表層下凹深度、地表坡度等。本次研究選擇地表層和土壤層的模型參數(shù)進(jìn)行敏感性分析（而其設(shè)計(jì)值參數(shù)不予考慮）。根據(jù)篩選結(jié)果，確定4 種LID 設(shè)施的25 個(gè)參數(shù)（表1）。各參數(shù)的取值范圍參考SWMM 手冊，其中SWMM 規(guī)定土壤的孔隙率應(yīng)大于田間持水量，田間持水量應(yīng)大于凋萎點(diǎn)。

表1 LID設(shè)施參數(shù)及取值范圍Tab.1 Parameters of LID facilities and their value ranges

1.2 敏感性分析方法

全局敏感性分析方法可分為兩種，定性分析和定量分析。定性分析是研究不同參數(shù)對輸出變量響應(yīng)的相對大小，如Morris 法、Delta Test 法等；而定量分析可獲得不同參數(shù)對輸出變量的貢獻(xiàn)率大小，如Sobol 法、EFAST 法等。本研究同時(shí)選用Morris、Sobol 和EFAST 三種敏感性分析方法開展研究，分別從定性和定量兩個(gè)方面揭示SWMM 模型結(jié)果對LID 設(shè)施參數(shù)的敏感性。

1.2.1 Morris法

Morris 法是目前廣泛應(yīng)用的定性全局敏感性分析法［15］，由Max D.Morris 于1991年提出，之后被Campolonge 等［16］改進(jìn)。Morris 法可以快速篩選出敏感的參數(shù)，但無法量化各參數(shù)的敏感性。該方法由單獨(dú)隨機(jī)的“一次變化法”實(shí)驗(yàn)組成，在這些實(shí)驗(yàn)中，依次評估改變每個(gè)被選擇因素的值所產(chǎn)生的影響［17］。

Morris法計(jì)算公式為：

式中：Si（X，Δx）為參數(shù)i的敏感性指數(shù)；y（x）為模型輸出結(jié)果；X=（x1，x2，…，xn）為參數(shù)的n維向量；Δx為x的變化量。

Morris法首先把各參數(shù)的取值范圍映射到［0，1］之間，構(gòu)成n維采樣空間。在采樣空間里，應(yīng)用一次變化法的試驗(yàn)設(shè)計(jì)對參數(shù)進(jìn)行隨機(jī)采樣。由于采樣方法具有隨機(jī)性，需進(jìn)行t次重復(fù)，得到更好的采樣分布。通過對采樣樣本進(jìn)行計(jì)算，可得到各參數(shù)Si的平均值μ和標(biāo)準(zhǔn)差σ。μ表示輸入變量對輸出變量影響的強(qiáng)弱，μ值越大則表明輸入變量的敏感性越強(qiáng)。σ表示輸入變量間相互作用的強(qiáng)弱，σ值越大則表明參數(shù)間的相互作用越強(qiáng)。

1.2.2 Sobol法

Sobol 法由Sobol 于1993年提出［18］，其核心是將模型輸入的總方差分解成各參數(shù)的方差和參數(shù)相互作用的方差之和，通過計(jì)算各參數(shù)對輸出方差的貢獻(xiàn)比來識別各參數(shù)的敏感性［19］。

假設(shè)模型為Y=f（X）=f（x1，x2，…，xn），X=［x1，x2，…，xn］是參數(shù)向量，n為參數(shù)的數(shù)量，方差分解公式表示為：

式中：V為模型的總方差；Vi為第i個(gè)參數(shù)作用的方差項(xiàng)；Vi，j為第i個(gè)和第j個(gè)參數(shù)共同作用的方差項(xiàng)；V1，2，…，n為所有參數(shù)共同作用的方差項(xiàng)；Si為第i個(gè)參數(shù)的敏感性指數(shù)。

1.2.3 EFAST法

EFAST（Extended Fourier Amplitude Sensitivity Test）稱為擴(kuò)展傅里葉振幅靈敏度檢驗(yàn)法［20］，是在FAST 法的基礎(chǔ)上結(jié)合Sobol 法改進(jìn)而來的全局敏感性分析方法。此方法采取了模型方差分析的思想，計(jì)算獲得模型預(yù)測的期望值和方差，來表達(dá)各個(gè)輸入?yún)?shù)間的相互作用［21］。因此，通過模型方差的分解得到各參數(shù)間的相互作用對總方差的貢獻(xiàn)比，進(jìn)而得到各參數(shù)的敏感性指數(shù)。

式中：V為模型的總方差；Vi為第i個(gè)參數(shù)作用的方差項(xiàng)；Vi，j為第i個(gè)和第j個(gè)參數(shù)共同作用的方差項(xiàng)；Si為第i個(gè)參數(shù)的敏感性指數(shù)。

1.3 敏感性分析過程

采用SimLab2.2軟件對LID設(shè)施的參數(shù)敏感性進(jìn)行分析（圖4），具體過程如下：

圖4 敏感性分析的流程圖Fig.4 Flow chart for sensitivity analysis

（1）定義LID 設(shè)施各參數(shù)的取值范圍及分布形式。根據(jù)SimLab 軟件的要求，將表1里各參數(shù)的取值范圍轉(zhuǎn)換到［0，1］區(qū)間內(nèi)，各參數(shù)分布形式均假設(shè)為均勻分布。

（2）使用SimLab 軟件在多維參數(shù)空間內(nèi)進(jìn)行隨機(jī)采集，生成參數(shù)樣本集。

（3）編寫MATLAB 程序，將生成的參數(shù)樣本集寫入對應(yīng)的SWMM模型文件中，運(yùn)行模型。

（4）使用MATLAB 程序讀取SWMM 模型結(jié)果，并將模型結(jié)果編寫成SimLab 軟件可讀取的文件格式。通過SimLab 軟件的Monte Carlo 方法進(jìn)行敏感性計(jì)算，最終得到各種方法的分析結(jié)果。

2 結(jié)果與分析

2.1 Morris法分析結(jié)果

采用Morris 法對25 個(gè)LID 設(shè)施模型參數(shù)抽樣10 次，對得到的260 組樣本量進(jìn)行敏感性分析計(jì)算，獲得4 種不同輸出變量的敏感性指數(shù)（μ，σ）散點(diǎn)圖（圖5）。圖5（a）至5（d）分別顯示了以蒸發(fā)量、雨水下滲、地表徑流量和峰值流量為響應(yīng)變量的敏感性分析結(jié)果。

對蒸發(fā)［圖5（a）］而言，敏感性最高的參數(shù)為生物滯留池土壤導(dǎo)水率（B3），其次為植被草溝植物覆蓋率（V1）、綠色屋頂導(dǎo)水率（G3）、植被草溝粗糙系數(shù)（V2）、綠色屋頂滲流率（G6）；而植被覆蓋率（V1）與各輸入?yún)?shù)的相互作用最強(qiáng)。雨水下滲［圖5（b）］對生物滯留池相關(guān)參數(shù)的敏感性最高，而對其他LID設(shè)施的參數(shù)敏感性較低，其中最為敏感及參數(shù)相互作用最顯著的3個(gè)參數(shù)均為生物滯留池的導(dǎo)水率（B3）、植物覆蓋率（B1）和滲流率（B6）。地表徑流量［圖5（c）］對生物滯留池導(dǎo)水率（B3）的敏感性最高，其次為透水鋪裝導(dǎo)水率（P3）和生物滯留池植物覆蓋率（B1）。峰值流量［圖5（d）］對生物滯留池植物覆蓋率（B1）最為敏感，其次為生物滯留池粗糙系數(shù)（B2）和綠色屋頂植物覆蓋率（G1），參數(shù)的敏感性與參數(shù)相互作用的排序一致?？傊?，從4組不同輸出變量敏感性分析得出，生物滯留池導(dǎo)水率對模型輸出變量的影響較大，其次是生物滯留池的植物覆蓋率和粗糙系數(shù)。

圖5 Morris法敏感性分析結(jié)果Fig.5 Sensitivity analysis results of the Morris method

2.2 Sobol法分析結(jié)果

選用Sobol 法進(jìn)行定量的敏感性分析，同樣模型參數(shù)為25個(gè)，生成1 508 組參數(shù)樣本，得到敏感性分析結(jié)果（圖6）。與Morris 法的結(jié)果相比，蒸發(fā)、雨水下滲、地表徑流量［圖6（a）、圖6（b）、圖6（c）］均對生物滯留池導(dǎo)水率（B3）最為敏感。而其他參數(shù)對蒸發(fā)、雨水下滲、地表徑流量的響應(yīng)較小。對峰值流量［圖6（d）］而言，Sobol 法中最敏感的參數(shù)是生物滯留池導(dǎo)水率（B3），Morris 法中最敏感的參數(shù)為生物滯留池植物覆蓋率（B1），兩者的分析結(jié)果不同。這可能由于Morris 法的采樣次數(shù)較少，或者Sobol 法考慮參數(shù)之間的相互作用影響較大，引起兩者分析結(jié)果之間的差異，這種差異原因在其他文獻(xiàn)中也有所討論［22］。但是從Sobol定量分析的結(jié)果來看，B1和B3對峰值流量影響均較大。Sobol 法的計(jì)算可作為Morris 法的一種補(bǔ)充與印證。

圖6 Sobol法敏感性分析結(jié)果Fig.6 Sensitivity analysis results of the Sobol method

2.3 EFAST法分析結(jié)果

為了與Sobol 法進(jìn)行定量比較，本研究同時(shí)采用EFAST 法進(jìn)行敏感性分析。樣本數(shù)量為2 150組，得到不同LID設(shè)施參數(shù)的敏感性分析結(jié)果（圖7）。在EFAST 法結(jié)果分析中，蒸發(fā)［圖7（a）］，最敏感的參數(shù)是生物滯留池導(dǎo)水率（B3），植被草溝植物覆蓋率（V1）、綠色屋頂田間持水量（G5）以及透水鋪裝粗糙系數(shù)（P2）也比較敏感。雨水下滲［圖7（b）］和地表徑流量［圖7（c）］，最敏感的參數(shù)是生物滯留池植物覆蓋率（B1），與Morris 法和Sobol 法的分析結(jié)果相同。而峰值流量［圖7（d）］，最敏感的參數(shù)是生物滯留池植物覆蓋率（B1），次之是生物滯留池粗糙系數(shù)（B2），此結(jié)果與Morris 法的結(jié)果一致，與Sobol 法的結(jié)果不同。造成兩種結(jié)果差異的原因，一方面由于采樣數(shù)量的影響，另一方面由于參數(shù)間相互耦合造成的影響。

圖7 EFAST法敏感性分析結(jié)果Fig.7 Sensitivity analysis results of the EFAST method

2.4 EFAST法、Morris法和Sobol法比較

表2比較了3 種全局敏感性分析方法篩選出的LID 設(shè)施對輸出變量影響最大的4 個(gè)參數(shù)。從表2可以看出，三者分析的參數(shù)敏感性極為相似，但是排序有所不同。原因一是分析方法的樣本數(shù)不同，Morris法的樣本量較小，而兩種定量分析法的樣本量，計(jì)算量均較大。原因二是各輸出變量的敏感性分析結(jié)果受各參數(shù)之間的相互作用的強(qiáng)弱影響，在此方面EFAST 方法優(yōu)勢性更強(qiáng)。與Morris法相比，Sobol法和EFAST法可以進(jìn)一步的定量分析參數(shù)對輸出變量的貢獻(xiàn)率。EFAST 法的計(jì)算更穩(wěn)健高效，因此在模型參數(shù)較多時(shí)，首選EFAST 法。而參數(shù)少，樣本數(shù)量少的情況下，首選Morris法進(jìn)行參數(shù)的定性分析。

3 結(jié) 論

通過SWMM 模型結(jié)合3 種全局敏感性分析方法，從蒸發(fā)、雨水下滲、地表徑流量及峰值流量4方面綜合評估4種LID設(shè)施中25個(gè)模型參數(shù)的敏感性。研究表明，生物滯留池導(dǎo)水率對蒸發(fā)量、雨水下滲和地表徑流量影響最顯著，生物滯留池植物覆

蓋率和導(dǎo)水率對峰值流量影響最顯著。本研究采用不同全局敏感性分析方法的結(jié)果具有一致性，因此證明對LID 設(shè)施參數(shù)的篩選具有可行性，對LID 設(shè)施的配置和海綿城市的優(yōu)化設(shè)計(jì)可起到指導(dǎo)性作用。3 種敏感性分析方法中，定性法（Morris）、定量法（Sobol/EFAST）具有不同的優(yōu)缺點(diǎn)，Morris 計(jì)算樣本少，計(jì)算速度快；Sobol 法和EFAST 法分析精準(zhǔn)，但是樣本需求量大，計(jì)算時(shí)間長。多種方法的應(yīng)用可使敏感性分析更具有可靠性和準(zhǔn)確性。