王渤權,王麗萍,唐 勇,劉明浩,李傳剛
(1.華北電力大學 可再生能源學院,北京 102206;2.國家電網公司國家電力調度控制中心,北京 100031)
我國具有十分豐富的水能資源,位居世界的第一位,多利用水能資源從而使水電廠多發電,以節約煤炭資源的消耗,這一舉動在當今資源短缺的時代有著重要的意義,這也是電力系統對水電部門提出的要求[1]。而水電站發電量的多少主要取決于徑流過程,目前水電站工作人員多采用徑流預報值來做發電計劃的制定,然而,由于河川徑流的隨機性、地區性,常常出現預報值不夠精準,預報誤差較大的情況,由此便產生一定的出力誤差,影響著發電計劃的制定,從而導致被迫修改發電計劃。由于目前國內外尚未有對由于預報誤差而引起的出力誤差進行相關研究,為方便水電站工作人員充分掌握由于預報級別以及預報誤差產生的出力誤差信息,以便能夠及時地對發電計劃進行制定與修改,本文對預報值、預報誤差以及出力誤差之間的關系進行探討,建立三者的三維聯合分布模型。對于聯合分布的求解,國內外學者進行了大量的工作,費永法[2]基于事件積原理,采用一種解析算法計算出不同河流洪水遭遇的概率問題,一定程度上解決了多元隨機變量的聯合分布問題。戴昌軍、梁忠民[3]通過將正態化變換的Moran方法、傳統經驗頻率方法以及多維轉化為一維的方法應用在水文領域中進行比較分析,結果表明Moran方法具有較好的擬合度。王文圣、丁晶[4]采用核估計理論來構造一種多變量的非參數模型,并且將其應用于水文學中,為隨機水文學的發展提供了新的思路。鄭紅星、劉昌明[5]采用了經驗頻率的方法對南水北調中不同的水文區域降水遭遇問題進行了計算與分析。以上專家學者們對多維隨機變量的聯合分布的研究均做出了一定的貢獻,但這些方法對于求解預報值、預報誤差以及出力誤差的聯合分布均存在一些不足之處:以上方法對于邊緣分布函數要求多為正態分布并且有些方法對于兩者的相關性為線性要求,而對于預報誤差與出力誤差而言,其分布函數不一定為正態分布,尤其對于出力誤差分布的研究甚少,其分布未知,并且二者存在著復雜的數學關系,呈現出非線性的關系,故不適合采用上述方法進行求解計算。最大熵法是依據樣本信息來推求未知分布的一種推斷方法,能夠客觀地反映事物的本質,可以擬合任意形式的分布,文獻[6]證明對于預報值以及預報誤差的分布特征均可應用最大熵法進行描述,并且文獻[7]指出運用最大熵方法擬合出來的最大熵分布函數較正態分布函數能夠更好地反映預報誤差的特點,而Copula函數[8-10]是通過邊緣分布函數與相關性結構兩方面構造的函數,并且適用于任何邊緣分布函數,對其變量間相關性的要求也比較廣泛,線性非線性均可,具有更強的靈活性與適應性。
基于此,本文將最大熵理論與Copula函數理論結合,建立三維聯合分布模型來求解預報值、預報誤差以及出力誤差的三維聯合分布,在此基礎上結合貝葉斯理論進而計算出在一定預報徑流值與預報誤差的情況下其出力誤差的條件概率的大小,為水電站工作人員調度管理提供理論依據。
針對出力誤差分布未知的情況,利用最大熵理論中可以求解出任意形式的分布,以及Copula函數中適用于任何形式邊緣函數的特點,將最大熵理論與Copula函數理論進行結合,從而建立三維聯合分布模型對預報值、預報誤差以及出力誤差的聯合分布進行求解。
本文建立的三維聯合分布模型形式如下所示:
P(X1≤x1,X2≤x2,X3≤x3)=C(u1,u2,u3)=
φ-1[φ(u1;θ)+φ(u2;θ)+φ(u3;θ)]
(1)

