基于SUFI-2參數(shù)最適置信區(qū)間的參數(shù)優(yōu)化方法結(jié)果對(duì)比分析

陳海濤，曹向真

（華北水利水電大學(xué)水資源學(xué)院，鄭州450045）

0 引言

SWAT模型作為分布式水文模型，由于其可很好的模擬和反應(yīng)不同水文要素和下墊面分布的不均勻性而在流域水文平衡、長(zhǎng)期的地表徑流模擬等方面得到了廣泛的應(yīng)用[1-3]。許多研究發(fā)現(xiàn)氣候和土地變化對(duì)地表徑流的形成具有很大的影響[4,5]。近年來基于SWAT 模型對(duì)研究區(qū)域的地表徑流模擬演變規(guī)律、大尺度農(nóng)業(yè)面源污染負(fù)荷以及城市水文響應(yīng)的分析逐漸增多[6-9]。

SWAT 模型是在國(guó)外的數(shù)據(jù)和水文試驗(yàn)基礎(chǔ)上建立的模型，直接運(yùn)用在國(guó)內(nèi)流域會(huì)出現(xiàn)不理想的情況[10,11]，因此要將模型運(yùn)用到中國(guó)流域，需要進(jìn)行完整的模型驗(yàn)證評(píng)估，水文模型參數(shù)的不確定性研究及參數(shù)率定、校準(zhǔn)問題一直是需要關(guān)注的重點(diǎn)[12,13]。水文參數(shù)的不確定性研究常采用SWAT 模型自帶的參數(shù)敏感性分析或是SWAT-CUP 中的LH-OAT 方法[14,15]。在參數(shù)率定方面，白薇等通過對(duì)SWAT 模型自動(dòng)率定模塊的改進(jìn)和應(yīng)用[16]，提出了一種快速率定SWAT 模型的方法，該方法只對(duì)部分子流域進(jìn)行模擬，大大減少了自動(dòng)率定所需要的時(shí)間。周帥等以黃河上游流域?yàn)槔骄苛怂哪Ｐ蛥?shù)獨(dú)立以及交互作用對(duì)徑流模擬的影響[17]，結(jié)果表明，參數(shù)間的相互作用在徑流模擬中具有很大的不確定性，對(duì)汛期的影響尤為顯著。參數(shù)校準(zhǔn)方法大多采用自動(dòng)校準(zhǔn)方法，即SWAT-CUP 自帶的四種優(yōu)化參數(shù)算法（SUFI-2、GLUE、PSO、Parasol），目前這些算法已被大量使用在水文模型的研究中。魏丹等構(gòu)建分布式水文模型并通過改進(jìn)的SUFI-2 算法對(duì)模型參數(shù)進(jìn)行了不確定分析，探討了SWAT模型在干旱和半干旱流域的適應(yīng)性[18]。董天奧以河西支流域?yàn)檠芯繀^(qū)域，引入了PSO 優(yōu)化算法，實(shí)現(xiàn)了SWAT 參數(shù)的批量?jī)?yōu)化[19]。薛晨等[20]則是采用SUFI-2和GLUE兩種算法對(duì)我國(guó)霍林河流域的SWAT模型參數(shù)進(jìn)行了不確定分析，比較了不同參數(shù)取值對(duì)徑流模擬的影響。樊琨[21]等運(yùn)用手動(dòng)和自動(dòng)兩種方法對(duì)三水河流域進(jìn)行徑流模擬和調(diào)參，結(jié)果表明SWAT-CUP 方法效率高但穩(wěn)定性差，結(jié)合人工干預(yù)后可獲得更好的結(jié)果。但四種算法在同一流域的適用性比較研究還比較少，以上研究也大多沒有考慮參數(shù)置信區(qū)間差異對(duì)四種校準(zhǔn)方法所得精度結(jié)果的影響。

本文以湖北麻城舉水河流域?yàn)槔?，通過SWAT建立月尺度徑流模型，并采用SUFI-2 方法得到參數(shù)適宜置信區(qū)間，比較參數(shù)置信區(qū)間影響前后的3 種優(yōu)化算法（GLUE、PSO、Parasol）徑流結(jié)果并進(jìn)行適應(yīng)性分析。通過探究參數(shù)置信區(qū)間不確定性對(duì)徑流模擬的影響，為模型參數(shù)率定及優(yōu)化算法的選擇提供了一種有效的方法。

