任東風(fēng)，胡馨予，齊 歡，董華飛

（1.遼寧工程技術(shù)大學(xué) 測繪與地理科學(xué)學(xué)院，遼寧 阜新 123000；2.沈陽市勘察測繪研究院 軟件開發(fā)部，遼寧 沈陽 110000）

0 引言

洪水預(yù)報是通過計算模型參數(shù)建立準確的水文模型.根據(jù)收到的實測降雨資料，工作人員利用建立的水文模型對流域內(nèi)的洪水流量進行預(yù)測.但由于受降雨、氣候、徑流及下墊面時空分布變異性的影響，洪水預(yù)報具有高度的不確定性.暴雨不僅在任何時候都有空間分布的變化，而且流域內(nèi)土壤、植被、地形、地貌、地質(zhì)、水文地質(zhì)等條件的空間分布是不均勻的.人類活動的影響一般是時空的[1]，故而研究洪水預(yù)報根本問題在于對降雨、氣候、產(chǎn)匯流以及下墊面等洪水基本因子的確定以及研究.朱穎元[2]等提出特定地域內(nèi)洪水發(fā)生的時間和量級具有周期變化的特征，因此歷史洪水一定程度上包含當(dāng)?shù)貧夂蚝拖聣|面變化等因素的信息，在此基礎(chǔ)上將各流域水文站洪水預(yù)報方案將降雨產(chǎn)匯流過程視為整體黑箱系統(tǒng)，建立預(yù)報根據(jù)雨量站雨量輸入-流域水文站流量輸出的關(guān)系模型進行洪水預(yù)報，利用灰色關(guān)聯(lián)分析法與多元線性回歸方程模型尋找流域范圍附近雨量站與流量之間存在的聯(lián)系.殷峻暹[3]基于GIS 的城市雨洪預(yù)報模型提出將數(shù)學(xué)計算模型與GIS 系統(tǒng)有機的結(jié)合在一起模擬雨洪模型.在灰色關(guān)聯(lián)分析方法與多元線性回歸方程在洪水預(yù)報中的應(yīng)用方面，李正最[4]認為灰色關(guān)聯(lián)分析與多元線性回歸模型在模型形式、參數(shù)辨識、方法處理等方面是一致的，這2 種方法在短期洪水預(yù)報中有效.夏軍[5]基于時間序列多信息利用的可拓原理和灰色系統(tǒng)理論的相關(guān)分析思想，提出一種適用于缺乏輸入因子數(shù)據(jù)或影響因子選取困難的中長期水文預(yù)報方法.陳意平[6]等認為灰色關(guān)聯(lián)分析為水利系統(tǒng)的中長期預(yù)報提供了一種新方法.

本文以遼河分支繞陽河流域的東白城子站為例，根據(jù)歷年洪水?dāng)?shù)據(jù)，對流量和各個雨量站進行灰色關(guān)聯(lián)性分析，判斷流域流量與雨量站之間的關(guān)聯(lián)強度，建立多元線性回歸方程，進行洪水洪峰流量方程驗證和預(yù)測.

1 研究區(qū)概況

遼西地區(qū)主要分布3 大水系，即大凌河、小凌河、遼河.遼河水系包括繞陽河、柳河、養(yǎng)息牧河、秀水河[7].繞陽河是遼河水系重要的河流之一，是遼河最下游右側(cè)的一條大支流，發(fā)源于遼寧省阜新縣扎蘭營子查喀爾山（海拔590 m）[8]，在盤錦市盤山縣萬金灘匯入遼河.東白城子站以上有韓家杖子水文站，15 處雨量站，流域面積2 138 km2，主要支流包括葦塘河和二道河，植被生長不良，土壤質(zhì)量差，水土流失嚴重，近些年來發(fā)生較大洪水年份有1969年、1975 年、1984 年、1994 年、1998 年、2005 年、2009 年、2013 年等.年降水量460～520 mm，從6月到9 月是一個汛期.這個地區(qū)有許多暴雨，其中大部分是在長期干旱之后.暴雨時空分布不均，降雨強度變化較大.這條河洪水主要是由大雨引起的.洪峰多出現(xiàn)在夜間，匯流時間短，波動大，峰高小，基本上屬于高滲透產(chǎn)流[9].徑流系數(shù)較小，一般在0.05～0.31.東白城子水文站2013 年7 月16 日洪水過程線，見圖1.

圖1 東白城子水文站2013 年7 月16 日洪水過程Fig.1 flood hydrograph of dongbaichengzi hydrological station on July 16,2013

2 研究方法

2.1 灰色關(guān)聯(lián)分析

