廣義相加模型在烏江夏季徑流預報中的應用

榮艷淑，胡玉恒，馮瑞瑞，殷雨婷，李崇浩

(1.河海大學水文水資源學院，江蘇 南京 210098；2.河海大學水利學科專業實驗教學中心，江蘇 南京 210098；3. 南京恩瑞特實業有限公司，江蘇 南京 211106；4.中國南方電網電力調度控制中心，廣東 廣州 510623)

中長期徑流預測一直是水文預報的一大難題。傳統的統計水文預報模型，如時間序列分析和各種回歸模型等，由于計算簡便，模型參數有一定的物理意義而被廣泛地應用于水文預報工作[1]。如，基于多種氣象因子建立的逐步回歸模型、多元線性回歸模型和時間序列分解模型對月、旬徑流都有較好的預報效果[2-4]。但是傳統線性統計模型沒有考慮地球氣候系統中各要素與徑流之間復雜的非線性關系，預測精度很難有進一步提高[5]。隨著計算機技術的發展，人工神經網絡、改進人工神經網絡和小波分解等，不僅能夠處理變量之間的非線性映射關系[6]，而且預測結果明顯優于線性回歸模型，預報精度也有明顯提高[7-9]。但是這些方法存在模型物理意義不明確，模型通用性較差等缺陷[10]。近年來涌現了一些考慮了非線性關系，或者使用多種模型融合的技術方法，明顯提高了預測水平[11-15]。廣義相加模型(generalized additive model,GAM)就是綜合了線性關系和非線性關系的一種復合模型，它承認預測變量與某些物理因子存在線性關系，與另外的物理因子可能存在非線性關系，并基于線性和非線性因子建立綜合模型。這一模型在基于海溫因子的降水預測、基于氣溫和降水的大氣污染影響分析和基于NAO的阿爾卑斯山冬季極端氣溫的預測中都得到較好的應用[16-18]。

烏江流域位于長江上游，流域內水電站網密布，為當地的農業灌溉、生活用水以及防洪發電帶來了顯著效益。但是，流域內徑流預測準確性一直困擾著水電站網管理部門，特別是汛期季風盛行期間，降水變化顯著，徑流變化更難以預測。在這一背景下，基于廣義相加模型，筆者對烏江流域的洪家渡夏季徑流進行預測研究，并對GAM的建立、模擬和預測進行評估和探討，以期能為烏江流域的徑流預測提供一個新途徑。

1 研究數據與方法

1.1 研究數據

以烏江流域的洪家渡夏季徑流為研究對象，以多種海溫指數作為預測因子，包括太平洋、大西洋和印度洋的26個海溫指數等(表1)。徑流資料來自洪家渡水電站，海溫指數來自國家氣候中心。由于海溫與我國夏季降水之間存在較長的滯后相關性[19]，因此選取前期冬季海溫指數進行研究。冬季指數為12月、翌年1—2月3個月的平均值，夏季徑流為6—8月3個月的平均值。資料時間范圍均為1952—2014年。

表1 海溫指數與因子代碼

1.2 GAM模型

GAM是由Hastie等[20]在廣義線性模型(GLM)[21]的基礎上改進而來，該模型包括兩部分，分別描述了自變量與因變量的線性和非線性關系，因此GAM存在參數和非參數兩部分，屬于半參數模型。這一改進打破了線性模型中關于自變量與因變量線性相關的假定，使得模型應用有了進一步拓展。GAM可表示為

(1)

式中：g(yi)——因變量yi的連接函數；l——序列長度；β0——截距；s(xki)——平滑函數，本文采用薄板樣條函數作為平滑函數；n——非線性因子的個數；xki——非線性因子；βj——線性部分的回歸系數；m——線性因子的個數；xji——線性因子；εi——隨機誤差，服從高斯正態分布。

取l= 63，為避免平滑過程中出現過度擬合，將薄板樣條函數節點的最大個數固定為5。GAM參數部分檢驗使用t檢驗，非參數部分使用F檢驗。

1.3 模型模擬和預測的評估方法

1.3.1 模型模擬的評估方法

基于式(1)對洪家渡夏季徑流序列進行模擬，模擬效果選用最小信息準則(Akaike information criterion, AIC)、均方根誤差(root mean square error, RMSE)、平均絕對誤差(mean absolute error, MAE)、概率空間線性誤差(linear error in probability space, LEPS)[22]和線性相關系數(correlation,r2)5種評價指標進行評估。根據最小信息準則，AIC值越小表明模型模擬效果越好；RMSE用來衡量模擬值與觀測值之間的相對偏差，其值越小表明偏差越小；MAE能夠反映模擬值與實測值的絕對偏離程度，其值越小模擬效果越好；LEPS使用了Potts等[23]改進的LEPS，其值越接近于1，模擬效果越好。r2刻畫2個序列之間的線性相關，其值越接近1表明整個序列模擬值與實測值的相關性越高。

