基于小世界復雜網(wǎng)絡的水文診斷模型仿真

裴源博，趙雪花

(太原理工大學水利科學與工程學院，山西 太原 030024)

1 引言

自然環(huán)境在氣候與人類各種行為活動的影響下，始終處于不斷變化的狀態(tài)，其中水文系統(tǒng)也會出現(xiàn)明顯改變。現(xiàn)階段，水資源短缺、承載能力下降已經(jīng)嚴重威脅社會可持續(xù)發(fā)展。由于水文系統(tǒng)的平衡關系被破壞，水文序列表現(xiàn)出的一致性特征也逐漸消失，統(tǒng)計規(guī)律將無法用于對未來情況的預測。非一致性、變異性已經(jīng)成為水文序列的新特征。當水文系統(tǒng)出現(xiàn)異常時，也代表著水循環(huán)過程遭到破壞，這會導致較多的工程和科學問題。例如，因水文系統(tǒng)出現(xiàn)變異，致使已有關于水文分析的成果遭到質疑；增加水文序列演變的復雜性，進而造成水資源管理的不可控性。這些改變對傳統(tǒng)的水文預測工作敲響警鐘，一些科研人員已經(jīng)將如何準確診斷水文序列的變異情況提上研究日程。

文獻[1]提出基于滑動偏相關算法的水文序列變異診斷方法。選取塔里木河流域的某時間段徑流資料，利用Pearson相關系數(shù)分析徑流和天氣狀況之間的相關性；引入滑動窗口技術完成水文氣象要素的聯(lián)合序列異常診斷，并通過雙累積曲線方式完成驗證。文獻[2]提出基于自回歸模型的水文序列異常識別算法。將初始序列及其相依成分的相關系數(shù)當作擬合指標，利用信息熵函數(shù)形式作為自回歸模型準則，構建診斷模型。除上述方法外，還有學者利用Mann-Kendall秩次檢驗法對水文序列的突變性與非一致性進行診斷。通過該方法識別突變點，同時采用Morlet小波分析序列周期特征，結合突變點診斷、周期分析與驅動機制三個方面分析序列的非一致性變化特點。

水文序列作為觀測樣本，在氣候與人文因素影響下會出現(xiàn)突變現(xiàn)象，通常情況下，突變具有多尺度特征，不同原因會導致不同時間與不同等級的突變。而上述方法并沒有結合水文序列的非線性特征，雖然能夠診斷出突變情況，但是診斷結果并不全面。復雜性正是當前水文研究的主要特征，為此，本文利用復雜性理論實現(xiàn)非一致性水文序列突變診斷。水文系統(tǒng)的復雜性主要表現(xiàn)在它是由天文、大氣、生物與人類共同作用與交叉的結果，應結合多方面進行考量。分析在不斷變化的環(huán)境中，水文序列表現(xiàn)出的非一致性現(xiàn)象，為突變診斷提供更加可靠的方法。

2 水文序列突變的驅動因子與統(tǒng)計特征分析

2.1 突變驅動因素分析

在不斷變化的自然界中，水文變化一般由多種因素相互作用形成。這些因素之間的作用關系如圖1所示。

圖1 驅動因素關系圖

1)自然變化

在人類干擾程度較小的情況下，水文變化一般受到自然因素的驅動，呈現(xiàn)出較為平穩(wěn)的波動趨勢，且不斷循環(huán)。其時間序列便于預測與管理，不會出現(xiàn)較大幅度改變。

2)氣候變化

近年來全球變暖已經(jīng)引起全世界的關注，相關研究顯示，暴雨、干旱等極端天氣發(fā)生頻率正逐年提高。而氣候變化大多與人類活動相關，可利用區(qū)域水循環(huán)來改變水文系統(tǒng)。如果不考慮人類活動對其產(chǎn)生的影響，氣候驅動因素對水文系統(tǒng)的影響微乎其微。

3)人類活動

自首座水壩建成以來，水文就與人類活動有著密不可分的聯(lián)系。現(xiàn)階段，人類活動對世界上82%以上的水資源均存在直接影響。此外，大多數(shù)人口都面臨非常嚴重的水資源威脅，同時也影響著生物的多樣性。由于人口迅猛增長，人們對水資源的依賴更加強烈。工業(yè)與經(jīng)濟發(fā)展也對水資源的需求量進一步增大，供需矛盾愈演愈烈。

社會行為的特征與差異性，使人們在面對水文問題的觀點與做法上形成區(qū)別。所以，必須從根本上分析人類活動的影響，構建行為和用水規(guī)律之間的關系。針對當前情況，人類行為對水文產(chǎn)生的影響分為直接與間接兩種。其中，前者代表從河流取水，用于生活、灌溉與工業(yè)等方面，主要目的是提高生活水平；后者則是通過人類一系列活動導致水環(huán)境產(chǎn)生改變，例如林地與土地建設等均會造成產(chǎn)匯流變化[3，4]。

綜上所述，水文變異并非受到單個因子驅動，而是在自然與社會共同作用下產(chǎn)生的，且二者相互影響形成一個復雜反饋過程。

