Lomax分布形狀參數(shù)變點的貝葉斯估計

沙雪云， 周菊玲， 董翠玲

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

0 引言

變點問題是近幾年統(tǒng)計學的熱點研究問題，變點模型則是研究變點問題的一種非常重要的統(tǒng)計模型，其應用廣泛，研究方法多樣. 常用的研究方法有貝葉斯方法、極大似然法和最小二乘法等，其中貝葉斯方法在解決變點問題上綜合了先驗信息以及樣本信息，使得判斷更為準確. MCMC算法是一種高效的貝葉斯方法，將Gibbs抽樣與M-H抽樣相結合的算法，根據參數(shù)的滿條件分布形式來選取Gibbs抽樣和M-H抽樣，進而得到參數(shù)的Gibbs樣本，最終把Gibbs樣本的均值作為各參數(shù)的貝葉斯估計. Lomax分布是一種非常重要的壽命分布，具有單調的失效率，所以在壽命試驗的數(shù)據處理中起著至關重要的作用. 姚惠等人分別研究了基于沒有任何數(shù)據缺失的情況下，Lomax分布在熵損失函數(shù)、Linex損失函數(shù)等不同的損失函數(shù)下參數(shù)的bayes估計[1-4]. 龍兵等人針對在數(shù)據不完整的情況下對Lomax分布進行各參數(shù)估計[5-8]. 岑泰林討論了在完全數(shù)據和隨機刪失數(shù)據不同情況下Lomax分布的參數(shù)估計問題[9]. 但是，目前對Lomax分布單變點問題的研究較少. 本文將研究在尺度參數(shù)已知的情況下，對Lomax分布的形狀參數(shù)及變點位置進行貝葉斯估計，并運用MCMC算法進行隨機模擬，結果表明，各參數(shù)的估計值與真實值的MC誤差較小，說明各參數(shù)估計值在較高水平上是行之有效的.

1 Lomax分布及變點模型

設隨機變量X服從參數(shù)為λ,θ的Lomax分布，則分布函數(shù)和密度函數(shù)分別為

其中尺度參數(shù)λ>0，形狀參數(shù)θ>0.

假設樣本yi(i=1,2,…,n)值服從Lomax分布，若在序列{y1,y2,…,yn}存在變點，則在該序列中存在一個時間點，在該時間點的前后樣本值yi所服從的Lomax分布的尺度參數(shù)θ也將會發(fā)生變化，即在這一時間點之前的序列服從參數(shù)θ1為的Lomax分布，在這一時間點之后的序列服從參數(shù)為θ2的Lomax分布，稱該時間點為序列中的一個變點，當在序列中存在兩個及以上的變點時，則稱此模型為多變點模型. 含有k個未知變點的Lomax分布的模型為

(1)

其中m1,m2,…,mk是所要求得未知變點參數(shù)，θ1,θ2,…,θk分別為不含變點區(qū)間的尺度參數(shù). 針對含有多個未知變點的情況，可以采用二分分段法和貝葉斯因子進行逐段尋找，其中用貝葉斯因子來判斷序列是否存在單個變點可以借助于Kass和Raftery在1995年提出的貝葉斯因子模型選擇的一般準則[10]. 下面給出用二分分段法及貝葉斯因子解決多變點問題的具體思路：

在處理多變點問題時，首先需要確定變點是否存在，若存在變點，則要在給定變點個數(shù)的情況下求出變點的所在位置，借助二分分段法及貝葉斯因子先來確定變點的個數(shù)，然后判斷在一個序列集中是否存在一個變點及這個變點所在的位置.

第1步：在整個序列集S中，利用貝葉斯因子來檢測是否存在單個變點. 若貝葉斯因子檢測得不存在變點，就表示整個序列集不存在變點，則終止程序；否則就能夠檢測出這個序列里的第一個變點；

第2步：基于檢測出的第一個變點，將整個序列集S分為兩個子序列，對每個子序列，用第1步分別尋找出一個變點，一直持續(xù)這個過程直到在每個子序列中找不到變點為止.