對于上述模型的求解,主要是針對兩個未知參數λij和θ的估計,其值的好壞直接影響著模型的準確程度,此處采用兩階段估計法對其進行求解,即先計算λij,之后由其估計值再計算出θ值。其主要計算步驟如下:
(1)對于λij的估計采用非線性規劃法求解,計算公式如下:
(3)
其中:
f(xi)=exp(λi0+∑λijxji)
式中:f(xi)為最大熵概率密度分布函數;mij為第i個變量第j階原點矩。
(2)利用得到的λij的值采用極大似然法[11]進行對θ值的估計,計算公式如下:
(4)
由此得到的λij和θ值代入式(1)即得到預報值、預報誤差以及出力誤差的三維聯合分布函數。
模型能否很好的反映樣本數據的分布特征,描述各變量之間的相關關系主要取決于對模型的檢驗評價,為此,本文分別采用P-P圖和均方根誤差法(RMSE)[12-14 ]對其進行檢驗,其中P-P圖是用來檢驗數據是否服從某一指定分布的一種方法,將變量的累積比例與相應的指定分布的累積比例點繪在一張圖中,并且計算回歸系數,其值越接近于1,說明指定分布的擬合效果越好。而RMSE是用來衡量理論分布頻率與經驗頻率之間的差距,其計算公式為:
(5)
式中:MSE為均方差;N為樣本容量;pc為由模型得出的理論頻率;p0為經驗頻率;i為樣本序號。
由式(5)可知,其值越小,模型的擬合效果越好。
模型的整體計算流程圖如圖1所示。

圖1 三維聯合分布模型計算流程圖Fig.1 The calculation flow chart of three-dimensional joint distribution model
本文選取某水電站A為例進行實例分析,分別利用該流域10年(2004-2013年)的每日12個時段的預報流量值以及實際值計算出預報出力值和實際出力值,進而得到與預報誤差相應的出力誤差,以此為數據樣本進行計算。
由于預報值、預報誤差與水電站出力誤差存在著復雜的數學聯系,出力誤差不僅和預報值大小以及預報誤差有關,還和水位因素有關。由于同一時期預報值的流量級別大致為一個范圍,在整個調度期水位范圍也大致相同,而不同時期水位范圍卻不盡相同,并且對于不同時期預報誤差與出力誤差也具有不同的相關關系,故采用時期的劃分來反映水位的特點。因此,本文分為枯水期、過渡期和汛期3個時期進行探討,按照當地流域資料劃分的月份如表1所示。

表1 時期劃分表Tab.1 Period division table
將劃分好的時期對預報值、預報誤差以及出力誤差進行相關性分析,進而確定三維聯合分布模型中生成元的類型,其中,弱相關采用AMH Copula[15]生成元:
φ(u)=ln{[1-θ(1-u)]/u}
正相關采用Gumbel-Hougaard Copula生成元:
φ(u)=-(lnu)θ
當呈現出負相關時,采用Frank Copula生成元:
φ(u)=-ln[(e-θu-1)/(e-θ-1)]
其不同時期三者之間的秩相關系數如表2所示。

表2 各變量間秩相關性系數Tab.2 The rank correlation coefficient among variables
由表2計算結果可知,在三個時期預報值與預報誤差、預報值與出力誤差的秩相關系數在0.1左右,均呈現出弱相關性。而預報誤差與出力誤差三個時期的秩相關系數在枯水期最大,為0.623;過渡期次之,為0.420;而在汛期二者的秩相關系數較弱,為0.251。通過分析上述變量間相關性可知,預報值與預報誤差,預報值與出力誤差呈現弱相關性,而預報誤差與出力誤差呈現出正相關性,故本文分別選取AMH Copula生成元以及Gumbel-Hougaard Copula生成元代入式(1)建立兩種形式的模型:模型1與模型2,并且對兩種模型進行檢驗,從而選出最優模型。
對模型1與模型2中的參數λi和θ利用式(2)~(4)進行估計,得到結果如表3~表6所示。
(1)清萬樹《詞律》中載《卜算子》調9體[注]柳永與張先的名為《卜算子》的詞作,實則是《卜算子慢》,且柳永之作為《卜算子慢》正體。[1]118-121

表3 各時期預報值拉格朗日因子Tab.3 Lagrangian multipliers of forecast value in each period

表4 各時期預報誤差拉格朗日因子Tab.4 Lagrangian multipliers of forecast error in each period

表5 各時期出力誤差拉格朗日因子Tab.5 Lagrangian multipliers of output error in each period
由計算出的參數λi和θ即可得到三維聯合分布模型1和模型2,為了驗證模型的可行性,采用χ2檢驗法檢驗模型的邊緣分布的擬合程度,判斷其是否通過檢驗。由于模型1與模型2所求邊緣分布相同,故只針對一種模型進行檢驗,顯著性水平取0.05,其檢驗結果如表7所示。

表6 兩種模型中的θ值Tab.6 Theta values of two models

