基于馬氏距離判別的丹江口水庫長期徑流分級預報

程忠良，劉 勇，高 成，胡 健

(1. 河海大學水文水資源學院，南京 210098；2. 南京水利科學研究院 水文水資源與水利工程科學國家重點實驗室，南京 210029；3. 中國科學院 陸地水循環及地表過程重點實驗室，北京 100101)

0 引 言

徑流預報是水文學的重要應用領域之一，也是水庫運行調度和防汛、抗旱以及水資源應急調度等的前提[1]。長期徑流預報因其預見期較長、預報精度不高，難以滿足人們對社會經濟生產安排的需求，一直以來都是水文領域研究的難點之一。水文預測方法眾多，通常大致可分為過程驅動和數據驅動兩大類[2]。其中過程驅動法是基于產匯流機制建立的模型，是徑流預測的一個發展方向。但由于徑流受到氣候氣象、下墊面、人類活動等多方面不確定性因素的影響[3]，其形成機理和規律尚未被完全掌握，預報精度不高，使得該方法的應用十分困難。隨著數據獲取能力和計算能力的發展，數據驅動模型在水文預報中的應用也愈加廣泛，神經網絡[4]、支持向量機[5]等預測方法在實際應用中也取得了一定的成果。近年來，通過發掘未來徑流與前期降雨、海表溫度、大氣環流指數等大尺度氣候-水文變量間的關聯關系，建立基于物理成因背景的徑流預報模型是目前研究的重要方向[6-9]。馮小沖[6]、劉勇[7]等通過對影響丹江口水庫入庫徑流的氣候和水文氣象物理因子的分析，篩選出與徑流密切相關的物理因子作為預報因子，分別采用逐步回歸和神經網絡的方法構建月徑流預報模型，取得了一定的成果。在氣候變化和高強度人類活動的背景下，全球海洋和大氣環流均發生了顯著的年代際變化[10]，水庫來水量的變化也隨之加劇，單獨采用定性或定量的回歸預報模型存在預報誤差偏大和精度偏低的問題。為進一步提高預報的精度和可靠性，本文擬在前人研究的基礎上，采用定性定量相結合的預報方法，從影響丹江口水庫長期徑流的物理成因出發選擇預報因子，采用AIC準則篩選關鍵因子，劃分徑流級別，利用馬氏距離判別分析的分類判別優勢，構建長期徑流分級預報模型，并與逐步回歸模型作比較。

1 研究思路及方法

1.1 研究思路

首先從物理成因出發，考慮影響丹江口水庫入庫徑流的海溫、環流等因素，采用單相關系數法，計算預報因子與預報徑流之間的相關系數r，并作0.05信度水平下的檢驗，初選出基礎預報因子集。然后，將資料序列劃分為模擬期和檢驗期，模擬期的徑流按照距平要素值劃分為豐、平、枯三個級別，采用AIC準則進一步篩選出關鍵因子集，構建分段多元線性回歸方程，對預報對象進行模擬；檢驗期時，采用馬氏距離判別法根據預報因子對預報徑流先做出定性判別，再依據定性判別類別代入相應的線性回歸方程，做出定量預報，并將該預報模型與逐步回歸模型的模擬和預報結果進行對比，從預報精度和穩定性兩個方面去比較兩種模型。

1.2 研究方法

1.2.1 單相關系數法

在水文中長期預測中常采用相關系數來考察預報因子與預報對象之間是否線性相關，并以此作為因子挑選的依據[3]。單相關系數的公式為：

(1)

1.2.2 AIC準則

因子個數p的選取對模型的預報精度和穩定性具有很大的影響。取較小的p，擬合程度較差；取較大的p，則受偶然變化影響太大，甚至出現過度擬合的現象，預報效果會受到影響。AIC是衡量統計模型擬合優良性的一種標準，由日本統計學家赤池弘次在1974年提出，它建立在熵的概念上，提供了權衡估計模型復雜度和擬合數據優良性的標準。現定義一個準則函數為：

(2)