2.2 水文突變統(tǒng)計特征表現(xiàn)

現(xiàn)有水文設計多數(shù)基于一致性假設，但是受到各類因素影響，一些地區(qū)表現(xiàn)出的水文一致性被破壞，這對水文設計產(chǎn)生嚴重影響。因此，變異診斷已經(jīng)成為水文設計的基礎，診斷結果的準確性關乎水利工程建設方案是否可行。如果方案出現(xiàn)問題，會有嚴重事故發(fā)生，對人們生命與財產(chǎn)安全帶來威脅。為準確判斷出水文序列是否出現(xiàn)變異現(xiàn)象，必須對水文系統(tǒng)的特征有一定了解。

將水文序列描述為x1，x2，…，xt，通常記為{X(t)}。其中，x1，x2，…，xt代表序列的t個實際值。結合時間序列相關理論，可利用下述統(tǒng)計特征對水文序列進行描述。

1)概率分布

通過分布函數(shù)[5]表示某水文過程的所有統(tǒng)計特征。若水文序列為任意變量，則一維分布函數(shù)表示為

F(x，t)=P[X(t)≤x]

(1)

與其對應的一維密度函數(shù)為

f(x，t)=?F(x，t)/?x

(2)

同理，X(t)的n維分布函數(shù)公式如下

F(x1，…，xn，t1，…，tn)=P[X(t1)≤x1，…，X(tn)≤xn]

(3)

2)均值

均值能夠體現(xiàn)序列平均水平，屬于水文序列的主要位置特征，其表達式如下

(4)

3)方差與變異系數(shù)

方差與變異系數(shù)分別表示水文序列的絕對、相對離散情況，計算公式如下

σ(t)2=E[X(t)-μ(t)]2

(5)

Cv=σ/E(X)

(6)

3 基于復雜性理論的水文序列突變診斷

3.1 復雜網(wǎng)絡模型構建

由于水文系統(tǒng)受到較多復雜因素影響，本文結合復雜性理論，建立復雜網(wǎng)絡。該網(wǎng)絡的拓撲無需考慮節(jié)點位置與邊的形狀。

將網(wǎng)絡抽象為網(wǎng)絡模型圖G(V，E)，其中，V代表節(jié)點集合，網(wǎng)絡中節(jié)點數(shù)量為|V|；E(1e，2e，…，me)表示邊集合，記為|E|，節(jié)點集合中的點與邊互相對應。

假設該網(wǎng)絡中共存在8個節(jié)點與9條邊，則從下述幾個維度對該網(wǎng)絡進行分析。

1)最短路徑

節(jié)點i與j的最短路徑為dij，描述由i到j的全部路徑中包含邊數(shù)最少的路徑，計算公式如下

(7)

式中，N為節(jié)點總數(shù)量。

2)度及度分布

節(jié)點的度ki就是與節(jié)點i相連的全部鄰居數(shù)量，表達式為

(8)

3)集聚系數(shù)

在該網(wǎng)絡中，任何節(jié)點之間的連線都具有隨機性[6]，其集聚系數(shù)表示為C=p，其中，p代表兩點相連的概率。對于該系數(shù)的分析可以得出網(wǎng)絡復雜程度等有用信息。集聚屬于網(wǎng)絡基本特征屬性，Ci即為網(wǎng)絡中各點集聚系數(shù)C的平均值

(9)

式中，ki代表節(jié)點i的鄰居數(shù)量，Ei代表ki個點具有的真實邊數(shù)。

4)節(jié)點介數(shù)

在復雜網(wǎng)絡模型中節(jié)點介數(shù)是體現(xiàn)該節(jié)點重要性的指標，計算公式如下

(10)

式中，njk表示節(jié)點j和k的最短路徑總數(shù)。

綜合上述維度可知，復雜網(wǎng)絡模型并非完全規(guī)則或完全隨機。因此，本文建立一個介于規(guī)則與隨機之間的小世界復雜網(wǎng)絡模型，建模過程為：

1)建立包含N個節(jié)點的最近鄰耦合模型，確保節(jié)點和其K/2個鄰居節(jié)點相連。

2)做隨機重連，在保證一個節(jié)點不變的情況下，根據(jù)概率p實現(xiàn)邊的重新連接。將規(guī)則的模型轉換成小世界網(wǎng)絡模型，并通過p來調節(jié)大小。當p=0時，屬于規(guī)則模型，當0

(11)

3.2 變異診斷

利用上述構建的復雜網(wǎng)絡模型，通過最小二乘算法確定模型目標函數(shù)[7]，獲取全局最優(yōu)解，即為變異點位置。該方法最大優(yōu)勢是決策者能夠設計不同變異點間存在的最小距離，完成準確診斷。診斷過程如下：

1)假設水文時間序列表示為{y(ti)}i=1，…，n。

2)當使用最小二乘算法完成未知變異點時間序列擬合時，如果要使殘差平方和為零，則最優(yōu)擬合策略為將序列中全部點進行連接。所以，必須設置約束條件，控制分段擬合數(shù)量與步長。假設序列{y(ti)}分為m′+1段，m′為變異點數(shù)量，最優(yōu)分段線性方程表示為