灰色關(guān)聯(lián)分析是一種通過灰色關(guān)聯(lián)度順序來描述因素間關(guān)系的強弱、大小、次序的系統(tǒng)理論分析方法[10].該方法利用各變量之間的發(fā)展態(tài)勢的相似度或相異度來表示關(guān)聯(lián)程度[11]，相似度越高，代表變量之間的關(guān)聯(lián)程度越高.本文對雨量站和流域流量采用灰色關(guān)聯(lián)分析，計算各雨量站與流域流量的關(guān)聯(lián)值，分析每個雨量站對流量變化的影響程度，選擇關(guān)聯(lián)值超過閾值0.75 的雨量站為主要影響因素，將這些起主要作用的因素作為多元線性回歸模型的模型參數(shù).灰色關(guān)聯(lián)分析的基本原理如下：

（1）確定分析數(shù)據(jù)序列.設(shè)置1 個母序列x0(t)和若干子序列xi(t)，t=1,2,…,n，由此構(gòu)成數(shù)據(jù)矩陣[12]為

式中，子序列為自變量，母序列為因變量.

（2）數(shù)據(jù)序列的無量綱處理.利用無量綱解決父序列單位與子序列單位之間的差異，表示的數(shù)據(jù)含義不同的問題，增強數(shù)據(jù)序列間的可比性，與代表新形成的母序列和子序列[13].

（3）計算母序列與子序列的絕對差值，構(gòu)成各因子的差序列，計算式為

（4）計算母序列與子序列的關(guān)聯(lián)系數(shù)

式中，Δmin和Δmax為式（3）中最小值和最大值，k∈（0，1）為分辨系數(shù)，其值不會影響母序列和子序列的關(guān)聯(lián)度分析，本文取值為0.5[14].

（5）計算關(guān)聯(lián)度為

式中，ζ0i為母序列與子序列的關(guān)聯(lián)程度，其值越高，說明二者的關(guān)聯(lián)程度越高，擬合程度越高.

2.2 多元線性回歸模型

多元線性回歸模型是指用2 個或2 個以上因子作自變量來描述因變量的變化.根據(jù)灰色關(guān)聯(lián)分析模型計算各變量之間的關(guān)聯(lián)度，選取與分析變量關(guān)系緊密的因素，設(shè)置為因變量.通過設(shè)置自變量參數(shù)的回歸方程對因變量進行預(yù)測，模型公式為

式中，Y為流量；a與C為回歸系數(shù)；x為設(shè)定的自變量.

回歸系數(shù)可根據(jù)數(shù)據(jù)資料按最小二乘法確定.對式（6）應(yīng)用最小二乘法可導(dǎo)出

式中，xik為第i個預(yù)測因素的第k年數(shù)值（i=1,2,…,m；k=1,2,…,n）；yk為預(yù)測對象的第k年觀測值，i=1,2,…,n；n為年數(shù)；m為預(yù)測因素總個數(shù).

本文將灰色關(guān)聯(lián)模型分析挑選的主要雨量站作為自變量，使用SPSS20.0 軟件，采用自變量全部進入回歸模型方法，對因變量流域洪峰流量進行預(yù)測.該方法具有簡單便捷、推算速度快、外推性等優(yōu)點.

2.3 洪水預(yù)報誤差修正模型

實時洪水預(yù)報誤差修正，利用實時系統(tǒng)獲取的觀測信息和預(yù)報誤差可利用的其他信息，對水文模型中有一定影響的誤差因素進行實時修正是必要的[15].圖2 表示模型機與實時校正結(jié)合進行洪水預(yù)報的結(jié)構(gòu)框圖.圖2 中，I(t)和Q(t)為ε時間之前測得的模型輸入和輸出；QQ為可用于實時校正的其他信息；QC(T+L)為未校正模型計算結(jié)果；QC(T+L/T)為未校正模型計算結(jié)果[16].

圖2 實時校正預(yù)報Fig.2 real-time correction forecast

實時校正技術(shù)的研究方法很多，可分為模型誤差校正、模型參數(shù)校正、模型輸入校正、模型狀態(tài)校正和綜合校正.模型誤差修正以自回歸方法為典型例子，即根據(jù)誤差序列建立自回歸模型，然后用實時誤差預(yù)測未來的誤差[17].

自回歸修正（Auto Regression Updating）方法對模型殘差數(shù)列

采用殘差自回歸估計式

那么預(yù)報結(jié)果的校正式為

式中，

式中，et為t時刻的模型計算誤差；ξt+L為t+L時刻經(jīng)實時校正后的預(yù)報系統(tǒng)殘差；c1,c2,…,cp為模型常系數(shù)；p為模型回歸階數(shù)；為估計的t+L時刻誤差.

