謝 平,霍競群,桑燕芳,吳林倩,李雅晴,牛靜怡
(1.武漢大學水資源與水電工程科學國家重點實驗室,湖北武漢 430072;2.中國科學院地理科學與資源研究所陸地水循環(huán)與地表過程重點實驗室,北京 100101)
水文時間序列分析是深入了解水文現象特點與性質的重要途徑,以及探究水文過程復雜變化規(guī)律的關鍵方法[1-2]。長久以來,對水文時間序列的認識與研究都是基于一致性假設(即物理成因不變)的基本條件開展的[1]。然而,受自然因素與人為因素等的影響,許多流域和地區(qū)的水文水資源狀況發(fā)生了很大改變,導致水文時間序列不再全部滿足一致性假設[3-4],即發(fā)生了水文變異[5],表現為統(tǒng)計規(guī)律隨時間發(fā)生明顯改變,且任意時刻數值均與前期時間內的歷史值有關[6]。
依據水文變異的具體表現形式與水文序列參數的變化特征,又可分為跳躍、趨勢、周期以及相依變異等[7]。其中,相依變異是水文變異的一種常見形式,表現為水文序列后一數據與緊鄰的前一(幾)數據具有成因、數值等方面的繼承性[2]。當前,對水文序列中相依變異成分的研究和分析主要是基于數理統(tǒng)計方法進行[8-9],如Hurst系數法[10-11]、極差和輪次分析法[12]、自相關分析法和譜分析法[13]等。然而,這些方法主要針對相依變異形式的識別和檢驗,但對相依變異程度進行量化的研究則很少涉及。在實際的水文統(tǒng)計分析與計算中,不僅要明確水文序列是否發(fā)生了相依變異,還要判定相依變異的強弱程度,如果未對相依變異的顯著性進行定量評估,可能會對水資源安全評估以及水利工程建設等工作造成較大影響。相關系數可以描述序列間相關程度,其絕對值越大,則序列間的相關性越強,故可將相關系數作為度量水文變異程度的指標[14]。
常規(guī)的時間序列模型包括自回歸模型、滑動平均模型和自回歸滑動平均模型。因其具有時間相依的直觀形式和較強的適用性,故水文序列的相依性一般用它們來進行描述[15]。其中,自回歸滑動平均模型ARMA(Auto-Regressive Moving Average)是自回歸模型AR(Auto-Regressive)和滑動平均模型MA(Moving Average)的混合,且較AR和MA能更好地反映水文變量在時序變化上的統(tǒng)計特性,即具有更大的彈性,在水文學中應用更為廣泛[7]。國內外許多學者用其來描述水文現象在時間上的相依性。翟顥瑾等[16]使用ARMA模型模擬水質污染的相依隨機過程,實現對長江未來水質情況的分析預測;栗現文和楊佳等[17-18]采用ARMA模型分析地下水位埋深時間序列中去除確定性組分之后的相依隨機和純隨機組分,從而達到對未來短期地下水位動態(tài)進行準確預測的目的;Dwivedi等[19]運用ARMA等時間序列模型對降雨和徑流的隨機相依水文過程建模,從而實現對降雨的預測和徑流量的估計;Atan 等[20]利用ARMA 模型的相依性特征模擬洪水時間序列并對流域洪災風險進行了評估與分析等。然而,目前的研究主要是運用ARMA模型表征水文序列的相依特性,對于水文序列中相依變異強弱的量化研究則鮮有提及。趙羽西等[21]的研究雖涉及到對序列相依變異程度的分級量化,卻也僅適用于AR 模型,至于其它模型尚需進一步補充和完善。此外,在運用AIC、BIC 等定階準則對含有AR?MA相依成分的水文序列進行定階時,只能確定模型的自回歸階數p與滑動平均階數q的和(即p+q的值),而在AIC和BIC函數最小值不唯一或階數p和q有多種組合方式時,并不能對其做出明確判斷。
為此,本文以適用性更為廣泛的ARMA模型為例,在已有研究[21]的基礎上,以原序列與其相依成分間的相關系數為衡量標準,提出了描述水文相依變異強弱的一種方法。通過相關系數表達式的推求說明選用該指標進行變異分級的基本原理,并設計統(tǒng)計試驗說明基本原理的合理性;以ARMA(1,1)和ARMA(1,2)模型為例,指出相關系數指標與自相關系數的聯系;最后通過模擬時間序列和實測水文序列對所提方法進行驗證,并結合水文序列的相依變異程度和相依成分與原序列的擬合效率系數為具有ARMA(p,q)相依成分的水文序列提供一種定階的新思路。
2.1 相依變異分級原理一般情況下,在扣除確定性成分的影響之后,剩余水文序列xt中的相依成分可用自回歸滑動平均模型ARMA(p,q)表示為:
式中:u為序列xt的均值;p為模型的自回歸階數;q為模型的滑動平均階數;參數φ1,φ2,…,φp是待定系數中的自回歸系數;參數θ1,θ2,…,θq是待定系數中的滑動平均系數;εt,εt-1,εt-2,…,εt-q為獨立隨機變量,其標準差為σε、均值為0,且與xt-1,xt-2,…,xt-q獨立無關。
式(1)中的φ1(xt-1-u)+φ2(xt-2-u)+…+φp(xt-p-u)-θ1εt-1-θ2εt-2-…-θqεt-q為相依成分,用ηt表示,故xt能夠表示為線性疊加的形式:
式中:ut=u+εt,其均值與標準差分別為u、σε,為水文序列中的純隨機成分,在我國水文統(tǒng)計分析與計算中通常認為其服從P-Ⅲ型分布[7]。
又因為ηt與ut相互獨立,而xt的期望E(xt)=u、ut的期望E(ut)=u,故可得相依成分序列的期望值E(ηt)=0。
xt與ηt間的相關系數表達式為:
式中:
因此,由式(4)可得:
式中ση2與σx2分別為ηt和xt的方差,它們之間具有以下關系:
式中σu2為純隨機序列的方差,根據式(6)和式(7)可得:
將式(1)進行移項并整理得:
將式(9)兩邊同乘以xt-u并取數學期望得:
上式兩端除以σx2后考慮到自相關系數的定義:
可得:
式(12)兩端再同乘以σx2可得:
將式(9)兩邊同乘以εt并取數學期望得:
將式(9)兩邊同乘以εt-1并取數學期望得:
同理,依此類推,直至將式(9)兩邊同乘以εt-q再取數學期望得:
為此,式(13)可以化簡為:
將上式進行整理得:
因為σε2=σu2,結合式(8)和(18)可得:
對式(19)進行初步驗證:
當p=0時,即為滑動平均模型MA(q),式(19)可化簡為:
當q=0時,即為自回歸模型AR(p),式(19)可化簡為:
式(20)和(21)與現有結論相符合[22],故初步認為式(19)合理可靠。
ARMA 模型中的滑動平均系數θ1,θ2,…,θq和自回歸系數φ1,φ2,…,φp的參數解與樣本自相關系數存在函數關系,而樣本自相關系數能夠描述水文序列的相依變異程度[7],故相關系數r能夠表示樣本序列的相依關系并可描述其相依程度。因此,可以用相關系數對水文序列的相依程度進行分級。
2.2 相依變異分級方法基于以上推導證明的分級原理,借鑒水文變異分級方法[21,23-25],提出水文時間序列的相依變異分級方法。其中,選取rα、rβ、0.6和0.8作為相依變異分級界限,相應地把相關系數r分為5段區(qū)間,對應描述不同強弱的相依變異程度,具體分級標準見表1。rα和rβ分別是在已知模型階數和序列長度的情況下,顯著性水平為α和β(α>β)下的分級閾值,通常取α=0.05,β=0.01。
表1 相依變異程度分級
以ARMA(p,q)為例,利用所提方法對水文序列的相依變異程度量化分級的詳細步驟如下:
(1)去除水文序列中的確定性成分,得到含有獨立隨機或相依成分的剩余序列xt。
(2)判斷xt中是否存在相依成分。畫出自相關與偏相關系數圖,如果自相關系數均在容許限以內,認為剩余序列為純隨機成分;否則,認為剩余序列含有相依成分,且當自相關與偏相關系數均拖尾時,認為相依成分符合ARMA(p,q)模型。
(3)確定模型的階數。不同于AR(p)與MA(q),據自相關系數圖和偏相關系數圖并不能判斷出ARMA(p,q)模型的階數,故選用兩種定階準則(AIC和BIC準則)來判定其階數p、q。
(4)估計模型參數。用矩估計法[26]對AR?MA(p,q)的參數φi(i=1,2,...,p)和θi(i=1,2,...,q)進行估計,由此得到相依成分序列ηt。
(5)計算相關系數r。利用式(3)求得序列xt與ηt間的r值。
(6)判定變異等級。根據得出的相關系數r值對照表1找到對應的等級區(qū)間,由此判定出水文序列相依變異程度的強弱。
3.1 分級原理合理性驗證為驗證水文相依變異分級方法對ARMA(p,q)的適用性,而實際中基于時間序列模型的水文序列預測與延展工作通常以低階為主,故以較低階數(p+q)的ARMA(1,1)和ARMA(1,2)模型為例,分別設計統(tǒng)計試驗,以驗證分級原理中式(19)推導過程的合理性。
3.1.1ARMA(1,1)模型ARMA(1,1)模型可以表示為xt=φ1(xt-1-u)-θ1εt-1+ut,據式(19)可得:
用Yule-Walker方程估計自回歸系數φi(i=1,2,…,p)[27]:
由Yule-Walker方程可得ρ1=φ1,將其代入式(22)得到:
為滿足平穩(wěn)性和可逆性,各系數需符合:
統(tǒng)計試驗詳細步驟如下:
(1)利用P-Ⅲ型分布生成長為n=1000、均值u=100、變差系數Cvu=0.2、偏態(tài)系數Csu=0.4 的序列ut,然后隨機生成20組自回歸系數φ1值和滑動平均系數θ1值,代入式(1)得到20個ARMA(1,1)序列xi(i=1,2,3,…,20);
(2)相依成分ηi由所得ARMA(1,1)序列xi扣除ut求取;
(3)每組系數值內試驗1000 次,ηij為第i組第j次統(tǒng)計試驗得到的相依成分序列(i=1,2,…,20;j=1,2,…,1000);
(4)依據式(3)計算得到每組原始序列與其相依成分序列的相關系數rij(i=1,2,···,20;j=1,2,…,1000),再求出每組的平均相關系數
表2 ARMA(1,1)模型中不同參數下的ra、ri 及相對誤差(δ)
3.1.2ARMA(1,2)模型ARMA(1,2)模型可以表示為xt=φ1(xt-1-u)-θ1εt-1-θ2εt-2+ut,將p=1,q=2 帶入式(19)得:
因為由Yule-Walker方程得到ρ1=φ1,代入式(26)得:
為滿足平穩(wěn)性和可逆性,各系數需符合如下條件:
ARMA(1,2)模型的統(tǒng)計試驗步驟參照ARMA(1,1)模型,并保持參數不變。在符合式(28)要求的情況下,隨機生成20個不同φ1、θ1、θ2的系數值組合,依據模型表達式得到ARMA(1,2)序列,對比統(tǒng)計試驗得到的平均相關系數ri與由公式求出的相關系數ra,結果見表3??梢钥闯?,20組統(tǒng)計試驗中只有兩組結果的δ值超出5%,但不足6%,且參數個數越少,相對誤差δ越小,由此驗證了式(19)的合理性,說明本文所提分級原理對ARMA(1,2)模型也具有適用性。
表3 ARMA(1,2)模型中不同參數下的ra、ri 及相對誤差(δ)
此外,進行了針對ARMA(2,1)、ARMA(1,3)、ARMA(2,2)和ARMA(3,1)模型的統(tǒng)計試驗,且根據統(tǒng)計試驗結果,本文所提分級原理與方法對于其它階數p+q等于4階以內的ARMA模型也是合理的,限于篇幅,詳細分析結果不再展開贅述。因此,本文所提的分級原理與方法對較低階數的AR?MA(p,q)模型具有普適性。
3.2 相關系數與自相關系數的關系式(19)給出了r與自相關系數的關系,其成為本文所提分級方法的理論基礎。為使該分級方法的合理性更具說服力,現繼續(xù)探討相關系數與自相關系數在ARMA模型中的具體聯系,此處僅以ARMA(1,1)和ARMA(1,2)模型為例。
3.2.1ARMA(1,1)模型 利用統(tǒng)計試驗法生成具有ARMA(1,1)相依成分的序列,隨機序列ut服從P-Ⅲ型分布,其長度n=1000、均值u=100、變差系數Cvu=0.2、偏態(tài)系數Csu=0.4。為單獨研究兩種待定系數的影響,選取兩組φ1、θ1值。第一組滑動平均系數θ1取固定值0.2,自回歸系數分別取φ1=-0.9、-0.5和-0.2 ; 第二組自回歸系數φ1取固定值0.2, 滑動平均系數分別取θ1=0.2、0.5和0.9。圖1、圖2為兩組情況所得各階自相關系數。
圖1和圖2顯示, |ρ1|不僅隨著 |φ1|的增大而增大,也隨著 |θ1|增加而增加。而且,序列與相依成分的相關系數r也存在與 |ρ1|一致的變化趨勢,即隨著 |φ1|與 |θ1|的增大而增大。由此可知,對含有ARMA(1,1)相依成分的序列而言,相關系數可以刻畫其相依變異程度。
圖1 滑動平均系數θ1 固定時的各階自相關系數
圖2 自回歸系數φ1固定時的各階自相關系數
3.2.2ARMA(1,2)模型 自相關系數可以描述水文序列的相依程度[7],在ARMA(1,2)模型中,式(27)顯示r2與參數φ1、θ1、θ2有關,而滑動平均系數和自回歸系數的參數解與樣本自相關系數存在函數關系,所以r與自相關系數具有函數關系。根據式(28)以及自相關系數、自回歸系數和滑動平均系數的關系,選取多組符合條件的ρ1、ρ2的值,并由矩估計法[26]得到ARMA(1,2)模型的參數φ1、θ1、θ2值;純隨機成分的生成同3.2.1 節(jié)。ARMA(1,2)序列減去純隨機成分得到相依成分序列。求出每種ρ1和ρ2值下的相關系數,結果如表4所示。能夠看出:一階自回歸系數φ1的大小與ρ1的值無關[26,28],且當二階自相關系數ρ2固定不變時,相關系數r隨著一階自相關系數ρ1的增加而增加;反之,當ρ1固定時,相關系數r隨著ρ2的增加而減小。說明r與 |ρ1|呈正相關關系、與 |ρ2|呈負相關關系,故對于含有ARMA(1,2)相依成分的序列,其相依性強弱也能用相關系數進行描述。
表4 不同ρ1、 ρ2 下序列的相關系數
3.3 模擬序列分析驗證將所提分級方法應用于具有相依成分的模擬序列進行驗證。首先從常用的顯著性水平中選取第一顯著性水平α=0.05、第二顯著性水平β=0.01;然后根據相依變異程度分級表(即表1)中各個變異等級對應的相關系數r的取值范圍,并結合所提水文相依變異分級方法的理論基礎式(19),求取各相依變異程度所對應的整體變量的取值區(qū)間;最后選取滿足整體變量H在取值區(qū)間內的系數值組合,生成具有相依變異成分的模擬序列。下面仍以AR?MA(1,1)和ARMA(1,2)模型為例進行驗證。純隨機成分在模擬序列中的長度n=100、均值u=100、變差系數Cvu=0.2、偏態(tài)系數Csu=0.4,且服從P-Ⅲ型分布。
3.3.1ARMA(1,1)模型 在ARMA(1,1)模型中,它的自回歸部分和滑動平均部分的階數分別為p=1、q=1,根據式(22)—(24)可知此時整體變量H可表示為。各相依變異程度對應的整體變量H的取值區(qū)間見表5。
表5 ARMA(1,1)模型各相依變異程度區(qū)間所對應的整體變量H的取值范圍
在上述條件下,分別選取一階自回歸與一階滑動平均φ1=0.12、θ1=0.10,φ1=0.20、θ1=0.16,φ1=0.25、θ1=0.24,φ1=0.50、θ1=0.70 和φ1=0.80、θ1=0.95 共5 組系數值,由這5 組系數值求得的整體變量H分別對應相依變異程度的5個等級。模擬生成這5種情形下的序列,繪于圖3,并計算出所生成序列與其相依成分之間的相關系數。
圖3為相依無、弱、中、強、巨變異下模擬序列與其相依成分示意圖。從圖3可以直觀看出,各相依成分與原序列的擬合程度依次加強,表明模擬序列的相依變異程度依次加強;求得各圖的模擬序列與相依成分間的相關系數值分別為0.1556、0.2530、0.3370、0.7047、0.9004,滿足表5所列出的各級相依變異程度對應的相關系數的取值區(qū)間,也從定量角度很好地驗證了所提方法的合理性。
圖3 ARMA(1,1)模型不同H值下相依變異分級序列
3.3.2ARMA(1,2)模型ARMA(1,2)模型的階數為p=1、q=2,由式(26)、(27)可知整體變量H可表示為,各相依變異程度所對應的整體變量H的取值范圍見表6。
表6 ARMA(1,2)模型各相依變異程度區(qū)間所對應的整體變量H的取值范圍
此種條件下,分別在相依變異的5 個等級對應的整體變量H區(qū)間內各取一組滿足條件的φ1、θ1、θ2組合值,同時一階自回歸系數φ1、一階滑動平均系數θ1以及二階滑動平均系數θ2還需滿足式(28)的限制條件。選取的5 組系數值分別為φ1=0.02、θ1=0.07、θ2=0.08,φ1=0.15、θ1=0.18、θ2=0.20,φ1=0.30、θ1=0.28、θ2=0.30,φ1=0.50、θ1=-0.50、θ2=-0.40,φ1=0.80、θ1=-0.80、θ2=-0.85。所生成序列與其相依成分的示意圖如圖4所示,圖4依次為相依變異從弱到強的5 個等級,所對應的相關系數r值分別為0.1077、0.2975、0.4703、0.6842、0.9207,均滿足表6中各變異等級對應的r值區(qū)間。圖4直觀地顯示了模擬序列相依變異的程度依次增強,而相關系數r的值也依次增加,故所提分級方法的合理性得到了驗證。
圖4 ARMA(1,2)模型不同H值下相依變異分級序列
為進一步說明所提分級方法的合理性與適用性,選取實測水文序列對其進行驗證,包括黃河花園口站的年徑流序列(1919—2000年)、長江漢口(武漢關)站的月徑流序列(2010—2015年)及長江沙市(二郎磯)站的月徑流序列(1991—1998年)。按照相依變異分級方法的步驟順序,要先扣除水文序列中的確定性成分,得到含有獨立隨機或相依成分的剩余序列,可以采用綜合加權法[1]、滑動相關系數法[29]、Brown-Forsythe 法[30]等檢測序列的跳躍點;采用分段趨勢相關系數識別法[31]、過程線法、Kendall檢驗法[32]等鑒別趨勢成分;采用滑動周期相關系數識別法[33]、周期圖法、最大熵譜分析法[34-35]等檢驗序列的周期成分。
顯著性水平分別選取α=0.05、β=0.01,對以上序列進行變異診斷分析[36],其中結果如表7所示。結果表明,花園口站的年徑流序列在1932年發(fā)生了明顯的跳躍變異,這主要是氣候變異直接導致的[37];漢口站的月徑流序列和沙市站的月徑流序列都存在12 a 的周期變異,這種周期波動現象主要和區(qū)域的氣候條件及水文地理特征密切相關[24]。
表7 實測徑流序列變異診斷結果
圖5為各站點實測徑流去除確定性成分之后的自相關系數和偏相關系數圖,藍線為容許度為95%時的界限。由圖可知,花園口站的自相關系數均在容許限以內,說明花園口站年徑流序列不含有相依成分;而漢口站和沙市站序列的自相關系數和偏相關系數均呈現出拖尾的特性,表明漢口站和沙市站的月徑流序列都存在相依成分,且為ARMA(p,q)相依成分。故采用ARMA(p,q)模型對漢口站和沙市站的月徑流的剩余序列進行建模,ARMA(p,q)中的階數p,q運用AIC與BIC準則進行判定:
圖5 實測水文序列的自相關系數及偏相關系數
式中:n與σε2分別為序列長、殘差方差。AIC準則和BIC準則類似,都是確定在有最小函數值時所對應的p和q為最佳階數。
可得在p+q等于3和2時,漢口站和沙市站有最小的AIC與BIC值,由此可確定沙市站月徑流的剩余序列符合ARMA(1,1)模型,因在ARMA(p,q)模型中q的取值一般小于p[28],故初步判斷漢口站月徑流的剩余序列符合ARMA(2,1)模型,還需要根據相依成分序列與其原始水文序列的擬合效果進一步進行驗證。采用ARMA(2,1)和ARMA(1,1)模型分別對漢口站和沙市站的剩余序列建模,由矩估計法[26]估計模型未知參數,并求出相關系數,具體結果如表8所示。可以看出:花園口站年徑流剩余序列由AIC和BIC準則判定的階數為1階,但與其相依成分之間的相關系數為0.0240,比rα=0.2172小,因此該序列無相依變異,階數也最終確定為0;漢口站序列的相關系數r=0.3418,在區(qū)間[rα,rβ]內,故其相依變異程度為弱變異;沙市站序列的相關系數r=0.4022,在區(qū)間[rβ,0.6] 內,因此其相依變異程度為中變異。
表8 實測序列相依變異強弱分級
圖6為各站點實測徑流序列與其相依成分的示意圖。圖6(b)與圖6(c)為分別使用ARMA(2,1)模型和ARMA(1,2)模型對漢口站月徑流序列進行擬合得到的原始水文序列及其相依成分示意圖,其中ARMA(2,1)模型的相關系數為0.3418、ARMA(1,2)模型的相關系數為0.3343,從相關系數的大小初步判斷ARMA(2,1)模型優(yōu)于ARMA(1,2)。進一步引入一個能夠驗證水文模型模擬結果優(yōu)劣的評價指標即納什效率系數(NSE)[38],其表達式為:
圖6 各站點原始水文序列及其相依成分示意
式中:T為序列長度;Q0t為實測序列的第t個值;Qmt為模擬序列的第t個值;表示實測序列的均值。E的取值范圍為(-∞,1],E越接近1表示擬合效果越好,模型可信度越高;E接近0,表示模擬結果總體可信,但擬合過程存在一定的誤差;E遠小于0,則模型是不可信的。經計算,使用AR?MA(2,1)模型和ARMA(1,2)模型對漢口站月徑流序列進行擬合所得的納什效率系數分別為E1=-0.6104 和E2=-41.5755,考慮到漢口站序列的相依變異程度為弱變異,而相依變異程度越高,擬合效果越好,即擬合誤差越小,所以漢口站的相依成分序列與原序列間的擬合效果相比于相依中、強、巨變異要差,此為擬合過程存在誤差的原因,從而使用ARMA(2,1)模型對漢口站月徑流序列進行擬合所得的納什效率系數為E1=-0.6104 得到了很好的解釋。故根據擬合效果與漢口站月徑流序列的相依變異程度可以判定漢口站月徑流的剩余序列符合ARMA(2,1)模型,這為AIC和BIC函數最小值不唯一或階數p和q有多種組合方式時,提供了一種對ARMA(p,q)相依成分進行定階的新思路。圖6直觀地顯示,實測徑流序列與其相依成分的擬合程度逐漸加強,說明實測序列的相依變異等級越來越高,故該分級方法的適用性在實測水文序列的應用中得到了驗證。
受自然和人為等因素的影響,水文要素常常呈現出一定的相依特性,故從氣候和人為兩方面對各站徑流序列的相依變異做出物理成因分析[21]?;▓@口站位于黃河下游,它的年徑流序列無相依變異主要是其總體受氣候變化與人類活動的影響較?。?3];相較于氣候變化,人類活動對長江中游徑流變化的影響更加顯著[39],沙市水文站地處湖北省荊州市,作為長江中游荊江河段的關鍵控制站,它上游葛洲壩和三峽水利樞紐的運行和建設是導致此站月徑流中變異的重要原因,其次還受到沿程護岸工程以及洲灘守護工程的影響[40-41];漢口站地處湖北武漢,是對漢江入匯長江后的水情進行控制和監(jiān)測的重要站點,雖然受到長江上游水利工程的影響,但漢江的水利工程設施相對較少,匯入長江后,水文相依變異程度減輕為弱變異[13]。
為定量研究水文序列中相依性的強弱,本文提出了一種基于自回歸滑動平均模型(ARMA)的水文序列相依變異分級方法,分析驗證后得出結論如下:
(1)由分級原理可知,ARMA模型中的自回歸系數和滑動平均系數的大小影響著相依成分與水文時間序列間相關系數r的取值,而滑動平均系數和自回歸系數的參數解與樣本自相關系數存在函數關系,基于此建立了相關系數r與能夠表征水文序列相依程度的樣本自相關系數之間的聯系,即相關系數r可以表示樣本序列的相依關系并描述其相依程度。此外,本文所推導的相關系數r與ARMA(p,q)模型參數間的關系(即式(19))在以往的文獻中并未提及,該關系可為后續(xù)對自回歸滑動平均模型與相關系數的研究提供參考。
(2)運用統(tǒng)計試驗證明了分級原理與方法對較低階數的ARMA模型的合理性。①公式計算和統(tǒng)計試驗所得r值差別很小,證明了式(19)的合理性;②在ARMA(1,1)模型中,相關系數r和自相關系數|ρ1|均隨自回歸系數的絕對值 |φ1|和滑動平均系數的絕對值 |θ1|的增大而增大;在ARMA(1,2)模型中,r與 |ρ1|呈正相關關系、與 |ρ2|呈負相關關系,說明以相關系數作為指標衡量序列的相依性強弱是合理的;③根據整體變量H的取值范圍,生成模擬水文序列,所得相關系數處于對應變異程度的取值區(qū)間,進一步說明此方法合理可靠;④而本文所提的分級原理與方法對于更高階數的ARMA(p,q)模型是否適用尚需進一步研究與探討。
(3)實測序列分析表明:黃河花園口站年徑流序列無相依變異、長江漢口(武漢關)站月徑流序列為相依弱變異、長江沙市(二郎磯)站月徑流序列為中等程度的相依變異。并結合物理成因從氣候變化和人類活動兩個方面很好地解釋了各站點實測水文序列相依性依次加強的原因,驗證了所提分級方法的合理性及對實測水文資料的適用性。
(4)提供了一種具有ARMA(p,q)相依成分的水文序列定階的新思路:在AIC和BIC函數最小值不唯一或階數p和q有多種組合方式時,可以先根據原始水文序列的相依變異程度初步判定模型階數是否合適,再根據其與相依成分的擬合效率系數最終確定模型階數。