表7 χ2值檢驗結果Tab.7 The test results of χ2
由表7可以看出,預報值、預報誤差以及出力誤差分布的χ2值分別為3.21、5.51和6.17,均小于檢驗臨界值11.07、9.48、9.48。通過檢驗要求,說明采用最大熵法擬合出的預報值分布、預報誤差分布以及出力誤差分布滿足檢驗要求,是可行的、合理的。
針對兩種形式的模型,為進一步選擇出擬合度較好的一種,采用P-P圖以及均方根誤差法來評價兩種模型優劣,從而確定最終的三維聯合分布模型形式。
本文利用樣本數據計算經驗頻率,再分別用模型1與模型2計算相應的理論頻率,其中經驗頻率計算公式如式(15)所示:
P(x1,x2,x3)=P(X1≤x1,X2≤x2,X3≤x3)=
(6)
式中:∑Ni為樣本系列中滿足條件的總的數據對數;N為樣本系列的數據對數。
得到的經驗頻率與模型1、模型2計算出相應的理論頻率,共同點繪到圖中,如圖2~圖7所示。

圖2 模型1枯水期P-P圖Fig.2 The P-P figure of model 1 in dry period

圖3 模型2枯水期P-P圖Fig.3 The P-P figure of model 2 in dry period

圖4 模型1過渡期P-P圖Fig.4 The P-P figure of model 1 in transitional period

圖5 模型2過渡期P-P圖Fig.5 The P-P figure of model 2 in transitional period

圖6 模型1汛期P-P圖Fig.6 The P-P figure of model 1 in flood period

圖7 模型2汛期P-P圖Fig.7 The P-P figure of model 2 in flood period
為了進一步評價兩個模型的優劣性,分別按照式(5)計算二者的均方根誤差值進行比較,所得結果如表8所示。

表8 兩種模型在不同時期的均方根誤差值Tab.8 The RSME values of two models in different periods
由表8可以看出,采用模型1函數計算得到的三個時期的均方根誤差分別為0.024 1、0.030 7、0.022 5,而模型2的均方根誤差分別為0.019 8、0.024 6、0.025 1,在枯水期以及過渡期模型2均方根誤差較小,模型2更優;在汛期模型1的均方根誤差較小,模型1更優。故在枯水期以及過渡期采用模型2來描述其三維聯合分布,在汛期采用模型1來描述。
由此便得到在不同時期預報值、預報誤差以及出力誤差三者的聯合分布,對于不同時期而言,其預報誤差以及出力誤差的二維聯合分布如圖8~圖10。

圖8 枯水期二維聯合分布圖Fig.8 Two-dimensional joint distribution figure in dry period

圖9 過渡期二維聯合分布圖Fig.9 Two-dimensional joint distribution figure in transitional period

圖10 汛期二維聯合分布圖Fig.10 Two-dimensional joint distribution figure in flood period
由圖8~圖10可以看出,在不同時期,預報誤差與出力誤差的二維聯合分布呈現出不同的形狀,說明預報誤差與出力誤差隨著時期的不同而呈現出不同的關系,分時期進行研究是有必要并且合理的。在此基礎上,運用所求出的三維聯合分布模型結合貝葉斯理論便可以求得在預報值與預報誤差一定時,水電站出力誤差的條件概率大小,使得水電站工作人員方便制定與修改發電計劃,更好地進行規劃與管理工作。
由貝葉斯理論可知:
P(A|B)=P(AB)/P(B)
(7)
式中:P(A|B)為在事件B發生的情況下發生A的條件概率;P(AB)為事件A和B同時發生的概率;P(B)為事件B發生的概率;此處事件A表示出力誤差;事件B表示預報值與預報誤差。
對于連續型的概率密度函數來說,一個點的概率值為0,一個范圍的概率才有意義,故首先對各個時期的預報值、預報誤差進行級別劃分,以此來進行概率計算,其分類結果如表9、表10所示。

表9 預報值分類結果 m3/sTab.9 Classification results of forecast value

表10 預報誤差與出力誤差級別劃分 %Tab.10 Level classification of forecast error and output error
由上述分類結果結合貝葉斯統計理論即可得到出力誤差的條件概率分布,其計算公式如式(17)所示:
(8)
式中:X、Y、Z分別表示為預報值,預報誤差以及出力誤差。
將三維聯合分布模型代入式(8)便可計算出各個時期的出力誤差條件概率大小,如:當預報值為6級,預報誤差為3級時,計算得到各個時期的出力誤差條件概率分布,其結果如表11所示。