1.3.2 模型預測的評估方法

一般的模型都是取序列長度的一部分數據建立模型，剩余數據作為預測樣本，用以評估模型的模擬能力[24]，但是這種方法存在一定的主觀性以及測試樣本長度的局限性，特別是當預測樣本過長時，模型易出現衰退，導致預測效果不好。

交叉驗證是一個可以克服上述問題的模型預測檢驗方法，其主要思想是將數據樣本分成兩部分，一部分作為訓練集進行模型訓練，另一部分作為測試集進行模型預測，但是該方法是將數據循環分組，對每個原數據都進行模擬，這樣通過多次重復試驗能得到較好的模型估計結果[23]。為了充分利用已有的樣本數據得到更可靠的預測結果，本文使用了廣義交叉驗證[25]確定薄板樣條函數平滑自由度的取值，使用留一法交叉驗證[26]進行模型預測結果的評估。

留一法交叉驗證的主要思想是從l個樣本數據中取一個作為測試集，剩下的l-1個作為訓練集，例如，對于某時間序列的第一個樣本，用i= 2，3，…，l個樣本建立模型，獲得模型參數，再對第一個樣本進行模擬；對第二個樣本，則是用i= 1和i= 3，4，…，l個樣本建立模型，再對第二個樣本進行模擬，如此往復，共進行l次訓練，得到l個預測結果。最后，通過模擬序列和實測序列的相關性分析和預報相對誤差頻數統計，對GAM預報效果進行評價。

表2 GAM中各預報因子的線性和非線性屬性

2 結果與討論

2.1 模型因子挑選

本文選擇前期冬季作為海溫影響的關鍵時期，將前期冬季26個海溫指數作為洪家渡夏季徑流預測的潛在預報因子。同時，為了驗證GAM模擬的效果，使用GLM進行對比分析。根據逐步回歸挑選因子的最小信息準則[27](AIC)分別篩選出GAM和GLM模型的預報因子。GAM挑選出13個預報因子，分別為x1、x2、x3、 x4、 x6、 x7、 x10、 x14、 x17、 x20、 x22、 x24、 x26；GLM挑選出14個預報因子，分別為x4、 x5、 x7、 x8、 x9、 x10、x11、 x12、 x13、x14、 x17、 x20、x22、 x24。2個模型共有的預報因子是8個，這是線性模型和非線性模型對自變量影響因變量的物理機制存在差異的原因。

用薄板樣條函數分別對GAM中的13個因子進行平滑處理，平滑自由度(用edf表示)的值為1時，表明自變量與因變量為線性關系，超過1時，表明自變量與因變量為非線性關系，其值越大，表明非線性關系越明顯。經過平滑處理，有9個因子與洪家渡徑流是線性關系，有4個因子與洪家渡徑流是非線性關系(表2)。由表2可見，NINO1+2區海表溫度距平指數(x1)、黑潮海溫指數(x20)、熱帶印度洋全區一致海溫指數(x24)和副熱帶南印度洋偶極子指數(x26)均對洪家渡徑流有非線性影響，而其他指數，包括NINO3區(x2)、NINO4區(x3)、NINO3.4區(x4)、NINO C區(x6)和NINO A區(x7)海表溫度距平指數、熱帶北大西洋海溫指數(x10)、印度洋暖池強度指數(x14)、大西洋多年代振蕩指數(x17)和暖池型ENSO指數(x22)與洪家渡徑流的關系可用線性關系描述。

圖1 x20和x24與洪家渡夏季徑流偏差的相關關系Fig.1 Relations between summer runoff deviation at Hongjiadu and factors, x20 and for x24

以x20和x24為例，圖1給出它們與模擬徑流量偏差的相關關系，圖中實線為平滑函數值，點線之間的區域為信度95%水平的范圍?？梢钥闯?，x20和x24與模擬值偏差的散點基本分布在信度95%水平之間，表明對x20和x24作為非線性因子的判定是合理的。

2.2 模型模擬效果評估

由圖2可以看出GAM和GLM對徑流量模擬的趨勢相似，且與實測序列基本相一致。從一些異常年份看，比如1954年、1957年、1967年、1998年、2001年、2008年、2012年等幾個流量異常偏多的年份，以及2005年之后幾個流量異常偏少的年份，均得到了較為準確的模擬。通過GAM與GLM的對比分析，在1954年、1957年、1998年、2001年、2008年等幾個流量異常偏多的年份，GAM的模擬效果更接近實測值。在2個模型模擬值明顯不一致的年份，如1967年、1985年、1994年，GAM的模擬值與實測值更為接近，模擬效果更好。

圖2 GAM和GLM模擬值與實測值比較Fig.2 Comparison between simulated values of GAM and GLM and observed values

由表3可知， GAM的AIC值為733.04，而GLM的AIC值為752.57，表明GAM的模擬擬合效果比GLM好。GAM的RMSE與MAE值均低于GLM，說明GAM在整體序列上的模擬精度高于GLM。GAM模擬值與實測值序列的相關性(r2)與GLM相比也有很大的提升。綜合分析，GAM模擬效果比GLM模擬效果有了顯著提升。