1 研究區(qū)域應(yīng)用

1.1 研究區(qū)域概況

舉水河流域位于湖北省東北部，大別山南麓，長(zhǎng)江中游下段北岸，北與河南接壤，東與巴水毗鄰，西與倒水交界，南臨長(zhǎng)江，是鄂東北地區(qū)跨行政區(qū)域的一條較大水系。位于東經(jīng)114°40'～115°28'，北緯30°52'～31°36'。全長(zhǎng)170.4 km，流域面積4367.6 km2，上游地區(qū)是著名的大別山區(qū)，重巒疊嶂，地勢(shì)陡峭。中游地區(qū)主要是丘陵，自新洲城關(guān)以下進(jìn)入鄂東濱江平原，地勢(shì)平坦。全流域中山地約43.5%，丘陵占46.2%，平原10.3%。流域地勢(shì)坡向與陸水相反，北高南低，河水自大別山南麓，順著地勢(shì)，直瀉長(zhǎng)江。舉水兩岸支流眾多，流程短促，集水面積不大。河長(zhǎng)超過2.5 km 的支流共有136條。舉水流域年平均降水量1 224 mm，比全省平均值高出近百毫米，且降水季節(jié)分配不均并多暴雨，年際變化大，如上游麻城全年70%的降水量集中于4-8月。最大年降水量超過2 000 mm，最小年降水量不足700 mm。

1.2 研究數(shù)據(jù)來源

DEM 數(shù)據(jù)和坡度高程數(shù)據(jù)來源于中科院地理空間數(shù)據(jù)云，使用ArcGIS 進(jìn)行剪裁變換后得到所需區(qū)域的DEM 圖（圖1）。土地利用數(shù)據(jù)來源于中國(guó)科學(xué)院資源環(huán)境科學(xué)與數(shù)據(jù)中心的Landsat 數(shù)據(jù)，按照國(guó)家土地利用分類方法解譯得到研究區(qū)域土地利用類型圖（圖2）并與SWAT 中的代碼建立土地利用分類檢索表（表1）。土壤數(shù)據(jù)來源于世界土壤數(shù)據(jù)庫(kù)（Harmonized World Soil Database，HWSD）和1∶100 萬的中國(guó)土壤數(shù)據(jù)庫(kù)，經(jīng)GIS重分類處理后得研究區(qū)域土壤類型分布圖（圖3），研究區(qū)域內(nèi)的主要土壤類型為不飽和疏松巖性土，石灰性沖擊土，飽和沖擊土，人為堆積土，飽和潛育土，簡(jiǎn)育高活性淋溶土，不飽和雛形土，飽和黏磐土和水體。氣象數(shù)據(jù)來源于世界氣象數(shù)據(jù)庫(kù)和麻城氣象站的雨量資料（1980-2017年）。水文數(shù)據(jù)來源于柳子港水文站的月徑流實(shí)測(cè)資料（1953-2006年）。

圖1 麻城舉水河流域數(shù)字高程圖Fig.1 Digitai elevation map of Jushui River basin in Macheng City

圖2 麻城舉水河流域土地利用圖Fig.2 Land use map of Jushui River basin in Macheng City

表1 麻城舉水河流域土地重分類檢索表Tab.1 Land use classification index for the Jushui River basin，Macheng City

圖3 麻城舉水河流域土壤類型分布圖Fig.3 Soil map of Jushui River basin in Macheng City

結(jié)合以上數(shù)據(jù)，使用SWAT模塊中的自動(dòng)水系提取功能加載研究區(qū)域的DEM 圖，生成河網(wǎng)并以研究流域出口點(diǎn)作為界定進(jìn)行水文流域的劃分，將研究區(qū)劃分為39 個(gè)子流域，結(jié)果如圖4所示。隨后進(jìn)行HRU 分析與定義，最終輸入氣象數(shù)據(jù)進(jìn)行模擬。

