基于混合正態(tài)分布的降水量分布特征估計(jì) *

李嬌華

（廣西師范大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，廣西 桂林 541004）

0 引言

近年來，全球性氣候變化對(duì)環(huán)境、糧食生產(chǎn)的影響受到世界各國政府和專家學(xué)者越來越多的關(guān)注。我國是農(nóng)業(yè)大國，國民經(jīng)濟(jì)的基礎(chǔ)是農(nóng)業(yè)，農(nóng)業(yè)是受天氣和氣候影響最大的領(lǐng)域。氣候的變化將會(huì)直接影響我國的糧食安全及其可持續(xù)發(fā)展。因此，研究降水量及其分布規(guī)律不僅有重要的氣候?qū)W意義，而且對(duì)農(nóng)業(yè)生產(chǎn)、水資源管理以及防災(zāi)、減災(zāi)有著重要的現(xiàn)實(shí)意義。

在應(yīng)用中，單分布函數(shù)模型成為降水量概率密度估計(jì)中最常用的方法。正態(tài)分布模型、對(duì)數(shù)正態(tài)分布模型、Gamma分布模型是單分布函數(shù)模型中最常用的分布模型。在某些情況下，使用這些模型分析往往能取得不錯(cuò)的效果。吳慧在文獻(xiàn)［1］中使用海南省18個(gè)測(cè)站收集1966-2003年的資料，結(jié)果得出在0.01的信度檢驗(yàn)下春、夏、秋季降水?dāng)?shù)據(jù)基本符合正態(tài)分布。劉慧等人在文獻(xiàn)［2］中使用宜賓市1953-2010年的春季降水量數(shù)據(jù)資料，降水量分布分別采用正態(tài)分布、Gamma分布、對(duì)數(shù)正態(tài)分布和P-Ⅲ型分布對(duì)春季降水量分布進(jìn)行研究。通過檢驗(yàn)發(fā)現(xiàn)，對(duì)數(shù)正態(tài)分布對(duì)宜賓市春季降水量的擬合效果最好。Mooley在文獻(xiàn)［3］中使用Gamma分布對(duì)亞洲夏季月降水量進(jìn)行模擬，結(jié)果表明該模型適合亞洲夏季季風(fēng)的月降雨量。然而不少研究表明，氣候變化較大的地區(qū)，使用單分布模型并不能很好地估計(jì)某地區(qū)的降水量分布，需要新的模型來估計(jì)。Li等人在文獻(xiàn)［4］中使用指數(shù)分布、Gamma分布、Weibull分布等單分布模型以及混合概率分布擬合日降水量，考查其模擬中國黃土高原地區(qū)日降水分布特征的效果。

曹杰等人在文獻(xiàn)［5］中研究發(fā)現(xiàn)，我國降水量分布具有一定的區(qū)域性。事實(shí)上，降水量分布不會(huì)很好地服從正態(tài)分布。如果對(duì)降水量的研究均是假定服從正態(tài)分布，其預(yù)報(bào)區(qū)間估計(jì)和顯著性檢驗(yàn)結(jié)果將會(huì)產(chǎn)生偏差，導(dǎo)致估計(jì)精度降低。下面我們以南寧市2000-2019年的年降水量為例來分析其概率分布，原始數(shù)據(jù)來源于國家氣象科學(xué)數(shù)據(jù)共享服務(wù)平臺(tái)-中國地面氣候資料日值數(shù)據(jù)集（V3.0）。降水量的直方圖和混合正態(tài)分布擬合曲線（實(shí)線）如圖1所示。圖1的橫坐標(biāo)表示南寧市2000-2019年的年降水量（單位：100 mm），縱坐標(biāo)表示出現(xiàn)該降水量的頻率。由圖1看出，南寧市2000-2019年降水量數(shù)據(jù)出現(xiàn)雙峰狀，該數(shù)據(jù)顯然不是對(duì)稱型的正態(tài)分布，應(yīng)該為混合正態(tài)分布。為進(jìn)一步分析此情況，計(jì)算其相關(guān)數(shù)字特征，得到其均值為14.56，方差為3.95，偏度為0.50<0，峰度為-0.52<3，結(jié)果說明與正態(tài)分布的數(shù)字特征有一定偏差。Shapiro-Wilk正態(tài)性檢驗(yàn)顯示，顯著性概率p值為5.529×10-7<0.05，因此，可以說明該地區(qū)年降水量數(shù)據(jù)不服從正態(tài)分布。

針對(duì)南寧市2000-2019年的年降水量數(shù)據(jù)，我們使用混合正態(tài)分布和正態(tài)分布進(jìn)行擬合比較，其擬合分布的密度曲線見圖1和圖2中的虛線。結(jié)合兩個(gè)圖，可以直觀看出，正態(tài)分布的刻畫并沒有很好體現(xiàn)出降水量數(shù)據(jù)的分布，相比之下，混合正態(tài)分布是對(duì)降水量數(shù)據(jù)較好的一個(gè)刻畫。

