程貝麗,周菊玲

(新疆師范大學(xué) 數(shù)學(xué)科學(xué)學(xué)院,新疆 烏魯木齊 830017)

0 引言

對(duì)數(shù)伽瑪分布是一種重要的壽命分布,廣泛應(yīng)用在生物學(xué)、金融及新型農(nóng)村合作醫(yī)療住院費(fèi)用等可靠性領(lǐng)域[1-4].由于對(duì)數(shù)伽瑪分布具有向右偏斜的特點(diǎn),所以在非壽險(xiǎn)精算中也得到了一定的應(yīng)用.Zhang等[1]從生物學(xué)角度研究了對(duì)數(shù)伽瑪分布的模型;Amini等[5]研究了對(duì)數(shù)伽瑪分布族的性質(zhì);王成元等[6-8]研究了不同損失函數(shù)下對(duì)數(shù)伽瑪分布尺度參數(shù)的貝葉斯估計(jì),并對(duì)不同的貝葉斯估計(jì)做出比較;金秀巖[9-10]分別探討了對(duì)數(shù)伽瑪分布和負(fù)對(duì)數(shù)伽瑪分布條件概率的封閉性、Pearson-χ2距離和漸進(jìn)性;熊福生[11]研究了對(duì)數(shù)伽瑪分布和負(fù)對(duì)數(shù)伽瑪分布的再生性,得到了一些與再生性、準(zhǔn)再生性相關(guān)的結(jié)果.在近幾年的研究中,變點(diǎn)問(wèn)題一直是統(tǒng)計(jì)學(xué)中的熱點(diǎn)問(wèn)題.變點(diǎn)問(wèn)題一般的研究方法有最小二乘法、貝葉斯方法、極大似然方法、累次計(jì)數(shù)法和非參數(shù)法.陳希孺[12]結(jié)合概率變點(diǎn)問(wèn)題,介紹了極大似然音樂(lè)方法、累計(jì)次數(shù)法和貝葉斯方法;黃月蘭[13]分別用極大似然方法、貝葉斯方法和最小二乘法討論了單參數(shù)指數(shù)分布的變點(diǎn)問(wèn)題.沙雪云等[14]用極大似然方法和貝葉斯方法討論了單參數(shù)Parato分布的變點(diǎn)問(wèn)題;陳杰[15]分別在非貝葉斯場(chǎng)景下和貝葉斯場(chǎng)景下對(duì)變點(diǎn)進(jìn)行檢測(cè).但是，用不同的方法研究對(duì)數(shù)伽瑪分布變點(diǎn)問(wèn)題的文章還比較少，本文用極大似然估計(jì)和貝葉斯估計(jì)研究單變點(diǎn)對(duì)數(shù)伽瑪分布的參數(shù)估計(jì)問(wèn)題.

1 對(duì)數(shù)伽瑪分布單變點(diǎn)模型

下面分別用極大似然估計(jì)方法和貝葉斯估計(jì)方法對(duì)單變點(diǎn)對(duì)數(shù)伽瑪分布的參數(shù)進(jìn)行估計(jì).

2 對(duì)數(shù)伽瑪分布的極大似然估計(jì)

設(shè)隨機(jī)變量X1,X2,…,Xn為來(lái)自對(duì)數(shù)伽瑪分布模型中的一個(gè)樣本,則對(duì)其n個(gè)獨(dú)立觀測(cè)值x1,x2,…,xn得到的似然函數(shù)為

(1)

令γ=(θ1,θ2,k),則單變點(diǎn)對(duì)數(shù)伽瑪分布的似然函數(shù)為

(2)

固定k,對(duì)θ1,θ2求極大,得

(3)

3 對(duì)數(shù)伽瑪分布的貝葉斯估計(jì)

在統(tǒng)計(jì)推斷中,樣本的狀態(tài)固然重要,但歷史的經(jīng)驗(yàn)也同樣舉足輕重.故為了提高統(tǒng)計(jì)推斷的質(zhì)量,應(yīng)注意先驗(yàn)信息對(duì)整個(gè)統(tǒng)計(jì)推斷的影響,沒(méi)有先驗(yàn)分布或取不恰當(dāng)?shù)南闰?yàn)分布都可能使統(tǒng)計(jì)推斷的結(jié)果產(chǎn)生誤差.