圖4 麻城舉水河流域子流域劃分圖Fig.4 Subarea division map of Jushui River basin in Macheng City

2 模型校準(zhǔn)與驗(yàn)證

2.1 模型校準(zhǔn)方法

模型的參數(shù)優(yōu)化方法選擇SWAT-CUP 自帶的四種參數(shù)優(yōu)化方法，即SUFI-2 算法、GLUE 算法、PSO 算法、Parasol 算法。模擬結(jié)果好壞程度可以通過一系列的統(tǒng)計(jì)指標(biāo)來表示，本文采用決定系數(shù)（R2）和效率系數(shù)（NSE）來評(píng)價(jià)模型的適用性。

SUFI-2 算法是利用拉丁超立方隨機(jī)采樣方法來進(jìn)行水文模型不確定分析的一種常用優(yōu)化方法。首先確定目標(biāo)函數(shù)和假設(shè)初始參數(shù)的范圍，然后進(jìn)行參數(shù)校準(zhǔn)，使得其落在95%的置信區(qū)間內(nèi)（即為95PPU，意為不考慮極壞的5%的情況），再通過拉丁超立方方法（LH-OAT）進(jìn)行參數(shù)結(jié)果的輸出，能得到多種參數(shù)結(jié)果，對(duì)每一種參數(shù)結(jié)果進(jìn)行分析，不斷模擬從而減少不確定的范圍。SUFI-2 算法結(jié)果即參數(shù)的不確定性程度是用不確定性程度因子P-factor（范圍是0～1）和校準(zhǔn)效果因子R-factor（范圍是0～∞）表示。當(dāng)P-factor 越來越接近1和R-factor 越來越接近0 時(shí)說明模擬結(jié)果與實(shí)測(cè)結(jié)果越吻合，模擬效果越好。SUFI-2 算法中的決定系數(shù)R2和效率系數(shù)NSE計(jì)算公式如下：

式中：QS,i為第i時(shí)段實(shí)測(cè)流量值；為實(shí)測(cè)徑流的平均值；QM,i為第i時(shí)段模擬流量值；為模擬流量的平均值；n'表示實(shí)測(cè)時(shí)間序列長(zhǎng)度；n表示樣本個(gè)數(shù)。

GLUE 算法（Generalized Likelihood Uncertainty Estimation）即普適的似然不確定性分析方法。GLUE 算法主要是進(jìn)行參數(shù)識(shí)別，避免在不確定條件下出現(xiàn)不同參數(shù)產(chǎn)生同樣結(jié)果的情況，解決最優(yōu)參數(shù)不唯一問題。GLUE 算法的核心在于似然函數(shù)的選擇，但選擇會(huì)受到人為主觀性地影響。該方法通過Monte-Carlo 隨機(jī)抽樣方法獲得多組參數(shù)，然后計(jì)算每組參數(shù)的模擬結(jié)果與測(cè)量值之間的似然函數(shù)以及權(quán)重，在似然函數(shù)值中選擇一個(gè)分界值，將所有參數(shù)組歸納為兩類，低于似然函數(shù)值的類別說明其不滿足模擬結(jié)果，不足以表征模型的特點(diǎn)。高于似然函數(shù)值的類別則與之相反，選擇高于似然函數(shù)值的類別，并將它們做歸一處理?？芍Ｐ驮谥眯艆^(qū)間的不確定性范圍與似然函數(shù)值的大小有關(guān)。該數(shù)學(xué)表達(dá)式如下：

式中：L(Y|θa)表示后驗(yàn)似然值；L(Qa|Y)表示觀測(cè)到的變量值；L0θa表示先驗(yàn)似然值；C表示歸一加權(quán)因子。

PSO 算法（Particle Swarm Optimization）即粒子群優(yōu)化算法。該方法指將系統(tǒng)初始化為一組隨機(jī)的粒子（有一個(gè)矢量速度決定它們飛行地路線）然后通過迭代的方法搜尋最優(yōu)值。粒子需要通過粒子本身和整個(gè)種群的移動(dòng)所找到的最優(yōu)解來調(diào)整自己下一步的搜尋方向。在SWAT模型中PSO算法的應(yīng)用是將需要校準(zhǔn)的每組參數(shù)作為每個(gè)粒子的位置坐標(biāo)進(jìn)行迭代運(yùn)算，最終根據(jù)模型目標(biāo)函數(shù)進(jìn)行適應(yīng)性評(píng)價(jià)。