圖1 南寧市2000-2019年的年降水量直方圖和混合正態(tài)分布擬合曲線

圖2 南寧市2000-2019年的年降水量直方圖和正態(tài)分布擬合曲線

本文以南寧市為例，借助混合正態(tài)分布對(duì)降水量分布的總體參數(shù)進(jìn)行估計(jì)，從而研究其分布規(guī)律。

1 有限混合正態(tài)分布及其參數(shù)估計(jì)

混合模型在計(jì)量經(jīng)濟(jì)學(xué)和社會(huì)科學(xué)中得到廣泛的應(yīng)用，其理論也得到廣泛的研究。其中，在混合模型中應(yīng)用較多的是混合正態(tài)分布模型。更多的混合正態(tài)分布模型理論研究可以參考文獻(xiàn)［6-9］。首先給出混合正態(tài)分布的概率分布模型：

其中，αk是系數(shù)是正態(tài)分布的密度函數(shù)，，則

稱為第K個(gè)分模型。

下面簡單介紹混合正態(tài)分布模型參數(shù)估計(jì)的EM算法，具體推導(dǎo)過程可參見文獻(xiàn)［9］。假設(shè)觀測(cè)數(shù)據(jù)y1，y2，…，yN由混合正態(tài)分布生成，

其中，θ=（α1，α2，…，αK；θ1，θ2，…，θK)，我們用EM算法來估計(jì)混合正態(tài)分布模型的參數(shù)θ。經(jīng)過計(jì)算，可以給出第i次迭代的參數(shù)為：

由此，估計(jì)混合正態(tài)分布參數(shù)的EM算法總結(jié)如下：

（1）給出參數(shù)的初始值進(jìn)行迭代；

（2）E步：根據(jù)當(dāng)前模型參數(shù)，計(jì)算分模型k對(duì)觀測(cè)數(shù)據(jù)yj的響應(yīng)度

（3）M步：計(jì)算第i步迭代的模型參數(shù)：

（4）重復(fù)第（2）和（3）步驟，直到對(duì)數(shù)似然函數(shù)值不再有明顯變化為止，最終可以得到模型的參數(shù)估計(jì)為

加入的廢液偏少，BT值會(huì)偏高，溶液呈乳白色，易形成“堿式絮狀物”，液相與固相比重差縮小后造成鎘渣漂移、出口含固量升高、首槽鎘渣品位低、除鎘反應(yīng)器含鎘梯度不明顯、迫使鋅粉更換加快等一系列化“惡性循環(huán)”。廢液加入過多，BT值偏低，消耗鋅粉量增加，槽內(nèi)產(chǎn)生氣泡，同樣影響沸騰層的穩(wěn)定。通過長時(shí)間試驗(yàn)發(fā)現(xiàn)：除鎘 BT值適宜控制在0.7～1.2之間，宜低不宜高。

本文主要探討k=2的情形，即兩個(gè)混合正態(tài)分布，計(jì)算該分布的均值和方差為：

其中α1+α2=1，假設(shè)mα表示兩個(gè)混合正態(tài)分布的α分位點(diǎn)，可以得到公式

其中Φ（·)表示標(biāo)準(zhǔn)正態(tài)分布的分布函數(shù)，按照此公式可以通過隱函數(shù)算出mα。利用α分位點(diǎn)和混合正態(tài)分布的參數(shù)可以掌握樣本的重要信息，并推出其潛在規(guī)律以及變化趨勢(shì)。利用前面收集到的樣本數(shù)據(jù)并使用EM迭代算法，可以獲得參數(shù)估計(jì)

本文接下來將通過模擬對(duì)模型參數(shù)及其變化規(guī)律做簡單探討。

2 數(shù)值模擬

基于混合正態(tài)分布模型的參數(shù)迭代公式，并借助R軟件進(jìn)行數(shù)值模擬實(shí)驗(yàn)，對(duì)EM算法求解混合正態(tài)分布模型中參數(shù)估計(jì)問題的可行性和準(zhǔn)確性進(jìn)行驗(yàn)證。并考慮當(dāng)數(shù)據(jù)來自混合正態(tài)分布模型時(shí)，使用混合正態(tài)分布估計(jì)的參數(shù)和使用正態(tài)分布估計(jì)的參數(shù)結(jié)果進(jìn)行比較。假設(shè)數(shù)據(jù)y1，…，yn是來自總體分布

我們考慮以下幾種情形，主要考慮混合正態(tài)分布的參數(shù)變化情形，包括混合比例變化、均值變化、方差變化等方面。情形1和情形2主要觀察混合比例的變化；情形2和情形3主要觀察均值的變化；情形3和情形4主要觀察方差的變化。