校正模型假設(shè)t+L時刻的模型誤差與t時刻以前的模型誤差有關(guān).誤差的預(yù)測估計式，依賴于回歸系數(shù)的確定.設(shè)已知觀測數(shù)列

模型計算數(shù)列

可得模型誤差數(shù)列

分別代入式（8）有

令

則式（11）的向量矩陣形式

式（12）的參數(shù)向量不隨時間改變，那么可用最小二乘法來確定如下

對式（13）求導(dǎo)得

2.4 繞陽河流域洪水預(yù)報模型

本文用到灰色關(guān)聯(lián)分析方法、多元線性回歸模型以及自回歸誤差修正模型，首先利用灰色關(guān)聯(lián)分析方法篩選出對河流流量影響較大的雨量站；再通過多元線性回歸模型建立多元線性回歸方程，利用自回歸誤差修正模型對多元線性回歸方程進行誤差校正；最后得出該流域降雨雨量與河流流量的回歸方程.繞陽河洪水預(yù)報模型結(jié)構(gòu)見圖3.

圖3 繞陽河洪水預(yù)報模型結(jié)構(gòu)Fig.3 flood forecast model structure of Raoyang River

3 結(jié)果與分析

3.1 數(shù)據(jù)準備

通過水文自動預(yù)報系統(tǒng)獲取降水等資料.水情數(shù)據(jù)自動采集上報系統(tǒng)是一個以遙測終端為主的實時水情數(shù)據(jù)采集、傳輸和處理系統(tǒng)、通訊網(wǎng)絡(luò)、信息中心3 部分組成[18]，見圖4.

圖4 水情自動測報系統(tǒng)架構(gòu)Fig.4 architecture of automatic hydrological forecast system

系統(tǒng)實時采集和處理阜新市各站的雨量、流量、水位等水文、水文、氣象資料，通過無線方式傳輸?shù)剿木直O(jiān)測中心的關(guān)系數(shù)據(jù)庫，再利用阜新市水情防汛值班應(yīng)用系統(tǒng)通過GIS 手段對區(qū)域內(nèi)雨量數(shù)據(jù)進行圖形顯示、查詢、分析.阜新市水情防汛值班應(yīng)用系統(tǒng)是阜新水文局聯(lián)合遼寧工程技術(shù)大學(xué)自主研發(fā)的水情防汛系統(tǒng).該系統(tǒng)主要實現(xiàn)雨情信息的查詢、顯示、統(tǒng)計分析、預(yù)警等功能，通過實時的降雨量數(shù)據(jù)生成等值面、等值線，將雨量信息清晰的顯示在地圖上，并當(dāng)降雨，流量等信息超過閾值時，發(fā)送預(yù)警短信，方便工作者分析阜新市雨情信息的空間分布情況，系統(tǒng)功能見圖5.

圖5 二維系統(tǒng)功能Fig.5 system function

本文研究東白城子站洪水預(yù)報，采用許家、雙廟、上招束溝、大興莊、東白城子、韓家杖子6 個雨量站，流域見圖6.

圖6 東白城子站流域Fig.6 drainage basin of dongbaicheng sub station

3.2 流量與影響因子灰色關(guān)聯(lián)計算

根據(jù)2013 年7 月13 日和2005 年8 月10 日2次洪水過程線和數(shù)據(jù)庫中水情數(shù)據(jù)為依據(jù)，提取降雨量和洪水流量數(shù)據(jù)，洪水過程線見圖7.圖7 中Y軸的降雨量代表時段范圍內(nèi)的面雨量，指東白城子站流域內(nèi)各個雨量站在時段內(nèi)測得的降雨量推求出來的區(qū)域內(nèi)點平均降雨量；流量是指流域在某一時刻實測的洪水流Δt量大小.降雨經(jīng)過復(fù)雜的產(chǎn)流和匯流兩個過程注入流域從而影響流域流量，但降雨產(chǎn)匯流過程受到距離、地形地貌以及下墊面等因素的影響，因而降雨量與流域流量之間具有一定的時間滯后性，洪水洪峰的形成相對于降雨也有一定的時間滯后性，即降雨匯流時間Δt.

圖7 2013 年繞陽河?xùn)|白城子站洪水過程Fig.7 flood hydrograph of East Baichengzi Station of Raoyang River in 2013