Parasol 算法的核心內(nèi)容是SCE-UA 算法，是一種全局性的優(yōu)化算法。Parasol 算法能成功解決非線性的復(fù)雜分布式水文模型問題。該法結(jié)合了空間復(fù)合點(diǎn)、系統(tǒng)演化、競(jìng)爭(zhēng)演化和混合復(fù)合形狀四個(gè)概念。算法的第一步是在模型需要率定的可行空間內(nèi)將隨機(jī)產(chǎn)生的a×b個(gè)點(diǎn)作為初始群體。第二步是按照目標(biāo)函數(shù)增長(zhǎng)的序列將初始群體分成a個(gè)種群，則每個(gè)種群包含了b個(gè)成員。第三步是進(jìn)行種群之間的若干代獨(dú)立競(jìng)爭(zhēng)進(jìn)化和定期交叉作用產(chǎn)生新的種群，這種情況下可以實(shí)現(xiàn)信息被種群共享。最后檢查是否滿足收斂的要求，不滿足則重新進(jìn)行第二步，算法特點(diǎn)是能夠很好地表達(dá)出不同參數(shù)之間的相關(guān)性。以下是Parasol 算法中的目標(biāo)函數(shù)式和NSE 效率系數(shù)計(jì)算式：

式中：SSQ表示殘差平方和；表示實(shí)測(cè)流量平均值；QS,i表示第i時(shí)段實(shí)測(cè)流量值；QM,i表示第i時(shí)段模擬流量值；n表示實(shí)測(cè)時(shí)間序列長(zhǎng)度。

2.2 參數(shù)敏感性分析

參數(shù)敏感性分析一般采用SWAT自帶的LH-OAT方法或者SWAT-CUP 程序中的LH-OAT 方法，本文選擇了第二種方法進(jìn)行敏感性分析。參數(shù)敏感性程度的識(shí)別則可以采用t檢驗(yàn)法，即當(dāng)t值越大，p值越小時(shí)，參數(shù)的敏感性越高。在此基礎(chǔ)上最終確定了對(duì)徑流參數(shù)敏感影響程度最高的4 個(gè)參數(shù)（表2）。

表2 舉水河流域參數(shù)敏感性排名Tab.2 Sensitivity ranking of parameters in Jushui River basin

3 模擬結(jié)果分析

決定系數(shù)R2和效率系數(shù)NSE被采用作為本研究區(qū)域模型的適應(yīng)性評(píng)價(jià)指標(biāo)。R2表示實(shí)測(cè)值與觀測(cè)值的擬合程度，R2越接近1，表明模型越吻合；NSE能判斷模型模擬程度的好壞，當(dāng)NSE≥0.75 時(shí)，表明模型的率定效果很好。一般認(rèn)為，在月尺度下同時(shí)滿足R2＞0.6，NSE＞0.5 時(shí)[22]，說明構(gòu)建的SWAT模型適用于本流域的水文模擬。