首先確定各參數(shù)的先驗(yàn)分布:

1)變點(diǎn)k的先驗(yàn)分布為無(wú)信息先驗(yàn)分布:π(k)∝1；

2)尺度參數(shù)θ1,θ2的先驗(yàn)分布分別取無(wú)信息先驗(yàn)分布和共軛先驗(yàn)分布.

討論尺度參數(shù)θ1,θ2的先驗(yàn)分布取無(wú)信息分布時(shí)參數(shù)的貝葉斯估計(jì),假定θ1,θ2,k之間相互獨(dú)立.

引理1 設(shè)隨機(jī)變量X服從對(duì)數(shù)伽瑪分布,參數(shù)θi的先驗(yàn)分布π(θi)為無(wú)信息先驗(yàn)分布,X1,X2,…,Xn為服從對(duì)數(shù)伽瑪分布的一個(gè)樣本,記X=(X1,X2,…,Xn),T=(t1,t2),則θi的后驗(yàn)分布服從Ga(nα,T).

證由式(2)可知,x1,x2,…,xn是樣本X1,X2,…,Xn的n個(gè)獨(dú)立觀測(cè)值的似然函數(shù),又已知參數(shù)θi的先驗(yàn)分布為無(wú)信息先驗(yàn)分布,即π(θi)∝1.則參數(shù)θi的后驗(yàn)概率密度函數(shù)為

于是θi|x～Ga(nα,ti),故θi|X～Ga(nα,T).

則當(dāng)參數(shù)θi的先驗(yàn)分布為無(wú)信息先驗(yàn)分布時(shí),對(duì)數(shù)伽瑪分布單變點(diǎn)模型的似然函數(shù)為

(4)

得到各參數(shù)的滿(mǎn)條件分布如下：

討論尺度參數(shù)θ1,θ2的共軛先驗(yàn)分布為伽瑪分布時(shí)參數(shù)的貝葉斯估計(jì),假定θ1,θ2,k之間相互獨(dú)立.

引理2[3]設(shè)隨機(jī)變量X服從對(duì)數(shù)伽瑪分布,參數(shù)θi的先驗(yàn)分布π(θi)為伽瑪分布Ga(ai,bi),X1,X2,…,Xn是服從對(duì)數(shù)伽瑪分布的一個(gè)樣本,記X=(X1,X2,…,Xn),T=(t1,t2),則θi的后驗(yàn)分布服從Ga(nα+ai,T+bi).

則參數(shù)θi的后驗(yàn)概率密度函數(shù)為

于是θi|x～Ga(nα+ai,bi+ti),故θi|X～Ga(nα+ai,bi+T).

則當(dāng)θi的先驗(yàn)分布為伽瑪分布時(shí),對(duì)數(shù)伽瑪分布單變點(diǎn)模型的似然函數(shù)為

(5)

得到各參數(shù)的滿(mǎn)條件分布如下：

尺度參數(shù)θ1,θ2滿(mǎn)條件分布的形式較為簡(jiǎn)單,可直接用Gibbs抽樣,而變點(diǎn)k的滿(mǎn)條件分布較為復(fù)雜,不再適合Gibbs抽樣方法,則利用MH算法進(jìn)行抽樣.

MCMC算法可分為如下幾個(gè)具體步驟.

為了確保抽樣的收斂性,在模擬過(guò)程中,先進(jìn)行10 000次的預(yù)迭代,再進(jìn)行20 000次的Gibbs迭代,即從10 001次開(kāi)始到30 000次結(jié)束.假設(shè)Gibbs樣本的容量為M,如果第B次以后逐漸收斂,則把后M-B個(gè)樣本的均值作為各未知參數(shù)的估計(jì)值.

4 隨機(jī)模擬

4.1 極大似然估計(jì)

根據(jù)對(duì)數(shù)伽瑪分布的概率密度函數(shù)和極大似然估計(jì)式(3),利用MATLAB軟件對(duì)未知參數(shù)θ1,θ2,k進(jìn)行估計(jì),結(jié)果如表1所示.

表1 參數(shù)θ1,θ2,k的極大似然估計(jì)

