因子服從指數(shù)分布的因子分析模型的參數(shù)估計(jì)研究

周?chē)?guó)瓊, 蔣文江

(1.昭通學(xué)院 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，云南 昭通 657000 ；2.云南師范大學(xué) 泛亞商學(xué)院，云南 昆明 650092)

在傳統(tǒng)的因子分析模型中，因子被假設(shè)為服從標(biāo)準(zhǔn)正態(tài)分布，可取實(shí)數(shù)域上的任何值；而在日常生活中，往往存在一些非負(fù)或非正的數(shù)據(jù)，如學(xué)生成績(jī)和壽命數(shù)據(jù)等，當(dāng)對(duì)這類數(shù)據(jù)做因子分析時(shí)，假設(shè)因子服從一個(gè)非負(fù)的概率分布(如指數(shù)分布)是較為合理的選擇.

傳統(tǒng)的因子分析模型假設(shè)模型中各個(gè)公共因子及特殊因子之間相互獨(dú)立且公共因子服從標(biāo)準(zhǔn)正態(tài)分布，因而可以直接利用極大似然估計(jì)[1-2]對(duì)模型中的參數(shù)做估計(jì)；但當(dāng)因子分析模型的因子服從指數(shù)分布時(shí)，會(huì)導(dǎo)致模型中的似然函數(shù)沒(méi)有顯式表達(dá)，不能直接使用極大似然估計(jì)法估計(jì)參數(shù)[3]；針對(duì)這一情況，本文采用蒙特卡洛方法來(lái)解決極大似然估計(jì)中似然函數(shù)沒(méi)有顯式表達(dá)的問(wèn)題，即利用基于EM算法的極大似然估計(jì)法對(duì)模型中的參數(shù)做估計(jì)[4-6]，其中EM算法中的E步采用馬爾科夫鏈蒙特卡洛方法中的M-H算法[7]從一個(gè)非常規(guī)的復(fù)雜分布中通過(guò)抽樣來(lái)完成積分計(jì)算[8].

1 模型及待估參數(shù)

設(shè)y=(y1,y2,…,yp)T是p維可觀測(cè)的隨機(jī)向量，μ=(μ1，μ2，…，μp)T是p維截距向量，Λ=(λij)p×q是p×q(p>q)維的因子載荷矩陣，x=(x1,x2,…,xq)T是q維潛在因子的隨機(jī)向量，ε=(ε1，ε2，…，εp)T是p維誤差隨機(jī)向量，本文所研究的模型為

y=μ+Λx+ε

(1)

該模型與傳統(tǒng)的因子分析模型[1]的不同之處是本文假設(shè)誤差隨機(jī)向量ε～N(0,∑)，其中∑=(σij)p×p為對(duì)角矩陣；潛在因子的隨機(jī)向量x中每一個(gè)隨機(jī)變量xk～exp(βk),k=1,2,…,q；誤差隨機(jī)向量的分量εj(j=1,2,…,p)與潛在因子隨機(jī)向量的分量xk之間相互獨(dú)立.模型中的待估參數(shù)有μ,Λ,Σ,βk,本文的核心工作就是對(duì)這些參數(shù)進(jìn)行估計(jì).

2 參數(shù)估計(jì)研究

2.1 模型參數(shù)極大似然估計(jì)的方法與原理

把模型(1)中的待估參數(shù)記為參數(shù)向量θ=(μ，Λ，∑，βk)，記Y=(y1,y2,…,yn)為觀測(cè)到的數(shù)據(jù)矩陣，X=(x1,x2,…,xn)是潛在的因子矩陣.

2.1.1 模型的似然函數(shù)

由于潛在因子隨機(jī)向量x的每一個(gè)分量之間相互獨(dú)立且xk～exp(βk),所以隨機(jī)向量x的聯(lián)合密度函數(shù)

(2)

由于誤差向量ε～N(0,∑)，所以當(dāng)潛在的因子隨機(jī)向量已知時(shí)，可觀測(cè)隨機(jī)向量y～N(μ+Λx,∑)，于是根據(jù)條件概率的定義可得y的條件概率密度函數(shù)[9]