本研究構(gòu)建舉水河流域的月尺度徑流模型。選取流域出口的柳子巷水文站1995-2004年的水文實(shí)測(cè)資料作為依據(jù)，其中1995-1996年為預(yù)熱期，1997-2000年為率定期，2001-2004年為驗(yàn)證期，應(yīng)用SUFI-2、GLUE、PSO、Parasol 四種算法進(jìn)行徑流模擬。研究區(qū)域SWAT 模型模擬結(jié)果評(píng)價(jià)見表3，月徑流尺度下實(shí)測(cè)徑流與模擬徑流對(duì)比分析結(jié)果見圖5～圖8。由表1可知，初次使用4種方法進(jìn)行校準(zhǔn)時(shí)，SUFI-2算法在率定期與驗(yàn)證期R2和NSE分別為0.73、0.66 和0.62、0.59，是唯一滿足精度要求的，GLUE 僅在驗(yàn)證期內(nèi)滿足要求。表明SWAT模型在麻城舉水河流域的建立具有一定的適應(yīng)性，可以進(jìn)行相應(yīng)的模擬與預(yù)測(cè)。模擬結(jié)果見圖5，在整個(gè)模擬階段內(nèi)率定期模擬效果比驗(yàn)證期好，率定期內(nèi)平水期模擬值偏低，但枯水期與豐水期實(shí)測(cè)值與模擬值的整體徑流趨勢(shì)線基本一致，徑流峰值誤差小；驗(yàn)證期年徑流量大，枯水期與平水期階段模擬效果較好，豐水期（4-8月）階段實(shí)測(cè)值遠(yuǎn)低于模擬值，月徑流峰值誤差大，模型模擬效果較差。GLUE、PSO、Parasol算法結(jié)果如表3所示，率定期與驗(yàn)證期的R2和NSE沒有同時(shí)滿足R2＞0.6，NSE＞0.5，不符合精度要求。對(duì)比圖6～圖8可知，GLUE算法與Parasol算法結(jié)果相近，兩種算法在枯水期階段模擬徑流曲線變化趨勢(shì)基本吻合，但在豐水期和平水期階段GLUE 模擬結(jié)果與實(shí)測(cè)值之間誤差更小，模擬效果更佳。PSO 算法模擬效果是最差的，在整個(gè)模擬時(shí)段內(nèi)模擬值都偏低，尤其在豐水期，不能體現(xiàn)本流域的徑流變化。

表3 舉水河流域SWAT模型模擬結(jié)果評(píng)價(jià)Tab.3 Evaluation of SWAT Model Simulation Results in Jushui River basin