(m′=0，…，m；Tm′≤ti≤Tm′+l)

(12)

式中，Tm′代表變異點位置，設置T0=0，Tm′=tn，b表示回歸系數(shù)矢量。

假設水文序列中出現(xiàn)3個變異點

y(ti)={y(t1)，y(tbp(2))，y(tbp(3))，…，y(tbp(m′))}

(13)

式中，tbp(2)、tbp(3)與tbp(4)描述變異點出現(xiàn)的位置，分別構建如下回歸方程

y(ti)=a1ti+c1，i=1，…，bp

(14)

y(ti)=a2ti+c2，i=bp(2)+1，…，bp(3)

(15)

y(ti)=a3ti+c3，i=bp(3)+1，…，bp(4)

(16)

y(ti)=a4ti+c4，i=bp(4)+1，…，m

(17)

根據(jù)變異點具有的連續(xù)性特征，獲得下述聯(lián)立方程

c2=c1+(a1-a2)tbp(2)

(18)

c3=c1+(a1-a2)tbp(2)+(a2-a3)tbp(3)

(19)

c4=c1+(a1-a2)tbp(2)+(a2-a3)tbp(3)+(a3-a4)tbp(4)

(20)

4)為避免變異點區(qū)間設置存在主觀性，利用構建的復雜網(wǎng)絡模型分別獲取變異點的置信區(qū)間[8]與變異點前后時間序列的顯著性。

假設變異點兩側的回歸方程表示為：

RLa=Aa(X-AvXa)+AvYa

(21)

RLb=Ab(X-AvXb)+AvYb

(22)

式中，Aa與Ab分別代表變異點前、后回歸系數(shù)，AvXa和AvXb分別描述變異點前、后全部數(shù)據(jù)的橫坐標均值。

由于

(23)

因此，變異點左右兩側置信的上、下限分別表示為

L上=(X，Yla+ts.StDevYc)

(24)

L下=(X，Yla-ts.StDevYc)

(25)

R上=(X，Ylb+ts.StDevYc)

(26)

R下=(X，Ylb-ts.StDevYc)

(27)

根據(jù)變異點的置信區(qū)間，準確確定序列變異點位置，完成序列突變診斷。

4 實驗數(shù)據(jù)分析與研究

為證明所提方法對于水文序列突變診斷的可行性，進行仿真。將渭河流域當作診斷目標，該段流域長度是734km，面積為13.11萬km2。此地地形差異較大，呈現(xiàn)出東低西高趨勢，最大高程差為2500m之上，流勢變化緩慢，河道淤泥較多。此外，水系發(fā)達，支流眾多，但近年來的徑流量出現(xiàn)整體下降趨勢。具體情況如表1所示，水文循環(huán)示意圖見圖2。

表1 渭河流域多尺度特征表

圖2 渭河流域水文循環(huán)特征示意圖

由于本文方法是基于復雜性理論的分析，因此結合該地區(qū)實際情況，選擇診斷指標，將該指標引入到復雜網(wǎng)絡模型中，再通過回歸變異診斷，輸出最終診斷結果。本文選擇的指標如表2所示。這些指標不但要體現(xiàn)流域特征，還要保證指標之間存在一定影響與制約。

表2 該地診斷指標表

結合上述診斷指標，利用文獻[1]、文獻[2]、文獻[3]方法與本文方法對淮河流域2011年-2019年的水文序列突變情況進行診斷，診斷結果如圖3-6所示。

圖3 文獻[1]方法診斷結果

圖4 文獻[2]方法診斷結果

圖5 文獻[3]方法診斷結果

圖6 本文方法診斷結果

相關水文資料表明，渭河流域在2011-2019年期間，共出現(xiàn)三次水文序列突變狀況，時間點分別是2014年、2016年與2017年。由圖6可知，本文方法在這三個時間點上的變異系數(shù)尤為突出，而文獻[1]和文獻[3]的方法都少診斷了一個變異點，文獻[2]方法的變異曲線始終較為平穩(wěn)，雖然變異系數(shù)有所改變，但趨勢較弱，并不明顯。資料顯示，在2014年，該流域煤炭開采大幅增加，而在2016與2017年均修建了多個水庫，洪峰減少。因此，本文診斷結果與歷史資料相符，表明該方法診斷結果更加可靠。

5 結論

為準確診斷出水文序列的變異情況，本文構建了復雜網(wǎng)絡模型，再結合回歸變異算法完成序列變異診斷。在實驗過程中，綜合當?shù)厮奶卣鳎O置多個診斷指標，給出可靠的診斷結果。現(xiàn)階段，基于復雜性理論的水文序列突變診斷還處于起步階段，尚未形成全面的診斷體系。在今后研究中，需引入混沌分析等方法，提高診斷速度，將多診斷方法相結合，全面體現(xiàn)序列的變異特征，避免單一方法得出的結果不夠科學合理。