式中：σ2(p)為誤差的方差；p為因子個數；n為樣本個數。

其中第一項代表的是擬合優度，第二項代表了增加因子后的懲罰，權衡兩者，選取最小的AIC(p)作為合理的因子個數p。利用AIC準則篩選出關鍵預報因子集，用于建立回歸方程。

1.2.3 馬氏距離判別分析

馬氏距離是由印度統計學家馬哈拉諾比斯(Mahalanobis)提出的，表示數據的協方差距離。馬氏距離判別分析法是根據觀測到的樣本的若干數量特征對新獲得的樣本進行歸類、識別，判別其所屬類型的一種有效計算樣本集相似度的統計分析方法。相比于歐氏距離判別方法，馬氏距離的優勢在于不受量綱的影響，兩點間的馬氏距離與原始數據的測量無關，可排除變量間相關性的干擾。該方法的主要思想是比較樣本到各個總體的馬氏距離，然后將其判給馬氏距離最近的那個總體。

設有k個m維總體G1，G2，…，Gk，均值向量分別為μ1，μ2，…，μk，協方差矩陣分別為∑1，∑2，…，∑k，則樣本X到各組的馬氏距離為：

D2(X,Gα)=(X-μα)T∑α-1(X-μα)

(3)

α=1，2，…，k

判別規則為：若D2(X，Gi)=min1≤α≤kD2(X，Gα)則X∈Gi。

針對實際問題，當μ1，μ2，…，μk和∑1，∑2，…，∑k均未知時，可以通過相應的樣本估計值來替代，馬氏距離的具體計算步驟如下。

(4)

(5)

(6)

采用馬氏距離判別法，計算樣本到各個總體的馬氏距離，利用馬氏距離最小準則對預報期徑流的豐、平、枯類別做出定性預報。

1.2.4 多元線性回歸

多元線性回歸是用于研究一個隨機變量和多個變量之間相關關系的方法，利用多個因子X1，X2，…，Xm與對象Y建立多元回歸方程：

Y=b0+b1X1+b2X2+…+bmXm

(7)

其中b0，b1，…，bm為回歸系數，多元線性回歸方程的回歸系數采用最小二乘法來確定。根據預報因子和預報對象之間的關系可以建立相應的回歸方程，做出定量預報。

1.3 模型精度評定方法

模型的精度和穩定性采用《水文情報預報規范》(GB/T22482-2008)[11]中的確定性系數和平均絕對誤差來評定，確定性系數DC和平均絕對誤差e的計算公式見式(8)和式(9)。

(8)

(9)

2 丹江口9月入庫徑流預報

2.1 研究區域概況

丹江口水庫位于漢江中上游，丹江與漢江干流匯合以下800 m處，其地理位置處于東經106°12′～111°26′，北緯31°24′～34°11′之間，是南水北調中線工程的水源地。漢江丹江口水庫具有防洪、發電、灌溉、航運和養殖等五大功能，水庫多年平均入庫水量為394.8 億m3，控制流域面積95 217 km2,是我國大型綜合利用水庫之一。

2.2 預報因子初選

利用收集到的丹江口水庫1952-2008年的9月入庫徑流資料以及北太平洋逐月平均海溫、北半球100 hPa和500 hPa逐月平均高度場以及74項逐月環流特征量資料，分別對各要素場的因子進行相關性普查，并通過0.05的顯著性檢驗。初步篩選出北半球100 hPa逐月平均高度場因子5項、北半球500 hPa逐月平均高度場因子7項、北太平洋海溫場因子5項、74項環流特征量因子6項，建立徑流預報基礎因子集。

2.3 基于馬氏距離判別分析的長期徑流分級預報

把丹江口水庫1952-2008年9月入庫徑流量資料分為模擬期(1952-2000年)和檢驗期(2001-2008年)，并根據《水文情報預報規范》(GB/T22482-2008)[11]，將丹江口水庫模擬期的9月入庫徑流按照距平值劃分為豐、平、枯三段，詳細劃分標準見表1。

表1 豐、平、枯級別劃分標準Tab.1 standard of division of high, flat and dry