(3)

根據(jù)式(2)和式(3)，基于觀測(cè)數(shù)據(jù)Y的似然函數(shù)

(4)

對(duì)應(yīng)的對(duì)數(shù)似然函數(shù)

(5)

(5)中的多重積分沒(méi)有顯式表達(dá)，所以很難通過(guò)極大化(5)來(lái)獲得參數(shù)θ的極大似然估計(jì)[10].

為解決上述積分問(wèn)題，考慮使用數(shù)據(jù)添加算法[11]，把潛在的因子矩陣X=(x1,x2,…,xn)看作缺失的數(shù)據(jù)添加到觀測(cè)到的數(shù)據(jù)矩陣Y=(y1,y2,…,yn)中，從而得到完全數(shù)據(jù)矩陣Z=(X,Y)，然后用EM算法來(lái)獲得參數(shù)θ的極大似然估計(jì).

2.1.2 EM算法估計(jì)參數(shù)

記Z=(X,Y)為添加數(shù)據(jù)后的完全數(shù)據(jù)矩陣，由于把潛在的不可直接觀測(cè)到的因子矩陣當(dāng)作是已知的能夠觀測(cè)到的數(shù)據(jù)矩陣添加到觀測(cè)數(shù)據(jù)Y中，式(4)和式(5)中對(duì)X的積分運(yùn)算失效，則添加數(shù)據(jù)后的對(duì)數(shù)似然函數(shù)

(6)

通過(guò)EM算法可求得式(6)中參數(shù)θ的極大似然估計(jì).

根據(jù)EM算法的原理[12]，基于θ的第t次迭代值對(duì)對(duì)數(shù)似然函數(shù)(6)中的X求期望，得到Q函數(shù)

(7)

需要計(jì)算的積分有E[xi|θ(t),Y]和E[xixiT|θ(t),Y].

利用蒙特卡洛方法，通過(guò)從密度函數(shù)p(x|y,θ)中抽取M個(gè)x的樣本x(l),l=1,2,…,M，然后根據(jù)樣本觀測(cè)值對(duì)x求期望，得

由貝葉斯公式[13]可得密度函數(shù)

(8)

該密度函數(shù)并不是某一常見(jiàn)分布[14]的密度函數(shù)，所以想要用常規(guī)的方法從(8)中抽取x的樣本非常困難[15].本文采用馬爾科夫鏈蒙特卡洛(MCMC)方法中的Metropolis-Hastings(M-H)算法來(lái)對(duì)x進(jìn)行抽樣.

根據(jù)Metropolis-Hastings算法的原理[16]，視密度函數(shù)p(x|y,θ)為產(chǎn)生樣本x的目標(biāo)密度，然后找一個(gè)容易產(chǎn)生樣本的分布作為建議分布，通常選取正態(tài)分布N(x(l-1),σ2Ω)作為建議分布[17]，其中x(l-1)是產(chǎn)生第l個(gè)x樣本的上一步第l-1步所產(chǎn)生的x樣本，Ω為目標(biāo)密度關(guān)于x的Fisher信息矩陣，σ2是任一給定的調(diào)節(jié)參數(shù)，通過(guò)調(diào)整不同的σ2值來(lái)得到樣本不同的接受概率，通常需要把樣本的接受概率控制在0.25到0.50之間[18]，而目標(biāo)分布與建議分布之間通過(guò)信息矩陣Ω進(jìn)行聯(lián)系.

根據(jù)極大似然估計(jì)原理，M步要對(duì)E步中所得的Q函數(shù)進(jìn)行極大化，從而求得參數(shù)θ的極大似然估計(jì)，即解決優(yōu)化問(wèn)題

也即求解如下方程

用矩陣微商的求導(dǎo)公式[19]分別對(duì)函數(shù)Q(θ|θ(t))中的參數(shù)μ、Λ、∑和βk求一階偏導(dǎo),得