表11 出力誤差條件概率分布Tab.11 Conditional probability distribution of output error
由表11可以看出,當預報誤差一定時,在不同時期呈現出不同的出力誤差概率值。其中枯水期出力誤差為3級時的概率為90%,在過渡期出力誤差為4級時的概率值最大,為93.7%,而對于汛期,出力誤差主要為6級,其概率為98.3%,造成這種不同主要是因為導致出力誤差的原因除了預報誤差外還與來流大小有著重要的關系,來流的大小主要從兩方面影響著出力誤差的結果,一是發電流量的多少,二是水庫水位的高低;在汛期,由于此級別下的來水已足夠大,導致水電站接近滿發,預報誤差對出力誤差的影響較弱,故此時對出力誤差基本沒有影響,出力誤差在6級處較為集中。此外,由于出力誤差除了和當前時段的預報誤差有關,還有當前水位與之前的預報誤差有著重要的關系,故導致同一時期相同預報值與預報誤差下出力誤差出現不同級別的概率情況。此結果可為水電站工作人員在預報值與實際值產生一定偏差時修改發電計劃提供科學的參考依據。
針對水電站運行中由于來水偏差產生出力誤差而被迫修改發電計劃的情況,本文對入庫流量預報值、預報誤差及其相應的出力誤差三者之間的相關性進行了分析,并且基于最大熵理論與Copula函數理論建立三維聯合分布模型,求解出三者之間的聯合分布。在此基礎上,利用該模型結合貝葉斯理論求解出不同時期不同預報級別與預報誤差情況下水電站出力誤差的條件概率分布。結果表明,在枯水期和過渡期,預報誤差與出力誤差相關性較強,并且成正相關,汛期各變量之間的相關性較弱;所建立的三維聯合分布模型與樣本數據擬合較好,能夠充分反映樣本數據的分布特征,各項指標均滿足檢驗要求,用該模型表述預報值、預報誤差以及出力誤差的聯合分布是可行的、合理的;此外,由該模型求解出的水電站出力誤差條件概率不同時期呈現出不同的分布狀態,符合實際情況,能夠有效反映不同時期不同預報值與預報誤差下其出力誤差的特點。該結果可為水電站工作人員由于預報來水偏差而及時修改發電計劃提供重要的指導信息,為水電站管理工作提供科學的參考依據。
由于影響水電站出力誤差的因素有很多,除了預報值大小,預報誤差以外,還有輸電線效率,意外事故,人為因素等相關內容,本文僅考慮了一般情況,即預報值、預報誤差對出力誤差的影響,在下一階段,將考慮對上述其他因素進行探討,進一步完善對出力誤差分布的研究。
□
[1] 宋萌勃,岳延兵,陳吉琴,等.水庫調度與管理[M]. 鄭州:黃河水利出版社,2013.
[2] 費永法. 一種計算洪水條件概率的方法[J]. 水文,1989,(1):18-23.
[3] 戴昌軍,梁忠民. 多維聯合分布計算方法及其在水文中的應用[J]. 水利學報,2006,(2):160-165.
[4] 王文圣,丁 晶.水文學確定性和不確定性方法的耦合初探[J].水文,2014,(2):3-7.
[5] 鄭紅星,劉昌明. 南水北調東中兩線不同水文區降水豐枯遭遇性分析[J]. 地理學報,2000,(5):523-532.
[6] 李憲東. 基于最大熵原理的確定概率分布的方法研究[D]. 北京:華北電力大學,2008.
[7] 刁艷芳,王本德,劉 冀.基于最大熵原理方法的洪水預報誤差分布研究[J].水利學報,2007,(5):591-595.
[8] 方 彬,郭生練,肖 義,等.年最大洪水兩變量聯合分布研究[J]. 水科學進展,2008,(4):505-511.
[9] 張 翔,冉啟香,夏 軍,等. 基于Copula函數的水量水質聯合分布函數[J].水利學報,2011,(4):483-489.
[10] CHOW DHARY H, ESCOBAR L A, SINGH V P. Identification of suitable copulas for bivariate frequency analysis of flood peak and volume data[J]. Hydrology Reasearch, 2011,42(2/3):193-215.
[11] 謝 華,羅 強,黃介生. 基于三維copula函數的多水文區豐枯遭遇分析[J].水科學進展,2012,(2):186-193.
[12] 康 玲,何小聰. 南水北調中線降水豐枯遭遇風險分析[J]. 水科學進展,2011,(1):44-50.
[13] Vandenberghe S,Verhoest N E C,Onof C, et al. A comparative copula-based bivariate frequency analysis of observed and simulated storm events: A case study on Barlett-Lewis modeled rainfall[J]. Water Resource Research, 2011,47:W07529.
[14] SONG S B, Vijay P Singh. Frequency analysis of droughts using the Plackett copula and parameter estimation by genetic algotirhm[J]. Stochastic Enviromental Research and Risk Assessment,2010,doi:10.1007/s00477-0101-0364-5.
[15] 李 計,李 毅,宋松柏,等.基于Copulas函數的多維干旱變量聯合分布[J].自然資源學報,2013,(2):312-320.