利用二分分段法，只需要檢測并估計出單個變點的位置，接下來，利用貝葉斯方法來討論Lomax分布的單變點問題.

假設變點個數(shù)有且只有一個，Lomax分布形狀參數(shù)單變點模型為

(2)

其中x>0,θ1>0,θ2>0,λ=c(c為常數(shù))，且m,θ1,θ2均為未知. 當θ1≠θ2時，m=2,3,…,n-1就是所要研究的變點，同時還需要對變點m以及參數(shù)θ1,θ2進行貝葉斯估計.

2 單變點模型參數(shù)的貝葉斯估計

當θ1≠θ2時，設m為變點，記α=(m,θ1,θ2)，則此變點問題的似然函數(shù)為

1)θ1,θ2通過用Fisher信息陣確定無信息先驗分布.

首先, 寫出樣本似然函數(shù)為

通過樣本似然函數(shù)可得到與其相關的樣本的信息矩陣,

進而得到θ1,θ2的無信息先驗矩陣為,

最后θ1,θ2的無信息先驗分布為

由貝葉斯公式求得聯(lián)合后驗分布為,

π(α|x)∝L(α|x)π(α)=

各參數(shù)滿條件分布為,

2) 取θ1,θ2的共軛先驗分布為Ga(bi,ci)，即

且m,θ1,θ2相互獨立，由貝葉斯公式得聯(lián)合后驗分布為,

π(α|x)∝L((α|x)π(m)π(θ1)π(θ2))∝

可得各參數(shù)滿條件分布為,

3 隨機模擬

下面進行隨機模擬估計，通過上述分別得到了在無信息先驗分布和共軛先驗分布下的Lomax分布變點位置及各參數(shù)的滿條件分布，可以利用MCMC算法對變點位置及各參數(shù)進行隨機模擬估計，由于參數(shù)θ1,θ2,m的滿條件分布都較為復雜，用Gibbs方法對樣本進行抽樣較為困難，所以利用M-H方法對各參數(shù)θ1,θ2,m滿條件分布進行抽樣.

取隨機產生n=183個數(shù)據作為樣本，參數(shù)(m,θ1,θ2,λ)的真實值取(100,3,18,5)，即，

(3)

利用得到的各參數(shù)滿條件分布，運用MATLAB軟件來進行MCMC模擬，在模擬過程中檢驗模型參數(shù)的收斂性，通過輸入多組初始值，形成多層Markov鏈，再用MATLAB軟件對參數(shù)進行多層Markov鏈迭代分析，當參數(shù)m,θ1,θ2收斂時，Markov鏈迭代圖將趨于重合. 為確保參數(shù)的收斂性，先進行10 000次的預迭代的基礎上在進行20 000次迭代，即從第10 001次開始至30 000次開始迭代，可得MATLAB的運行結果，如表1所示.

當θ1,θ2的先驗分布為無信息先驗分布時:

表1 無信息先驗分布下，參數(shù)m,θ1,θ2的貝葉斯估計

當θ1,θ2取共軛先驗分布時，θ1～Ga(15,3.2)，θ2～Ga(27,1.2)，

表2 共軛先驗分布下，參數(shù)m,θ1,θ2的貝葉斯估計

對模擬結果進行分析，首先，通過表1、表2可以觀察到對參數(shù)選取不同的先驗分布后，分別對它們進行隨機模擬，得到的估計值與真實值的MC誤差均不超過2%，說明各參數(shù)的估計值在較高水平上是有效的；其次，可取[2.5%分位數(shù)，97.5%分位數(shù)]作為參數(shù)置信水平為0.95的置信區(qū)間，從表1、表2可以看出置信區(qū)間較窄，區(qū)間估計效果良好；最后，通過圖2、圖4可以看出變點m的兩條Markov鏈圖像趨于重合，收斂性較好. 綜上分析，利用MCMC算法得到的參數(shù)估計及變點位置估計的效果都較為理想，因此使用該方法處理Lomax分布的變點問題是可行且有效的.