張響亮 許文華 楊明燦 王萬宣 梁寶元

(1.仙游縣氣象局，福建 莆田 351200; 2.莆田市氣象局，福建 莆田 351100)

準確預(yù)測內(nèi)澇點積水深度有利于防洪力量的合理調(diào)配，提高防汛救災(zāi)部門的工作效率，降低內(nèi)澇造成的人員傷亡和經(jīng)濟損失[1]。隨著全球氣候變化，強降水及其不同尺度極端事件頻發(fā)[2]，城市內(nèi)澇出現(xiàn)的頻率及危害將不斷增大[3]，面對嚴峻的防洪防澇挑戰(zhàn)，內(nèi)澇點積水深度預(yù)測模型研究成為國內(nèi)外熱點。城市水文模型經(jīng)歷了集總式水文模型向更具物理意義的分布式水文模型發(fā)展[4]，國內(nèi)外常用的內(nèi)澇模型有美國環(huán)境保護署的SWMM模型(暴雨洪水管理模型)、法國電力的TELEMAC-MASCARET(水動力學和水文學領(lǐng)域的高性能數(shù)值仿真開源軟件)等。國內(nèi)學者在SWMM模型的基礎(chǔ)上開發(fā)出了許多適用于我國的內(nèi)澇模型。隨著人工智能技術(shù)的發(fā)展，應(yīng)用機器學習方法建立的回歸模型因其預(yù)測效果較水文模型并不遜色，又具有計算快捷的優(yōu)勢，具有較好的實用性[5]。研究表明，降水是內(nèi)澇災(zāi)害最主要的誘發(fā)因素[6]。此外城市內(nèi)澇還受到地表匯流、管網(wǎng)排水和河道排澇等水文動力過程影響[7]。

綜上所述，水文模型受到地形數(shù)據(jù)、降水數(shù)據(jù)、運算能力的影響，存在較大局限性。本文基于內(nèi)澇點積水深度資料和臨近氣象觀測站小時降水量資料，設(shè)計滯后系數(shù)與堵塞系數(shù)的計算方法，滯后系數(shù)可表示地表匯流對積水深度的影響，堵塞系數(shù)可表示排水設(shè)施的運行情況。將滯后系數(shù)、堵塞系數(shù)、小時雨量和前一時次積水深度作為輸入層，采用機器學習中的嶺回歸方法建立了內(nèi)澇點積水深度預(yù)報方程。預(yù)測檢驗顯示了其較好的實用性。

1 資料與方法

1.1 資料選取

雨量數(shù)據(jù)選用莆田市國家氣象觀測站(以下簡稱莆田站)的逐小時降水數(shù)據(jù)，積水深度數(shù)據(jù)為莆田市荔園路才子路口(以下簡稱才子路口)的積水深度數(shù)據(jù)，莆田站與才子路口直線距離約3km。

選用2018年6月6日至8月31日莆田站逐小時降水量和才子路口整點積水深度建立回歸方程?；貧w方程檢驗資料選用2018年9月1日至2018年9月30日莆田站逐小時降水量和才子路口整點積水深度。

1.2 研究方法

1.2.1 回歸方程

通過特征值與目標值進行擬合，可建立回歸模型從而實現(xiàn)利用特征值對目標值進行預(yù)測。機器學習在進行回歸模型預(yù)測時，預(yù)測值和實際值之間存在一定誤差。找到最優(yōu)的估計參數(shù)使誤差最小，即可建立最優(yōu)的回歸模型，但在數(shù)據(jù)量較少或者數(shù)據(jù)有共線性時，會導致權(quán)重趨于無窮大，需要引入一個懲罰項來改善這一問題[5]。

1.2.2 嶺回歸方程

嶺回歸算法是基于最小二乘法的有偏估計回歸方法，加上了L2正則化表達式，考慮了各特征值的重要程度，能夠適應(yīng)共線性問題。雖然該方法損失了無偏性，但是能夠?qū)簿€性數(shù)據(jù)有較好的擬合效果。本文采用嶺回歸算法建立內(nèi)澇點積水深度預(yù)報方程。

1.2.3 積水深度預(yù)報方程建立

基于2018年6月6日至8月31日莆田站逐小時降水數(shù)據(jù)和才子路口積水深度數(shù)據(jù)，計算滯后系數(shù)與堵塞系數(shù)。將滯后系數(shù)、堵塞系數(shù)、當前小時雨量和前一時次積水深度作為輸入層，采用嶺回歸方法，得到才子路口積水深度預(yù)報方程。

2 數(shù)據(jù)優(yōu)化

