王 超

（河北省石家莊水文水資源勘測局，石家莊 050051）

水資源預(yù)測是區(qū)域水資源規(guī)劃管理的基礎(chǔ)，受氣候變化密切影響，水資源表現(xiàn)出時間波動性與突變性，這給預(yù)測研究帶來一定挑戰(zhàn)。本文基于時間序列數(shù)據(jù)集的新型建模方案 （empirical dynamic modelling，EDM），闡述模型原理與建模過程，為水資源預(yù)測開拓新思路。

1 EDM原理與數(shù)據(jù)處理

經(jīng)驗動態(tài)模型 （empirical dynamic modelling，EDM）是針對時間序列環(huán)境變量的動態(tài)建模分析方法，該方法擴展了傳統(tǒng)的復(fù)雜矩陣中定義為平衡系統(tǒng)到動態(tài)系統(tǒng)過程中處于平衡狀態(tài)［1-2］。 傳統(tǒng)雅克比矩陣是一個線性化方程組，矩陣之間部矢量交互通常是安裝一定模式設(shè)計的固定系數(shù)。 然而在實際生態(tài)系統(tǒng)中幾乎不存在靜態(tài)平衡，更為普遍的是表現(xiàn)為非線性行為。 EDM框架強調(diào)系統(tǒng)狀態(tài)隨生態(tài)系統(tǒng)向量狀態(tài)的變化而變化，并具有時間動態(tài)序列特征。EDM的內(nèi)涵豐富，主要包含非線性動力學(xué)識別、CCM動態(tài)歸因、生態(tài)序列預(yù)測等內(nèi)容。

本文主要應(yīng)用EDM原理中的S-map和CCM算法進行論述。

1.1 S-map算法

S-map是基于時間滯后信息的一種系統(tǒng)動力學(xué)回歸模型，其框架如式（1）：

S-map通過線性逼近已經(jīng)完成在對引子流形上的每個局部位置x（t*）的權(quán)重設(shè)計［3］。 在吸引子流形上，對于更接近目標點x（t*）的向量，權(quán)重更大，其中參數(shù)k是權(quán)重的控制參數(shù)之一，d為平均距離，權(quán)重計算如式（2）：

1.2 CCM算法

交叉收斂映射CCM （Convergent Cross Mapping，CCM）算法是EDM框架的重要內(nèi)容之一，由Sugihara等［4］于2012年在《Science》上公開發(fā)表的新型因果關(guān)系檢測方法。 區(qū)別于格蘭杰等因果關(guān)系檢驗的互相關(guān)函數(shù)與信息熵傳遞內(nèi)核。 CCM是基于Takens嵌入定理，Mx和My是同胚胎中M中的影子流形，它們均是通過對原始序列x，y的時滯信息重構(gòu)得到的，由此產(chǎn)生了Mx與My之間的連續(xù)映射，若Mx對My有因果關(guān)系，則能夠在My上找到Mx的過程標簽。 對于n維系統(tǒng)M中關(guān)聯(lián)變量Mx，My，構(gòu)建影子流形如式（3）：

式中 E為嵌入維數(shù)；τ為時滯步長。

若在Mx上找到系列與My，k0對應(yīng)的流形，My，k0流形平均距離函數(shù)如式（4）～式（5）：

式中 k為時間序列集的長度，i，j∈k。

若Mx在My上的映射收斂到0～1之間，表明Mx是My的歸因。

1.3 數(shù)據(jù)處理

通過《河北年鑒》中提取衡水市水資源總量年值（1970～2018年），并以此作為目標變量進行預(yù)測分析。水資源總量（WR）的變化與氣候變化密切相關(guān)，通過文獻分析［4-5］，選擇若干高頻變量作為預(yù)測因子，分別是年降雨量 （Pre）、 年平均溫度 （Tem）、 濕度（Mois）、日照時數(shù)（Hour）、蒸發(fā)量（ET），氣象因子數(shù)據(jù)來源于中國氣象中心，運用R語言實現(xiàn)水資源的預(yù)測分析研究。

2 結(jié)果與分析

2.1 基于S-map的非線性動力學(xué)

EDM方法的應(yīng)用性具有嚴格標準，不僅要求數(shù)據(jù)集為連續(xù)性的時間序列變量，還要求事物本身符合動力學(xué)特征。 rEDM開源程序包中設(shè)計了相關(guān)檢測程序，通過S-map算法中的權(quán)重參數(shù)k進行迭代，對比不同迭代結(jié)果中預(yù)測值與實際值的相關(guān)系數(shù)（rEDM中規(guī)定相關(guān)系數(shù)用rho標識）。 當k為0時，狀態(tài)系統(tǒng)中吸引子流形各點的權(quán)重相同，此時的S-map算法簡化為標準的線性模型，當k大于0時，則為非線性預(yù)測模型。 如圖1（b），rho在k大于0時的值均高于k為0時rho的值，說明線性動力學(xué)的存在，即運用EDM原理建立衡水市水資源總量的動態(tài)模型具有良好的適用性。另外，時滯步長表征了非線性預(yù)測的長度，如圖1（a），當步長為5時，rho達到最高水平，因此在后續(xù)的預(yù)測中選定該參數(shù)為5。