圖5 SUFI-2算法的舉水河流域?qū)崪y(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.5 Comparison and analysis of measured runoff and simulated runoff in Jushui River basin based on Sufi-2 algorithm

圖6 GLUE算法的舉水河流域?qū)崪y(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.6 Comparison and analysis of measured runoff and simulated runoff in Jushui River basin based on GLUE algorithm

圖7 PSO算法的舉水河流域?qū)崪y(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.7 Comparison and analysis of measured runoff and simulated runoff in Jushui River basin based on PSO algorithm

圖8 Parasol算法的舉水河流域?qū)崪y(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.8 Comparison and analysis of measured runoff and simulated runoff in Jushui River basin based on Parasol algorithm

在SUFI-2 算法中不斷調(diào)整置信區(qū)間直到達(dá)到最佳擬合效果，此時(shí)的參數(shù)置信區(qū)間為SUFI-2 最適置信區(qū)間見表4，將所得的參數(shù)置信區(qū)間作為已知信息，采用GLUE、PSO、Parasol 算法進(jìn)行第二次參數(shù)優(yōu)化，得到SWAT 模型模擬結(jié)果見表5，GLUE、PSO、Parasol 算法二次參數(shù)優(yōu)化后的實(shí)測(cè)徑流與模擬徑流的對(duì)比分析結(jié)果見圖9～圖11。

表4 基于SUFI-2算法的參數(shù)適宜置信區(qū)間Tab.4 Suitable confidence interval of parameters based on Sufi-2 algorithm

表5 舉水河流域SWAT模型二次模擬結(jié)果評(píng)價(jià)Tab.5 Evaluation of secondary simulation results of SWAT Model in Jushui River basin

圖9 GLUE算法二次優(yōu)化后實(shí)測(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.9 Comparison and analysis of measured runoff and simulated runoff by secondary simulation of GLUE algorithm

圖10 PSO算法二次優(yōu)化后實(shí)測(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.10 Comparison and analysis of measured runoff and simulated runoff by secondary simulation of Pso algorithm

圖11 Parasol算法二次優(yōu)化后實(shí)測(cè)徑流與模擬徑流結(jié)果對(duì)比分析圖Fig.11 Comparison and analysis of measured runoff and simulated runoff by secondary simulation of Parasol algorithm

由表5 可知，3 種算法的模擬精度均有明顯提高，其中GLUE 和Parasol 算法二次優(yōu)化后率定期與驗(yàn)證期的結(jié)果滿足R2＞0.7，NSE＞0.7，精度均高于SUFI-2 算法，模擬效果較好。PSO算法所得評(píng)價(jià)因子相較之前也有所提高，但在驗(yàn)證期內(nèi)仍不滿足精度要求，在本流域是不適用的。對(duì)比圖5 與圖9 可知，GLUE 算法優(yōu)化后相較于SUFI-2 算法，在枯水期與豐水期的模擬效果差別不大，但汛期模擬值與實(shí)測(cè)值之間誤差減小，月徑流模擬值與實(shí)測(cè)徑流趨勢(shì)更加吻合，徑流峰值也更相近，對(duì)于本流域參數(shù)優(yōu)化工作GLUE更合適。Parasol算法與GLUE 算法精度結(jié)果相近。Parasol算法在枯水期和豐水期的模擬效果比GLUE算法差，但在春汛出現(xiàn)單峰徑流的情況時(shí)模擬值是最接近實(shí)測(cè)值的，綜合來看本研究區(qū)域最佳優(yōu)化方法為GLUE算法或Parasol算法。

4 討論

由模擬評(píng)價(jià)結(jié)果可知，率定期的模擬效果均要高于驗(yàn)證期，這可能是由于驗(yàn)證期的年徑流量大，年際徑流變化大造成的。對(duì)比參數(shù)方法優(yōu)化前后算法的精度結(jié)果可知，三種算法的模擬精度都明顯提高，說明SUFI-2 算法所得的參數(shù)最適置信區(qū)間對(duì)算法徑流模擬結(jié)果具有很大影響。由以上模擬結(jié)果圖可知參數(shù)最適置信區(qū)間對(duì)模擬階段內(nèi)的枯水期徑流影響變化并不大，主要是針對(duì)汛期。這可能是由于汛期受到強(qiáng)降水和土壤植被覆蓋物、土壤類型的影響。汛期徑流量大，土壤層之間水量傳遞的變化特征顯著。4 種算法也是各有優(yōu)缺點(diǎn)，SUFI-2 算法耗時(shí)短在枯水期效果好但對(duì)于較為復(fù)雜的流域變化情況效果差；GLUE 算法在枯水期和平水期模擬效果好，常用于出現(xiàn)春汛和干旱缺水年份時(shí)；PSO算法精度低，不適用于干旱缺水地區(qū)的模擬；Parasol 算法在春汛出現(xiàn)單峰徑流時(shí)模擬效果最佳。

5 結(jié)論

在麻城舉水河流域建立SWAT模型對(duì)柳子港水文站的月徑流過程進(jìn)行模擬，通過分析結(jié)論如下：

（1）通過SWAT-CUP 中LH-OAT 進(jìn)行參數(shù)敏感性分析，并采用SUFI-2、GLUE、PSO、Parasol四種參數(shù)優(yōu)化方法進(jìn)行校準(zhǔn)?？芍狙芯繀^(qū)域敏感性影響最大的4個(gè)參數(shù)因子分別為CN2、ALPHA_BF、GW_DELAY、GWQMN。SUFI-2 算法在率定期與驗(yàn)證期的評(píng)價(jià)因子R2、NSE滿足精度要求，表明SWAT模型在本流域可以適用，但模擬精度不高。

（2）SUFI-2 算法不斷調(diào)參優(yōu)化后得到一個(gè)滿足最高精度要求的參數(shù)適宜置信區(qū)間，對(duì)比參數(shù)置信區(qū)間影響前后的3種優(yōu)化算法（GLUE、PSO、Parasol）徑流結(jié)果，GLUE、PSO、Parasol 算法的精度均有明顯提高，但PSO 算法仍不滿足精度要求，在本流域是不適用的。GLUE、Parasol算法滿足R2＞0.7，NSE＞0.7，精度均高于SUFI-2算法。本流域選擇GLUE 算法或Parasol算法更佳。

（3）根據(jù)徑流曲線變化選擇適宜的優(yōu)化參數(shù)方法前，先采用SUFI-2 方法獲取適宜參數(shù)置信區(qū)間再進(jìn)行校準(zhǔn)，有利于提高模擬精度，選取最佳校準(zhǔn)方法，可減少模擬過程中優(yōu)化方案選擇的任意性。