滯后系數(shù)與堵塞系數(shù)的計算需要連續(xù)性的雨量數(shù)據(jù)和積水深度數(shù)據(jù)，部分數(shù)據(jù)存在以下明顯錯誤：①降水過程停止后，低處的積水無法正常排出，導致內(nèi)澇點積水深度始終維持在較低的數(shù)值；②在小尺度的降水過程中，無法準確獲取內(nèi)澇點的降水情況，導致降水與積水深度變化不一致。為此，在建立回歸方程前，需要對內(nèi)澇點積水深度數(shù)據(jù)進行修訂。

2.1 低值積水數(shù)據(jù)處理

積水深度下降至低洼處時，積水無法匯入排水設(shè)施，排水速度異常緩慢。使用該數(shù)據(jù)建立的回歸方程整體排水速度偏慢，導致預(yù)測的積水深度偏高。因此需要判斷低洼處積水深度的臨界值，當積水深度低于臨界值時，計算出正常排水速度下的積水深度。

2.1.1 積水深度臨界值計算

原始數(shù)據(jù)中的積水深度下降至臨界值后，積水深度下降緩慢，導致積水深度低于臨界值的數(shù)據(jù)量較大。故認為積水深度低于某一值的數(shù)據(jù)量超過設(shè)定的閾值時，該深度為臨界值。本文設(shè)置的閾值為95%，計算的積水深度臨界值為2cm。

2.1.2 正常排水速度計算

積水深度大于臨界值時，積水可以順利匯入排水設(shè)施。篩選出積水深度大于臨界值、臨近2個時次無降水、積水深度比上一時次低的數(shù)據(jù)，計算積水深度降低的平均速度，作為正常排水速度，本文計算的正常排水速度為1.03cm/h。

2.2 積水深度異常上升

原始數(shù)據(jù)中，部分數(shù)據(jù)出現(xiàn)“降水停止后的積水深度依舊上升”的情況。使用該數(shù)據(jù)建立的回歸方程，降雨量的權(quán)重明顯偏低，甚至導致降雨量與積水深度呈反比。故使用正常排水速度計算出的積水深度替換臨近兩個時次無降水、且積水深度大于前一時次積水深度的數(shù)據(jù)。

3 特征值選取

內(nèi)澇點的積水深度主要受降水、前期累計降水、地表匯流、管網(wǎng)排水等因素影響。本文選用小時雨量、前一時次積水深度、滯后系數(shù)與堵塞系數(shù)作為影響積水深度的預(yù)報因子。以下介紹滯后系數(shù)和堵塞系數(shù)的計算方法。

3.1 滯后系數(shù)

姜曉岑等[8]指出，積水一般滯后降水約5～30min出現(xiàn)。由于小時雨量數(shù)據(jù)無法準確描述降水的時間分布情況，若較強的降水集中出現(xiàn)在中后時段，當前時次的積水深度沒有明顯變化，而是影響下一時次的積水深度。因此降水導致的積水深度若沒有達到預(yù)期的積水深度，則認為部分降雨量將滯后影響下一時次的積水深度。實際深度與預(yù)期深度相差越大，滯后系數(shù)也越大。積水深度預(yù)期值與滯后系數(shù)計算方法如下。

3.1.1 積水深度預(yù)期值

積水深度預(yù)期值計算公式如下：

hi=hi-1+(Δh×ri-v)

式中，hi為當前時次積水深度預(yù)期值，hi-1為前1小時的積水深度，Δh為單位雨量產(chǎn)生的積水深度，ri為當前時次降水量，v為正常排水速度。

計算單位降水量產(chǎn)生的積水深度，需要篩選降水超過排水設(shè)施承受臨界值的數(shù)據(jù)。由于降雨量超過臨界值的比重較大，在產(chǎn)生積水的數(shù)據(jù)中，分別計算降雨量大于50mm、40mm、30mm、20mm、10mm時產(chǎn)生積水的數(shù)據(jù)比重，當比重超過設(shè)置的閾值時，認為對應(yīng)的雨量為產(chǎn)生積水的雨量臨界值，本文閾值設(shè)置為95%，計算出的臨界值為10mm。由于產(chǎn)生積水的過程中，排水設(shè)施使積水深度下降，因此雨量產(chǎn)生的積水為積水深度實際變化值與排水速度之和，單位雨量對積水深度增長值的計算公式如下：

式中，Δh為單位雨量對積水深度增長值，hi為當前時次積水深度，hi-1為上一時次積水深度，v為積水下降速度(積水下降時為正值)，ri為當前時次降雨量，m為篩選后的樣本總數(shù)，本文計算的單位雨量對積水深度增長值為0.04cm/mm。

3.1.2 滯后系數(shù)計算

將上一時次積水深度的預(yù)期值與積水深度實際值的差值，作為當前時次的滯后系數(shù)，若滯后系數(shù)大于0，則考慮滯后降水影響當前時次的積水深度，若滯后系數(shù)小于0，則不考慮滯后降水對積水深度的影響，滯后系數(shù)修訂為0。