采用AIC準則進行關鍵預報因子的遴選。篩選出的關鍵因子集見表2。

從篩選的因子集來看，預報因子對于不同徑流級別的預報影響存在差異，海溫因子對偏豐徑流的預報影響較大，而100 hPa和500 hPa位勢高度對偏枯徑流的預報影響較大，未發現對豐、平、枯均有良好預報影響的因子。

表2 AIC準則篩選的關鍵因子集Tab.2 key factor sets selected by AIC

根據AIC準則篩選出的關鍵影響因子，利用多元線性回歸法，建立的分段線性回歸方程如下：

(10)

根據建立的回歸方程，可對模擬期1952-2000年的9月入庫徑流進行模擬預報。

根據篩選出的關鍵因子集，利用馬氏距離判別法進行判別分析，判別分析結果見表3。

表3 馬氏距離判別計算表Tab.3 Mahalanobis distance discriminate calculation table

訓練期1952-2000年的49年中有43年判定合格，合格率為87.8%；檢驗期2001-2008年的8年中有7年判定合格，合格率為87.5%，具有較高精度。對于檢驗期2001-2008年，根據馬氏距離判別分析的類別，代入相應的回歸方程，做出預報。

2.4 逐步回歸法預報模型

根據預報初選因子集，采用逐步回歸分析的方法進一步篩選和確定關鍵影響因子集，篩選出的關鍵影響因子集見表4。

利用篩選出的關鍵因子集構建回歸方程，預報方程見式(11)。

Y=99.96+0.206X1-0.180X2-0.142X3-0.397X4+

20.1X5+3.60X6+0.799X7+ 5.03X8-6.70X9-2.23X10

(11)

3 預報成果對比分析

利用基于馬氏距離判別的長期徑流預報模型(方法1)和逐步回歸預報模型(方法2)分別對丹江口水庫9月入庫徑流量進行模擬預報，模擬期預報成果見圖1。

表4 逐步回歸分析篩選的關鍵因子集Tab.4 key factor sets selected by stepwise regression analysis

圖1 丹江口水庫9月入庫徑流量實測值與模擬值Fig.1 The measured value and simulated value of the Danjiangkou reservoir in September

從模擬結果來看，根據《水文情報預報規范》(GB/T22482-2008)[11]中評定中長期預報精度的方案作為預報精度模型預報的評價標準：對于定量預報，水位(流量)按多年變幅的10%、其他要素按20%，要素極值的出現時間按多年變幅的30%作為許可誤差。采用多年變幅的10%(22.73 億m3)作為許可誤差，方法1的49年模擬期中有46年滿足要求，合格率為93.8%，方法2的49年模擬期中有44年滿足要求，合格率為89.8%。兩種方法的確定性系數DC和平均絕對誤差e見表5。

表5 模型確定性系數和平均絕對誤差統計表Tab.5 statistical table of deterministic coefficientand mean absolute error

在檢驗期2001-2008年的8年中，兩種方法的實測值與預報值的誤差見表6。

根據多年變幅的10%作為許可誤差，方法1中只有2003年的預報誤差超過允許值，預報的合格率為87.5%；方法2中的2003年和2007年的誤差超過允許值，合格率為75%。

表6 檢驗期實測值與預報值誤差統計表Tab.6 statistical table of measured and predictedvalues in the test period

4 結 語

(1)預報結果表明，基于馬氏距離判別的分級預報模型采用定性定量預報相結合的方法，預報的精度和穩定性優于采用單一定量預報的逐步回歸模型。

(2)利用AIC準則通過綜合考慮模型精度和復雜度來確定因子個數，對于提高預報精度和穩定性有一定的幫助，提高預報精度關鍵在于提高因子的代表性而不僅僅是增加因子的數量。

(3)基于物理成因出發篩選出的預報因子對提高預報精度有較好的效果，今后可加強預報因子之間的關聯性分析以及預報因子與預報對象之間的成因機理分析研究，進一步提高預報因子的相關性和代表性，從而減小預報的誤差。