以上4種情形分別產(chǎn)生20、50、100、300和500個(gè)樣本容量（n)進(jìn)行模擬重復(fù)10 000次。按照以上4種情形，分別計(jì)算相應(yīng)的參數(shù)估計(jì)。為了方便比較，我們主要 考慮不同 樣本下 的參數(shù)估 計(jì)μ?，σ?2，m?α及其 對(duì)應(yīng)標(biāo)準(zhǔn)差（即表格中括號(hào)內(nèi)的數(shù)字）；整體均值95%對(duì)稱置信區(qū)間（即模擬數(shù)據(jù)按照從小到大排序后首尾兩端各自去掉2.5%）以及對(duì)應(yīng)的區(qū)間長度。模擬實(shí)驗(yàn)的主要估計(jì)結(jié)果如表1和表2所示。

根據(jù)模擬結(jié)果表1和表2可以得到以下結(jié)論：

（1）從表1可以看出，混合正態(tài)分布的總體均值和總體方差估計(jì)具有相合性。隨著樣本容量不斷增加，估計(jì)值逐漸趨于真實(shí)值。估計(jì)值對(duì)應(yīng)的標(biāo)準(zhǔn)差也逐漸減小，即數(shù)據(jù)分布越穩(wěn)定。

（2）從表1中的情形1和情形2來看，混合比例的變化會(huì)影響總體均值和方差的數(shù)值改變。一般來看，均值增減的趨勢(shì)和方差增減的趨勢(shì)是相反的，即均值減小時(shí)，方差則增加。分位數(shù)的變化方向和均值變化方向是一致的。從情形3和情形4來看，方差的變化對(duì)均值和0.05分位數(shù)的影響不大。

表1 各情形下的參數(shù)估計(jì)

（3）從表2可以看出，當(dāng)真實(shí)總體服從混合正態(tài)分布時(shí)，在95%的置信水平下，總體均值的置信區(qū)間估計(jì)的區(qū)間長度小于正態(tài)分布估計(jì)的區(qū)間長度。

表2 各情形下整體均值的置信區(qū)間估計(jì)（置信水平為95%）

3 實(shí)證分析

在本節(jié)中，我們選取前面所述的南寧市2000-2019年的年降水量數(shù)據(jù)作為實(shí)例對(duì)前面介紹的方法進(jìn)行分析。年降水量數(shù)據(jù)共有240個(gè)數(shù)據(jù)。根據(jù)圖1可以看出，南寧市2000-2019年的年降水量數(shù)據(jù)直方圖呈雙峰狀，第一個(gè)峰值接近13，第二個(gè)峰值接近15，因此主要考慮兩個(gè)混合正態(tài)分布即k=2的情況來分析該降水量。

使用R軟件可以得到表3參數(shù)估計(jì)的結(jié)果，根據(jù)參數(shù)估計(jì)結(jié)果可知，年降水量數(shù)據(jù)服從混合正態(tài)分布，即：0.34N（12.87，0.22)+0.66N（15.43，3.65)。經(jīng)過計(jì)算，該數(shù)據(jù)的總體均值為14.55，總體方差為3.94，0.05分位數(shù)為12.17。其中，混合正態(tài)分布的0.05分位數(shù)為12.15。若降水量數(shù)據(jù)來自均值為14.55，方差為3.94的正態(tài)分布，則0.05分位數(shù)為11.28。因此，在混合分布假設(shè)下，混合正態(tài)分布的0.05分位數(shù)更接近樣本的0.05分位數(shù)。

表3 參數(shù)估計(jì)結(jié)果

整體均值置信區(qū)間估計(jì)結(jié)果如表4所示。從表4可以發(fā)現(xiàn)，在95%的置信水平下，降水量數(shù)據(jù)服從混合正態(tài)分布總體均值的置信區(qū)間估計(jì)的區(qū)間長度小于正態(tài)分布估計(jì)的區(qū)間長度，說明混合正態(tài)分布總體均值的置信區(qū)間估計(jì)更精確。因此，選擇混合正態(tài)分布擬合年降水量數(shù)據(jù)對(duì)降水徑流模擬、水庫調(diào)度及水資源規(guī)劃配置、農(nóng)業(yè)規(guī)劃等方面的研究均具有重要意義。

表4 降水量數(shù)據(jù)的整體均值置信區(qū)間估計(jì)（置信水平為95%）

4 結(jié)論

本文主要探討基于混合正態(tài)分布模型的降水量分布特征，利用南寧市2000-2019年的年降水量數(shù)據(jù)并通過模擬實(shí)驗(yàn)使用EM算法對(duì)分布參數(shù)進(jìn)行估計(jì)，找到最優(yōu)的參數(shù)進(jìn)行擬合。結(jié)果發(fā)現(xiàn)，在年降水量數(shù)據(jù)擬合時(shí)，混合正態(tài)分布模型的擬合效果要優(yōu)于正態(tài)分布模型的擬合效果。由此可見，混合正態(tài)分布可以應(yīng)用于呈多峰形的降水量數(shù)據(jù)分析，在氣候研究中具有很好的效果。