圖1 非線性動力學(xué)識別

2.2 相關(guān)參數(shù)設(shè)定

嵌入維數(shù)定義為重構(gòu)吸引子流形所需的變量維數(shù)，其不僅僅所指預(yù)測變量個數(shù)，包含相關(guān)變量的時滯子集，如圖2，圖中橫坐標為嵌入維數(shù)大?。‥），縱坐標為基于S-map 的預(yù)測值與實際值的相關(guān)系數(shù)rho，rho越大，說明嵌入維數(shù)最優(yōu)。 其中在WR，Pre兩變量中，最佳嵌入維數(shù)為8和2；在WR和Tem組合中最優(yōu)值為8和4；Mois在最佳嵌入維數(shù)為4；Hour和ET則均為7。 嵌入維數(shù)的大小決定模型的復(fù)雜程度，若其值過大，則會導(dǎo)致狀態(tài)空間過密和信息重疊，并損失動力結(jié)構(gòu)性規(guī)律； 反之過小則產(chǎn)生狀態(tài)空間依賴性不足，預(yù)測精度較低。 通過rEDM程序號中的simplex函數(shù)可給出最佳嵌入維數(shù)的平衡選擇。

圖2 嵌入維數(shù)選擇

2.3 基于CCM的水資源變化歸因

經(jīng)參數(shù)優(yōu)化后運用rEDM中的CCM_test函數(shù)對衡水市水資源總量WR與氣候環(huán)境因子之間的歸因識別。 如圖3，圖中橫坐標為時間序列數(shù)據(jù)集的時間變化； 縱坐標為各單變量通過時滯構(gòu)建吸引子流形后的點值與原來點值的相關(guān)性系數(shù)。

圖3 基于CCM的雙變量動態(tài)交互

圖3（a）中Pre在WR上的投影系數(shù)值為0.3，WR在Pre上的投影系數(shù)為0.2，均介于0～0.5，表明二者互為動力學(xué)因果關(guān)系，即降水量是水資源量的因，同時水資源量也是降水量的因，說明二者互為影響。 其中Tem map WR 0.01，說明溫度的變化是區(qū)域水資源變化的顯著性成因（P＜0.05）；同時水資源量的變化也會引起溫度的變化（WR map Tem 0.2），但這一過程并不顯著。 大氣濕度與水資源量之間互為成因，其中濕度對水資源量的投影值（Mois map WR）為0.03；水資源量對濕度的投影值（WR map Mois）為0.02；歸因均達到顯著水平。 日照時數(shù)對水資源量具有一定影響 （Hour map WR 0.3），前者是引起后者變化的原因；而后者不會引起前者的變化（WR map Hour 0.6，P＞0.5）。 圖3（e）中ET與WR交叉收斂特征十分明顯，二者互相投影值分別為0.01，0.02；通過5%水平信度檢驗，說明蒸發(fā)變化是衡水市水資源量變化的顯著性原因；同時后者也是前者的成因。

2.4 基于multuview的水資源預(yù)測

rEDM程序中基于EDM原理設(shè)計了圖4中非線性動力學(xué)預(yù)測方法，分別是Simplex、S-map、multivariate、multi-view等，Chun-Wei Chang等［6］研究表明以multiview的預(yù)測效果最佳，因此本文選擇該方法對1971～2018年衡水市水資源總量變化進行時間性預(yù)測，結(jié)果如圖4。

圖4 基于multi-view的衡水市水資源總量預(yù)測

從圖4（a）可知，multi-view模擬的水資源總量年值能夠較好地擬合其實際變化，雖然有誤差存在，但均較小；依據(jù)圖4（b）可知，實際值與預(yù)測值之間的相關(guān)系數(shù)高度0.95，Mae和RMSE為0.71億，1.01億m3，表明預(yù)測精度較高，預(yù)測結(jié)果可靠。

3 結(jié)語

（1）水資源是氣候系統(tǒng)耗散過程中的矢量之一，水資源的時空變化與氣候變化密切相關(guān)，通過建立氣候因子與水資源總量之間的關(guān)系模型能夠較好地開展水資源預(yù)測。

（2）介紹了應(yīng)用EDM原理對衡水市水資源進行預(yù)測的實證方案，建模過程表明水資源時序變化更加符合非線性特征，相對于傳統(tǒng)線性模型，EDM模型具有一定的適用性；EDM框架中的CCM算法檢測出了該市水資源變化的成因，從投影值來看，蒸發(fā)與溫度的影響其顯著成因（P＜0.05），雖然大氣濕度、降水量、日照時數(shù)對其具有一定影響，但并不顯著。