Fisher得分法與EM算法在隨機(jī)效應(yīng)模型中的應(yīng)用

趙鵬輝

(大慶師范學(xué)院 數(shù)學(xué)科學(xué)學(xué)院，黑龍江 大慶 163712)

0 引言

Harville[1]提出適合分析重復(fù)數(shù)據(jù)和增長曲線的隨機(jī)效應(yīng)模型的一般形式，Laird與Ware[2]，Jennrich與Schluchter[3]，Lindstrom與Bates[4]，Diggle[5]，Chi與Reinsel[6]等都致力于把此類模型應(yīng)用于縱向數(shù)據(jù)。Davidian與Giltinan[7]對(duì)隨機(jī)效應(yīng)模型用于重復(fù)測量數(shù)據(jù)給出一般性應(yīng)用。

隨機(jī)效應(yīng)模型在生物，醫(yī)學(xué)，經(jīng)濟(jì)，金融等領(lǐng)域具有廣泛應(yīng)用，因此，進(jìn)三十年來，關(guān)于隨機(jī)效應(yīng)模型的參數(shù)估計(jì)一直是線性模型的最活躍的研究方向之一。這方面已有一些專著，其中一些將隨機(jī)效應(yīng)模型應(yīng)用于處理縱向數(shù)據(jù)分析。Ware[8]指出，除了某些最簡單的情形，迭代方法是對(duì)重復(fù)測量數(shù)據(jù)使用隨機(jī)效應(yīng)模型而找出其中參數(shù)估計(jì)的必用方法。對(duì)于方差的研究，很多學(xué)者做了不少工作，當(dāng)設(shè)計(jì)是平衡的，但觀測值是隨機(jī)丟失時(shí)，多維模型可以對(duì)缺失觀測使用多維方法。

1 固定效應(yīng)和隨機(jī)效應(yīng)的估計(jì)

1.1 固定效應(yīng)的估計(jì)

在考慮固定效應(yīng)的估計(jì)時(shí)，我們將觀測向量y表達(dá)為如下形式：

yi=Xiα+Ziβi+εi，i=1,…n

1.2 隨機(jī)效應(yīng)的估計(jì)

1)α已知，隨機(jī)效應(yīng)被經(jīng)驗(yàn)Bayes估計(jì)，即用給定數(shù)據(jù)的后驗(yàn)估計(jì)：

2 方差分量估計(jì)

2.1 極大似然估計(jì)

Laird,Lange,andStram(1987)給出邊際似然(marginalloglikelihood)：

在用ML估計(jì)時(shí)注意了以下兩個(gè)問題：

1)用MLE估計(jì)方差成分時(shí)，沒有考慮到由于模型中固定效應(yīng)的估計(jì)引起的自由度的損失；

2)用MLE是在特定參數(shù)形式的假設(shè)下推出的。一般地，數(shù)據(jù)向量的分布是正態(tài)的。

2.2 限制極大似然估計(jì)

方差成分的限制極大似然估計(jì)(REML)是不含有固定效應(yīng)的那部分似然函數(shù)最大化。REML不僅對(duì)模型的固定效應(yīng)具有不變性，而且還不含有固定效應(yīng)。當(dāng)樣本量n遠(yuǎn)大于β的維數(shù)p時(shí)，ML估計(jì)與REML估計(jì)區(qū)別將非常小。

3 Fisher得分法與EM算法在方差成分中的應(yīng)用

3.1 EM算法用在方差成分的ML估計(jì)中

令D=diag(∑,∑,…,∑)，則

于是

(1)

(2)

同理可得

(3)

(4)

(2)和(4)為EM算法的E步，(1)和(3)為M步。給定σ2與∑的初值。由(1)-(4)得到迭代方程：

3.2 EM算法用在方差成分的REML估計(jì)中

對(duì)于限制極大似然估計(jì)：

上節(jié)中的(2)與(4)進(jìn)行相應(yīng)的改動(dòng)，從而得到迭代方程：

3.3 Fisher得分法用在方差成分的ML估計(jì)中

Fisher得分法用Hessian矩陣的期望替換Hessian矩陣，由EHαθ=0可知，可以通過分開解方程獲得更新方程，則

3.4 Fisher得分法用在方差成分的REML估計(jì)中

計(jì)算REML估計(jì)，log-likelihood變?yōu)?/p>

其中Iθθ=-EHθθ，