3.2 堵塞系數(shù)

排水設(shè)施在使用過程中存在被雜物堵塞等情況，導致排水速度降低?？紤]堵塞系數(shù)，可以解釋某些雨量較為均勻的降水過程中，積水深度變化不規(guī)律的情況，該系數(shù)亦可用于判斷各內(nèi)澇點排水設(shè)施是否存在堵塞。堵塞系數(shù)計算公式如下：

bi=(hi-1-hi-2)-(ri-1×Δh-v)

式中，bi為當前時次堵塞系數(shù)，hi-1、hi-2為前1小時和前2小時的積水深度，ri-1為當前時次的降雨量，Δh為單位雨量對積水深度的貢獻，v為積水深度下降速度。

計算得到的堵塞系數(shù)越大，則說明堵塞情況嚴重。例如降水停止后，計算的(ri-1×Δh-v)為負的定值，但若堵塞情況嚴重時，(hi-1-hi-2)應(yīng)為趨近于0的負值，計算的b值大。若堵塞情況較輕，則(hi-1-hi-2)項數(shù)值越小，計算的b值越小。

4 模型訓練及檢驗

4.1 算法與參數(shù)選取

scikit-learn庫(針對Python 編程語言的機器學習庫)中含有多種機器學習算法[9]，本文選擇scikit-learn庫中嶺回歸算法對內(nèi)澇方程進行訓練[10]。嶺回歸算法是基于最小二乘法的有偏估計回歸方法，加上了L2正則化表達式，考慮了各特征值的重要程度，能夠適應(yīng)共線性的問題。模型輸入變量為2018年6月至8月莆田站的小時雨量數(shù)據(jù)、才子路口的積水深度數(shù)據(jù)以及計算的滯后系數(shù)和堵塞系數(shù)。參數(shù)solver選擇“auto”(可以根據(jù)樣本量自動選擇回歸方程的優(yōu)化方案)、normalize選擇“False”(經(jīng)檢驗，關(guān)閉樣本特征歸一化后計算的均方跟誤差更小)，訓練出預(yù)報方程為：

yi=0.003ri+0.96hi-1+0.286li+0.264bi-0.208

式中，yi為預(yù)測的積水深度，ri當前降水量，hi-1為前一時次積水深度，li為滯后系數(shù)，bi為堵塞系數(shù)。

4.2 模型檢驗

通過9月1日00時至9月30日23時的積水深度數(shù)據(jù)對預(yù)報方程進行檢驗，預(yù)報結(jié)果與實際值的均方誤差為0.27cm，說明預(yù)報效果較好，但由于9月強降水的出現(xiàn)次數(shù)較少，積水深度為0的數(shù)據(jù)比例較高，因此檢驗結(jié)果無法代表強降水出現(xiàn)時積水深度的預(yù)報效果。故篩選9月出現(xiàn)降水過程的數(shù)據(jù)再次進行檢驗(9月降水數(shù)據(jù)中，僅 9月7日14時至9日07時、9月23日08時至24日18時出現(xiàn)較強降水)，降水過程的預(yù)報結(jié)果與實際值的均方誤差為1.09cm，總體的預(yù)報偏差在2cm以內(nèi)，能一定程度反映強降水過程的積水情況，但測試數(shù)據(jù)中出現(xiàn)雨量很小但積水深度突增的情況，可能是由于雨量站與積水點位置偏差導致雨量數(shù)據(jù)無法反映內(nèi)澇點的降水情況，因此對于無降水時積水深度上漲的時次無法準確預(yù)測。降水過程檢驗的部分結(jié)果見表1。

表1 9月降水過程降雨量、積水深度實際值與積水深度預(yù)測值

5 結(jié)論與討論

本文利用2018年6月至8月莆田站逐時降水資料、才子路口整點積水深度資料，應(yīng)用機器學習中的嶺回歸方法建立了內(nèi)澇點積水深度預(yù)報方程，并使用2018年9月1日00時至30日23時的降水過程進行檢驗，回歸方程的均方誤差為0.27cm，降水過程檢驗結(jié)果的均方誤差為1.09cm，預(yù)報效果較好。

回歸方程的效果依賴于實況數(shù)據(jù)的準確性，若使用更好的初始數(shù)據(jù)，并且加強內(nèi)澇點積水深度觀測站設(shè)備的維護和保障，建立的積水深度預(yù)報方程質(zhì)量更高。此外，本文的滯后系數(shù)和堵塞系數(shù)在其數(shù)值大小上能反映滯后和堵塞的情況，但其系數(shù)的具體算法還有待研究。