分別令上述一階偏導(dǎo)為零，即可得每一個(gè)參數(shù)的估計(jì)式

2.2 EM算法的收斂準(zhǔn)則

對(duì)于該模型參數(shù)的極大似然估計(jì)，本文在EM算法里采用的停止準(zhǔn)則[20]是當(dāng)估計(jì)值θ前后兩次迭代值的絕對(duì)差值|θ(t)-θ(t-1)|<10-2時(shí)停止迭代，重復(fù)進(jìn)行100次后計(jì)算參數(shù)極大似然估計(jì)的平均值；在EM算法的E步中，本文用M-H方法產(chǎn)生樣本量為15 000的樣本，然后使用后5 000個(gè)樣本來(lái)進(jìn)行統(tǒng)計(jì)推斷，以保證所使用的樣本來(lái)自目標(biāo)分布.

2.3 模擬研究

基于所述的參數(shù)估計(jì)的方法及原理，下面對(duì)該模型中的參數(shù)估計(jì)問(wèn)題進(jìn)行模擬研究，從下述定義的因子服從指數(shù)分布的因子分析模型中產(chǎn)生模擬數(shù)據(jù).

首先，假設(shè)模型中潛在因子隨機(jī)向量中包含四個(gè)潛變量，即x=(x1,x2,x3,x4)T；可觀測(cè)的隨機(jī)向量中包含8個(gè)可觀測(cè)變量，即y=(y1,y2,y3,y4,y5,y6,y7,y8)T，則模型為

y=μ+Λx+ε

待估參數(shù)有

上述矩陣中的常數(shù)表示不進(jìn)行參數(shù)估計(jì).根據(jù)模型的定義，共有31個(gè)待估參數(shù)，參數(shù)取如下真值進(jìn)行模擬：

分別在樣本量n為200、300和400的情況下進(jìn)行100次重復(fù)估計(jì)并取其均值作為參數(shù)估計(jì)值(結(jié)果見(jiàn)表1).

從表1的模擬結(jié)果可知，每種樣本量下參數(shù)估計(jì)值與真實(shí)值之間的偏差都較小，這表明用MCECM算法對(duì)該模型的參數(shù)做估計(jì)具有良好的效果，即估計(jì)值與真實(shí)值非常接近.特別地，當(dāng)樣本量為200時(shí)，真實(shí)值與估計(jì)值之間偏差絕對(duì)值的最大值約為0.06，而最小值約為0.001，所以在樣本量為200的情況下，參數(shù)的估計(jì)結(jié)果較為準(zhǔn)確；另一方面，隨著樣本量的增加，參數(shù)估計(jì)值的標(biāo)準(zhǔn)誤差整體在逐漸減小，可見(jiàn)隨著樣本量的增加，參數(shù)的估計(jì)結(jié)果越來(lái)越準(zhǔn)確，雖然其標(biāo)準(zhǔn)誤減小的幅度較小，但可認(rèn)為這種參數(shù)估計(jì)方法對(duì)樣本量的要求不高，故MCECM算法對(duì)于參數(shù)估計(jì)而言具有較強(qiáng)的實(shí)用性，針對(duì)需要解決的參數(shù)估計(jì)問(wèn)題，樣本量的要求不是很高.

表1 參數(shù)模擬結(jié)果

n代表樣本量；Bias代表真實(shí)值與估計(jì)值之間的偏差；SD代表標(biāo)準(zhǔn)差.

3 結(jié)語(yǔ)

針對(duì)諸多應(yīng)用領(lǐng)域中相關(guān)數(shù)據(jù)為非負(fù)的情形，提出了一個(gè)傳統(tǒng)因子模型的替代模型；用發(fā)展成熟的MCECM算法[21-22]來(lái)對(duì)模型中的參數(shù)進(jìn)行估計(jì).研究結(jié)果表明，用MCECM算法對(duì)因子分析模型中的參數(shù)進(jìn)行估計(jì)是一種有效的方法.