4 模擬和結(jié)論

4.1 模擬

這部分我們主要是用Monte Carlo方法比較EM算法和Fisher得分法，模擬使用軟件為S-plus統(tǒng)計(jì)軟件，由于機(jī)器設(shè)備條件有限，每種做1000次模擬，可以基本上看出各種方法的好壞。M與ni都不應(yīng)該太大，使它們在30之內(nèi)，選取不同樣本進(jìn)行比較，取p=q=1，p=q=2即可。

模擬一：m=10,20,30；ni=5,15,25；p=q=1。

假設(shè)模型為：yij=α+βi+εiji=1,…m;j=1,…ni。表1中為在該模型下，用EM算法和Fisher Score算法對(duì)∑估計(jì)的結(jié)果比較，其中“次數(shù)”是指平均迭代次數(shù)，“誤差”是指估計(jì)值與真值的平均誤差范數(shù)。

表1 ∑估計(jì)的EM算法和Fisher Score算法的比較

模擬二：m=10,20,30；ni=5,15,25；p=q=2。

假設(shè)模型為：yij=α1+xijα2+β1i+zijβ2i+εij，i=1,…,m，j=1,…,ni。

并且

Cov(βi,εi)=0。設(shè)置ρ=0，±0.1，±0.3，±0.5，±0.7，±0.9，分別計(jì)算。

4.2 結(jié)論

從模擬情況可以知道，p=q=1時(shí)，可以看出，EM算法收斂較慢，EM算法計(jì)算方差估計(jì)的誤差隨著數(shù)據(jù)的增多逐漸減小，平均迭代次數(shù)逐漸減少；而Fisher得分法迭代收斂速度很快，誤差隨著數(shù)據(jù)的增多逐漸減小，平均迭代次數(shù)并無明顯變化。當(dāng)p=q=2時(shí)，仍然如此。模擬中發(fā)現(xiàn)使用FisherScore算法比EM算法得到的結(jié)果迭代次數(shù)要少一些，但EM算法算法每步迭代所需時(shí)間和FisherScore相比較少。畢竟實(shí)際應(yīng)用中，花費(fèi)太多時(shí)間是沒有必要的。α估計(jì)在不同的方法中，即使改變樣本量及參數(shù)的個(gè)數(shù)，其結(jié)果沒有發(fā)生改變，歸根結(jié)底其估計(jì)是廣義最小二乘；σ2的估計(jì)無論使用哪種方法其兩種估計(jì)值無明顯差別?！乒烙?jì)用REML估計(jì)比ML估計(jì)隨著數(shù)據(jù)的增多，結(jié)果相差無幾；∑估計(jì)EM算法和FisherScore算法也相差不大。差別較大的就是運(yùn)算次數(shù)，并且數(shù)據(jù)越少，相關(guān)系數(shù)絕對(duì)值越大，EM算法所需運(yùn)算次數(shù)越高，而FisherScore算法很穩(wěn)定。

[參考文獻(xiàn)]

[1] David A. Harville. Maximum Likelihood Approaches to Variance Component Estimation and to Related Problems [J]. Journal of the American Statistical Association, 1977, 72(358):320-338.

[2] Nan M. Laird and James H. Ware. Random-Effects Models for Longitudinal Data [J]. Biometrics, 1982, 38(4):963-974.

[3] Robert I. Jennrich and Mark D. Schluchter.Unbalanced Repeated-Measures Models with Structured Covariance Matrices [J]. Biometrics, 1986, 42(4):805-820.

[4] Mary J. Lindstrom and Douglas M. Bates. Newton-Raphson and EM Algorithms for Linear Mixed-Effects Models for Repeated-Measures Data [J]. Journal of the American Statistical Association, 1988, 83(4):1014-1022.

[5] Peter J. Diggle. An approach to the analysis of repeated measurements [J]. Biometrics, 1988, 44(4):959-971.

[6] Eric M. Chi and Gregory C. Reinsel. Models for longitudinal data with random ecects and AR(1) errors [J]. Journal of the American Statistical Association, 1989, 84(6):452-459.

[7] Davidian.M. and Giltinan.D.M. Nonlinear Models for Repeated Measurement Data[M]. NewYork: Chapman and Hall,1995.

[8] James H. Ware. Linear models for the analysis of longitudinal studies [J]. The American Statistician, 1985,39(2):95-101.