模擬結(jié)果顯示,極大似然方法能夠估計(jì)未知參數(shù)的值,但參數(shù)θ1,θ2,k的極大似然估計(jì)值與真值相差較大,精度不高.

4.2 貝葉斯估計(jì)

取超參數(shù)ai=bi=1,i=1,2,取B=10 000,M=30 000,根據(jù)各參數(shù)的滿(mǎn)條件分布進(jìn)行MCMC模擬,將20 000個(gè)樣本的均值作為各未知參數(shù)的估計(jì)值.

1)當(dāng)先驗(yàn)分布取無(wú)信息先驗(yàn)分布時(shí),參數(shù)θ1,θ2,k的參數(shù)估計(jì)如表2所示.

表2 無(wú)信息先驗(yàn)分布情況下參數(shù)θ1,θ2,k的參數(shù)估計(jì)

2)當(dāng)先驗(yàn)分布取共軛先驗(yàn)分布時(shí)，參數(shù)θ1,θ2,k的參數(shù)估計(jì)如表3所示.

表3 共軛先驗(yàn)分布情況下參數(shù)θ1,θ2,k的參數(shù)估計(jì)

圖1、圖3是變點(diǎn)k抽樣迭代圖,是判斷Gibbs樣本是否收斂的重要依據(jù).用MCMC方法模擬時(shí),會(huì)同時(shí)產(chǎn)生兩條迭代鏈,若這兩條迭代鏈逐漸穩(wěn)定并趨于重合,則說(shuō)明Gibbs抽樣是收斂的,變點(diǎn)k的兩條迭代鏈如圖2、圖4所示.

圖1 變點(diǎn)k的抽樣迭代圖 圖2 變點(diǎn)k的兩條迭代鏈

圖3 變點(diǎn)k的抽樣迭代圖 圖4 變點(diǎn)k的兩條迭代鏈

下面對(duì)模擬的結(jié)果進(jìn)行分析,從表1、表2中可以看到θ1,θ2,k這3個(gè)參數(shù)的均值與真值都較為接近,且MC誤差都小于2%,即各參數(shù)的貝葉斯估計(jì)精度較高;各參數(shù)的標(biāo)準(zhǔn)差也比較小,說(shuō)明數(shù)據(jù)比較穩(wěn)定;各參數(shù)置信水平0.95的區(qū)間[2.5%分位數(shù),97.%分位數(shù)]的長(zhǎng)度也較短,即得到的區(qū)間估計(jì)效果較好.由圖1和圖3可以看出Gibbs抽樣基本在變點(diǎn)附近波動(dòng),呈現(xiàn)出一定的規(guī)律性.其次,觀察圖2和圖4,可以發(fā)現(xiàn)產(chǎn)生的兩條馬爾科夫鏈穩(wěn)定且趨于重合,說(shuō)明收斂性較好.綜上所述,對(duì)數(shù)伽瑪分布變點(diǎn)問(wèn)題可以用MCMC方法實(shí)現(xiàn),且在不同先驗(yàn)分布下,估計(jì)值的精度不同.

5 總結(jié)

本文主要研究對(duì)數(shù)伽瑪分布尺度參數(shù)存在變點(diǎn)時(shí)的參數(shù)估計(jì)問(wèn)題.隨機(jī)模擬的結(jié)果表明:兩種方法都能夠有效地對(duì)未知參數(shù)進(jìn)行估計(jì).在相同條件下,貝葉斯估計(jì)所需要的時(shí)間遠(yuǎn)遠(yuǎn)大于極大似然估計(jì),但是相比于極大似然估計(jì),貝葉斯估計(jì)的值更穩(wěn)定,也更接近真實(shí)值.通過(guò)比較可知,共軛先驗(yàn)分布下的貝葉斯估計(jì)值優(yōu)于無(wú)信息先驗(yàn)下的貝葉斯估計(jì)值.這說(shuō)明了選取恰當(dāng)?shù)南闰?yàn)分布可以提高統(tǒng)計(jì)推斷的質(zhì)量.因此,基于本文中單變點(diǎn)對(duì)數(shù)伽瑪分布的模型,先驗(yàn)分布取共軛先驗(yàn)分布時(shí),未知參數(shù)的貝葉斯估計(jì)是最優(yōu)的.