本文只研究各雨量站降雨量與流量之間的關(guān)系，不考慮地形地貌、下墊面等因素的影響，應(yīng)剔除時間滯后性的影響，因此將降雨量值平移Δt時長，使降雨量峰值與洪水-峰值對應(yīng).根據(jù)阜新市水文局提供的洪水資料數(shù)據(jù)，2013年洪水中Δt為9 h，2005 年洪水中Δt為9 h，圖7 中紅色線代表平移后的降水量，再依據(jù)降雨量與流量一一對應(yīng)關(guān)系，讀取各站點降雨數(shù)據(jù)及對應(yīng)流量數(shù)據(jù)，作為灰色關(guān)聯(lián)分析和多元線性回歸模型建立的數(shù)據(jù)基礎(chǔ)，見表1.

表1 降雨量-流量數(shù)據(jù)Tab.1 rainfall-flow data

東白城子站流量的變化與多個雨量站降雨量值相關(guān)，本文選取6 個與繞陽河流域流量相關(guān)的雨量站，通過這6 個雨量站的灰色關(guān)聯(lián)分析計算，最終確定各雨量站對繞陽河流域流量的影響程度，用于多元線性回歸方程設(shè)定參考系數(shù)，從而進行洪水流量預(yù)測.設(shè)定流量為隨時間變化的母序列x0(t)，將許家、雙廟、上招束溝、大興莊、東白城子、韓家杖子雨量站作為影響因素子序列xi(t).

將流量及其影響因素作灰色關(guān)聯(lián)分析的過程：首先設(shè)置數(shù)據(jù)序列，見表2；再利用表2 和式（2）對數(shù)據(jù)序列進行無量綱化處理，獲得初始值矩陣見表3；根據(jù)式（3），對初始值矩陣求絕對差值序列見表4，確定該矩陣中最大值和最小值（其中Δmin=0，Δmax=3.250 3），代入式（4），獲得關(guān)聯(lián)系數(shù)矩陣，結(jié)果見表5；最后，根據(jù)式（5），得到因素關(guān)聯(lián)度值，見圖8.

表2 流量及其影響因子Tab.2 flow and its influencing factors

表3 初始值Tab.3 initial value

表4 絕對差值Tab.4 absolute difference

表5 關(guān)聯(lián)系數(shù)Tab.5 relevance coefficient

圖8 洪水流量與影響因子的關(guān)聯(lián)度值Fig.8 relevance degree value of flood discharge and influencing factors

通過對繞陽河流域流量影響因子的灰關(guān)聯(lián)分析可知，關(guān)聯(lián)度大小依次為雙廟＞上招束溝＞大興莊＞韓家杖子＞東白城子，其中，雙廟、上招束溝對繞陽河流域流量的影響效果較大.流域流量與附近各個雨量站的關(guān)聯(lián)度值都超過閾值0.75，因此6 個雨量站都作為多元線性回歸模型自變量參數(shù)對繞陽河洪水流量進行預(yù)測.

3.3 洪水預(yù)報多元線性回歸模型

利用多元線性回歸模型預(yù)測阜新市東白城子站的洪水流量，選擇6 個雨量站因素作為輸入的模型參數(shù)，使用SPSS20.0 統(tǒng)計軟件，選取全部進入法計算回歸方程

式中，X1,X2,X3,X4,X5,X6分別為許家、雙廟、上招束溝、大興莊、東白城子、韓家杖子雨量站降雨量.由于該方程的回歸系數(shù)R為0.998，R2為0.996，說明影響因子之間有較強的相關(guān)性，回歸方程有意義，可以采用方程對東白城子站洪水流量進行預(yù)測.

3.4 回歸方程的檢驗

通過建立的洪水流量影響因子的回歸模型，根據(jù)東白城子站2003 年和1994 年的2 次洪水過程中的降雨量，推測兩場洪水的洪峰流量，并將預(yù)測結(jié)果與實際值進行對比與驗證，驗證結(jié)果見表6.由表6 可以看出，洪峰實際值與預(yù)測值的誤差在許可誤差范圍內(nèi)，實驗結(jié)果顯示預(yù)測效果良好.

表6 洪水洪峰擬合度分析Tab.6 analysis of flood peak fitting degree

4 結(jié)論

根據(jù)阜新水文局提供的各流域歷年洪水實測數(shù)據(jù)，分析雨量站降雨量和洪水流量，得到結(jié)論：

（1）采用灰色關(guān)聯(lián)模型計算各雨量站對流域流量波動的影響程度，選取關(guān)聯(lián)度高的雨量站作為主導(dǎo)因素，通過SPSS 軟件對主導(dǎo)因素和流量建立多元線性回歸模型，對洪水的洪峰流量進行有效預(yù)測.

（2）對比實測值和預(yù)測值的大小，確定回歸方程的實用性.驗證結(jié)果表明：多元線性回歸方程適用于遼西季節(jié)性河流的洪峰流量預(yù)測，該洪水預(yù)報模型可為遼西季節(jié)性河流的洪水預(yù)報提供一定的參考.