表3 GAM和GLM評價指標對比

由圖3(a)可以看出，模型殘差的分布大致呈正態分布的規律，而當值大于100以及小于-50時，部分散點并不完全處于1∶1直線上；圖3(b)也顯示出殘差分布出現了略呈偏態的特征，這表明殘差分布并不完全服從正態分布。采用Shapiro等[28]提出的W′正態檢驗方法，計算了統計量W′=0.991，達到了統計顯著性95%信度水平，認為殘差通過了正態性檢驗，服從正態分布。Montgomery[29]指出使用GAM時，樣本低于100時有可能出現類似圖3的現象。

圖3 模型殘差正態概率圖和殘差直方圖Fig.3 Normal probability plot and histogram of model residuals

為進一步驗證模型的適用性，探討了模型殘差散點分布以及實測流量與模擬流量散點的擬合情況。由圖4(a)可知，殘差大部分處于-50～ 50之間，分布基本合理。從圖4(b)可知，其值大致呈現1∶1線性分布的趨勢，流量存在一個異常大值(>500)也得到了較為準確的模擬，表明模型的擬合效果較好。

圖4 殘差散點圖及流量的實測值和模擬值散點分布Fig.4 Residual scatter plot and scatter plot of observed runoff versus predicted runoff

綜上所述，可以認定GAM用于本文的徑流模擬是合適的。

2.3 模型預測效果評估

采用留一法交叉驗證對GAM

和GLM的預測結果進行對比分析，其中樣本數據為1952—2014年。根據留一法交叉驗證的原理，每次均選取其中的62 a作為訓練集，剩下的1 a作為測試集驗證模型的優劣，循環往復，共進行63次模擬與預測。

圖5為實測值與2種模型預測值的對比，其中實測值與GLM預測序列的相關性為0.267，通過了信度95%的顯著性水平；實測值與GAM預測序列的相關系數為0.407，通過了信度99.99%的顯著性檢驗。因此，相較于GLM，GAM的預測能力也有了顯著提升。在一些流量異常年份的預測上，如1957年、1968年、2001年、2008年幾個流量偏多年份，GAM的預測效果相較于GLM均有明顯提升，而在1985年和1998年2個GLM預測值誤差顯著偏大的年份，GAM的預測結果也有了很大的改善。雖然也有部分年份GLM預測更接近實測值，但從整體趨勢及一些流量異常的年份看，GAM的預測能力均高于GLM。

圖5 2種模型的預測值與實測值對比Fig.5 Comparison between prediction values of GAM and GLM and observed values

從圖6可以看到在10%的相對誤差絕對值范圍內，GAM預報正確年份為13 a，與GLM相比少1 a，而在10%～20%以及20%～30%的誤差范圍內GAM預報準確的年份比GLM分別多3 a和5 a。而當相對誤差絕對值大于40%時，GLM預報的年份要明顯多于GAM。因此，在頻數分布上也可以看出，GAM模型的預報效果明顯優于GLM的預報效果。

圖6 GAM和 GLM預測值相對誤差頻數分布Fig.6 Frequency distribution of prediction relative error of GAM and GLM

如果以預測相對誤差絕對值小于等于30%為合格的預測結果，GAM和GLM的預測合格率分別為85.5%和79.0%。表明GAM的預測合格率高于GLM，預測效果更好。

文獻[15]使用了5種模型(時間序列方法、周期均值疊加法、逐步回歸、嶺回歸和主成分回歸)對2004年和2005年分別進行了預測分析，從6—8月平均相對誤差看，這2 a平均相對誤差分別是18.0%和142.0%，而GAM對這兩年預測的相對誤差分別是9.0%和68.8%，預測精度均明顯提高。這也印證了考慮前期海溫因子的非線性作用后，模型模擬及預測能力均有明顯提高。因此，GAM能夠對傳統線性回歸模型的預報效果加以改進。

3 結 語

在前期海溫因子中，挑選了4個非線性因子和9個線性因子，構建的GAM對洪家渡的夏季徑流預測非常有效，通過AIC、RMSE、MAE、LEPS、相關性5種評價指標評估，GAM都優于GLM。

對洪家渡夏季徑流的預測評估結果表明， GAM的預測效果比GLM明顯提高，從相對誤差絕對值的頻數、極值預報誤差和預測合格率等3個方面的分析表明，GAM的誤差范圍更小，極值預測的精度和預測合格率均比較高。

徑流的中長期預測對水庫調度和水資源管理有重要作用，研究高精度的徑流預測模型意義重大。本文的嘗試盡管提高了洪家渡夏季徑流的預測能力，但是，今后通過因子影響機理分析和消除因子之間的共線性等問題，選取更好的預測指標，或許可以進一步提高GAM的模擬精度。