基于Gibbs抽樣算法的貝葉斯動態面板數據模型分析

摘 要 針對現有動態面板數據分析中存在偶發參數和沒有考慮模型參數的不確定性風險問題，提出了基于Gibbs抽樣算法的貝葉斯隨機系數動態面板數據模型.假設初始值服從平穩分布，自回歸系數服從Logit正態分布的條件下，設計了Markov鏈Monte Carlo數值計算程序，得到了模型參數的貝葉斯估計值.實證研究結果表明：基于Gibbs抽樣方法的貝葉斯動態面板回歸模型能有效地揭示跨截面滯后變量對響應變量的位置、尺度和形狀的影響.

關鍵詞 動態面板數據;MCMC;Gibbs抽樣算法;貝葉斯推斷;后驗分布

中圖分類號 F064.1， O212.8 文獻標識碼 A

Bayesian Analysis for Dynamic Panel Data Models Using Gibbs Sampling Algorithm

ZHU Hui-ming1， ZHOU Shuai-wei1， LI Su-fang1， ZENG Zhao-fa2

（1. School of Business Administration， Hunan University， Changsha， Hunan 410082， China;

2. School of Finance and Statistics， Hunan University， Changsha， Hunan 410079， China）

Abstract We developed thehierarchical Bayesian approach for inference in random coefficient dynamic panel data models. Our approach allows for the initial values of each unit's process to be correlated with the unit-specific coefficients.A stationarity assumption was imposed on each unit's process by assuming that the unit-specific autoregressive coefficient is drawn from a logitnormal distribution. The research shows that our approach can efficiently reveal how the lagged variables affect the location， scale and shape of the response variable across section.

Key words dynamic panel data; MCMC;Gibbs sampling algorithm; Bayesian inference; Posterior distribution

1 引 言

自Mundlak（1961），Balwstra和Nerlove（1966）把面板數據引入到計量經濟學以來，面板數據的理論與應用研究得到了眾多學者的關注.經濟分析中經常會遇到具有三維（個體、時間、指標）信息的數據結構，即面板數據（panel data），是指在時間序列上取多個截面，在這些截面上同時選取樣本觀測值所構成的樣本數據，也就是把截面數據和時間序列數據融合在一起的數據.由于面板數據利用了更多數據的信息，提高了自由度和有效性，能得到更有效和更可靠的參數估計量，從而能夠更加精確地估計復雜的行為方程，檢測和度量單純使用橫截面數據或時間序列數據無法觀測到的影響.目前，面板數據在經濟管理學、社會學、心理學等領域都得到了廣泛應用.

動態面板數據模型常用于描述很多經濟變量的動態關系，因為經濟變量不僅受到本期因素的影響，而且受到非本期因素的影響.比如，物價往往受到上期物價的影響；各部門的投資額不僅影響本期生產量而且影響今后各期生產量.線性或非線性參數動態面板數據模型廣泛應用于當今經濟研究［1］，這些模型都是建立在各截面自回歸系數相同且為常數的假設之上.然而，在很多研究中，將各截面自回歸系數假設為跨截面變化更為符合實際.如購買力平價（Purchasing Power Parity， PPP）的回復速度（Speed of Reversion）［2］，個體對收入沖擊的調整速度［3］，以及各國儲蓄行為的差異等［4］.Yu，Jong和Lee使用偽極大似然估計方法研究了具有大N和大T的固定效應空間動態面板數據模型［5］，但沒有考慮外生變量對模型參數估計的影響.Pesaran和Smith研究發現，即使在具有大N和大T的面板數據中，忽略各截面的異質性會導致自回歸系數估計非一致［6］.因此，Pesaran和Smith提出了使用組平均估計，即對各截面分別回歸估計的系數進行平均.Hsiao，Pesaran和Tahmiscioglu指出，只要N/T→0，N→

，T→

，組平均估計是模型系數的一致漸進正態估計［7］.然而對于大N和小T的傳統面板數據來說，組平均估計為非一致估計.Hsiao，Pesaran和Tahmiscioglu提出了使用層次貝葉斯方法對隨機效應自回歸面板數據模型進行參數估計，指出對于大N和大T的面板數據，他們提出貝葉斯估計量漸進等價于組平均估計量，并且通過Monte Carlo實驗，說明其方法在小樣本下要優于組平均估計.Park，Sickles和Simar研究了僅含一個滯后因變量的動態面板數據模型的半參數有效估計量［8］.朱慧明基于VAR模型對中國居民消費水平進行貝葉斯單位根檢驗，解決了向量自回歸模型超參數估計的難題，克服了經典單位根檢驗在經濟時序小樣本下功效偏低的缺陷［9］.汪朝暉利用靜態面板數據模型研究了湖南省城鎮居民收入與消費的關系［10］；王立平基于空間動態面板數據模型研究了中國產業結構變遷對區域經濟增長的影響，發現產業結構變動對地區經濟增長存在顯著的促進作用，且外溢性顯著［11］.Ciarreta和Zarraga則利用1970~2007年的年度數據結合動態面板數據方法研究了歐洲12國電力消費和實際GDP的長期因果關系［12］，結果顯示電價與GDP之間存在雙向因果關系；Sarafidis，Yamagata和Robertson使用廣義矩估計方法（GMM）研究了英國企業勞動力雇傭情況［13］.

假定變量初始值服從平穩分布，相比將初始變量假定為固定值的傳統假設更為合理；同時，假定待估計隨機系數服從Logit正態分布，保證了自回歸過程的平穩性.在貝葉斯理論框架下，對含外生變量的動態隨機系數面板數據模型進行了貝葉斯推斷，并且使用層次貝葉斯超參數先驗設置，利用貝葉斯方法對模型進行參數估計，以我國東中西部六省的城鎮居民人均消費和可支配收入構成的面板數據進行實證研究，發現貝葉斯方法可以有效地對動態面板數據模型進行參數估計.

2 模型結構分析

考慮動態面板數據模型：

yi=δiyi，－1+Xiβi+ui，－1＜δi＜1，i=1，2，…，N. （1）

此處yi=(yi1，yi2，…，yiT)'為截面i各個時刻的觀察值序列，yi，－1=(yi0，yi1，…，yi，T－1)′為觀察值的一階滯后序列，Xi=(xi1，xi2，…，xiT)′為截面i的外生變量觀察值序列，隨機誤差項ui=(ui1，ui2，…，uiT)′，且ui相互獨立同分布.另外，式(1）可寫為：

yi=Ziωi+ui，Zi=(yi，－1，Xi)，ωi=(δi，βi)′.（2）

關于此模型，做出假設：

（ⅰ）假設各截面為平穩過程，即|δi|＜1，故令δi服從Logit正態分布，即

δi=2{exp (ωi)/［1+exp (ωi)］－0.5}.

（ⅱ）ωi=+εi，=(，′)′，εi=(ε1i，ε′2i)′；εi相互獨立同分布，均值為0，協方差為Δ，即

ωi~N(，Δ)，Cov(ωi，ωj)=0(i≠j).

（ⅲ）Xi為嚴格外生變量，即E(ui|Xi)=0；

（ⅳ）矩陣T－1(X'iXi)滿秩且隨著T→

，趨向于有限非奇異矩陣；

（ⅴ）隨機誤差項異方差且截面不相關，即

uit~i.i.d.N(0，σ2i)，E(uiu'j)=0(i≠j).

對于有限T，需要對被解釋變量的初始值進行設定.可以假設yi0固定且已知，即與截面系數不相關；或者假設其服從平穩分布.前者顯然與實際不符，若過程持續一段時間，沒有理由認為yi0不同于yit，因此假設初始值yi0產生于高斯分布N(βi/(1－δi)，σ2u/(1－δ2i)|βi，δ).主要考慮對參數和Δ的估計，并對所得估計的效果進行評價.

3 貝葉斯分析

3.1 參數的貝葉斯推斷

參數的先驗分布設置是貝葉斯統計分析的前提條件，在此動態面板數據模型中，參數的先驗分布可以層次先驗方法推導：y=(y′1，y′2，…，y′N)′的聯合密度函數為y(NT×1)~N((NT×Nm)ω(Nm×1)，Ω)，其中=diag(Z1，Z2，…，ZN)，ω=(ω′1，ω′2，…，ω′N)′，Ω為分塊對角矩陣，第i塊為Ωi=σ2iIT；ω~N(A，Γ)，A=eNIm，eN為元素全為1的列向量，Γ=INΔ；同樣可假設~N(μ，Ψ).

y的邊緣分布為正態分布N(Z，Φ)， Z=(Z1'，Z2'，…，ZN')'，Φ=Ω+Γ'，從而參數的后驗分布為：

p(|y，y0)∝p(y|)p()

∝exp {－(1/2)［(y－Z)′Φ－1(y－Z)+

(－μ)′Ψ－1(－μ)］}∝exp ［′(Z′Φ－1Z+

Ψ－1)－2(Z′Φ－1y+Ψ－1μ)′+C］∝exp {［－

(Z'Φ－1Z+Ψ－1)－1(Z'Φ－1y+Ψ－1μ)］′

(Z'Φ－1Z+Ψ－1)［－(Z′Φ－1Z+Ψ－1)－1

(Z'Φ－1y+Ψ－1μ)］+C}，

其中，C為常數；顯然

~N((Z′Φ－1Z+Ψ－1)－1(Z'Φ－1y+Ψ－1μ)，

(Z′Φ－1Z+Ψ－1)－1).

特別地，若Ψ－1→0，則

~N((Z'Φ－1Z)－1Z′Φ－1y，(Z′Φ－1Z)－1).

同時，經過代數變換，

ω*=(Z′Φ－1Z+Ψ－1)－1(Z′Φ－1y+Ψ－1μ)

=∑Ni=1Wii={［∑Nj=1(σ2j(Z′jZj)－1+

Δ)－1］－1［σ2i(Zi'Zi)－1+Δ］－1}i.

在實際應用中，由于上式的方差部分無法獲知，導致它一般只能用于經濟分析中的模型比較，從而具有一定的局限性.因此，考慮將超參數的先驗分布引入模型，從參數的聯合后驗密度中得到超參數的邊緣后驗密度.在上述動態面板數據模型中，參數ωi，，Δ，σ2i的聯合先驗分布為：

π(ωi，，Δ，σ2i)=π(ωi|，Δ)π(|Δ)

#8226;π(Δ)∏Ni=1π(σ2i) （3）

顯然，由假設2易見模型參數ωi關于(，Δ)的條件分布π(ωi|，Δ)為正態分布ωi~N(，Δ)；同時根據Hsiao［14］和Nandram［15］的觀點，各參數的先驗分布為：

σ2i~IG(σ2i)，Δ~IWυ0(Λ－10)，|Δ

~N(μ0，Δ/κ0).

此處IG(σ2i)表示逆伽瑪分布（Inverse-Gamma Distribution）；IWυ0(Λ－10)表示自由度為υ0，尺度矩陣為Λ0的逆Wishart分布（Inverse-Wishart Distribution）；μ0為先驗均值，κ0為Δ的先驗測量個數.

根據貝葉斯定理，參數的聯合后驗分布密度函數與模型似然函數的乘積成正比，即

π(ωi，，Δ，σ2iy，yi0)∝∏Ni=1π(y，yi0|ωi)

π(ωi|，Δ)π(|Δ)π(Δ)#8226;∏Ni=1π(σ2i)

∝∏Ni=1σ－Tiexp ［(－1/2)∑Ni=1σ－2i(yi－

ZiωI)'(yi－Ziωi)］×Δ－N2exp ［－

12∑Ni=1(ωi－)'Δ－1(ωi－)］×

Ψ－12exp ［(－12(－μ)'Ψ－1(－μ)］

×Δ－(υ0+N2+1)exp ［－12tr(Λ0Δ－1)

－κ02(－μ0)'Δ－1(－μ0)］×∏Ni=1σ－2i. （4）

顯然，參數的聯合后驗分布的形式比較復雜，下面分別討論各參數的后驗條件分布.

（ⅰ）參數ωi的后驗條件分布.根據條件概率定義，參數ωi關于(，Δ，σ21:N)的后驗條件分布密度函數為

π(ω|y1:N，，Δ－1，σ21:N)∝N［Ai(σ－2iZi'yi+

Δ－1)，Ai］，Ai=(σ－2iZ′iZi+Δ－1)－1.

（ⅱ）參數的后驗條件分布.類似地，參數關于(Δ，σ21:N)的后驗條件分布密度函數為

π(ω|y1:N，ω1:N，Δ－1，σ21:N)∝N［B(NΔ－1+

ψ－1μ)，B］，B=(NΔ－1+Ψ－1)－1，

=1N∑Niωi.

（ⅲ）參數Δ的后驗條件分布.類似地，參數Δ關于(ω1:N，，σi1:N)的后驗條件分布為Wishart分布

π(Δ－1|y1:N，ω1:N，，σ21:N)

∝W{［∑Ni=1(ωi－)(ωi－)'+Λ0］－1，υ0+N}.

（ⅳ）參數σ2i關于(ω1:N，，Δ－1)的后驗條件分布為逆伽瑪分布

π(σ2i|y1:N，ω1:N，，Δ－1)

∝IG{T/2，［(yi－Ziωi)'(yi－Ziωi)］/2}

將上述方法獲得的關于參數和Δ的估計稱為層次貝葉斯估計.

3.2 Gibbs抽樣算法

根據上述參數的后驗條件分布以及遍歷性定理，可以利用基于Gibbs抽樣的MCMC數值算法進行模擬仿真，以獲得模型參數的貝葉斯估計及其置信區間）.Gibbs抽樣算法是一種廣泛應用的Metropolis-Hastings抽樣方法，通過建立Markov鏈，對未知變量Uk進行模擬，當鏈達到平穩，收斂到目標分布，即得所求的后驗分布.Gibbs抽樣過程如下：首先，給定參數的初始值：(θ(0)1，θ(0)2，…，θ(0)k)，然后根據它們的條件后驗分布循環抽樣，抽樣的迭代步驟為：

1）從后驗條件分布f(θ1|θ(i)2，…，θ(i)k，y)抽取θ(i+1)1；

2）從后驗條件分布f(θ2|θ(i+1)1，…，θ(i)k，y)抽取θ(i+1)2；

…………

k）從后驗條件分布f(θk|θ(i+1)1，…，θ(i+1)k－1，y)抽取θ(i+1)k.θ(0)，θ(1)，…，θ(m)，…形成一條Markov鏈，θ′到θ的轉移概率為：

f(θ′，θ)=f(θ1|θ′2，…，θ′k，y)f(θ2|θ1，

θ′3，…，θ′k，y)…f(θk|θ1，…，θk－1，y)

當Markov鏈循環迭代m(m＜n)次收斂后，由Monte Carlo積分公式可以得到各個參數的后驗均值和方差分別為：

E(θk)≈1n－m∑nt=m+1θ(t)k，

Var(θk)≈1n－m－1∑nt=m+1［(θ(t)k)2－

1n－m∑nt=m+1θ(t)k)2］;

其中，θk表示表示(ωi，，Δ，σ2i)中的任一參數.

4 實證分析

選取我國東中西部六省山東、廣東、河南、湖南、甘肅、貴州的1997~2008年的居民收入、消費年度數據建立城鎮居民的消費模型，對各省的居民消費結構進行分析.模型中的被解釋變量CS為城鎮居民人均全年消費；同時，根據美國經濟學家杜森貝利的觀點：人們的消費具有慣性，前期消費水平高，會影響下一期的消費水平，故將前期城鎮居民人均全年消費和當期城鎮居民人均全年可支配收入YD作為模型的解釋變量.首先對各截面時間序列數據的統計特征進行分析，然后建立隨機系數動態面板數據消費模型，之后對模型進行貝葉斯參數估計.

1）數據特征分析

由表1可見，六省城鎮居民家庭平均消費序列中，均值廣東最高，河南最低；方差廣東最大；甘肅最低，貴州方差最小.Jarque-Bera檢驗統計量在1%的水平上不能拒絕正態分布的原假設，可以認為各省人均消費時間序列服從正態分布.從表2可以發現，六省城鎮居民可支配收入序列中，廣東的均值和方差最大，甘肅的均值和方差最小.Jarque-Bera檢驗統計量在1%的水平上也不能拒絕正態分布的原假設，說明各省人均可支配收入時間序列也服從正態分布.

2）模型形式設定

設CSit為第t年第i省城鎮居民家庭平均消費支出，CSit－1表示前期第i省城鎮居民家庭平均消費支出，YDit為第t年第i省城鎮居民人均可支配收入.針對人均消費和人均可支配收入的時間序列進行數據統計特征分析之后，建立隨機系數動態面板數據模型：

CSit=αi+γiCSit－1+βiYDit+uit，

i=1，2，…，6，t=1，2，…，T，

其中，CSit－1為消費的前一期值，綜合反映了前一期收入和消費對當前消費變量的影響，αi為各省的自發消費水平.

圖1 六省城鎮居民人均消費、前期人均消費和人均可支配收入構成的散布圖

3）參數估計

根據所構建的動態面板數據消費模型，基于層次貝葉斯方法，利用MCMC方法估計模型的參數.在模型運行中共進行了200 000次抽樣，為消除初始值的影響，將迭代的前100 000次舍棄，然后用100 001次到200 000次迭代得到的樣本進行參數估計.圖2~4給出了各截面參數α，β，γ的動態迭代軌跡.從迭代軌跡圖可以發現各參數的兩條馬爾可夫鏈是收斂的，說明MCMC仿真過程是平穩的.圖5~7給出了各截面參數α，β，γ的Gelman-Rubin統計量收斂性診斷，從圖中可以看出各參數的GR統計量隨著迭代次數增加趨近于1，表明抽樣方法的收斂性很好.

圖2 參數α的動態迭代軌跡

圖3 參數β的動態迭代軌跡

圖4參數γ的動態迭代軌跡

圖5 參數α的GR統計量收斂性診斷圖

圖6 參數β的GR統計量收斂性診斷圖

圖8~10給出了各截面參數α，β，γ的后驗分布的核密度估計曲線，其形狀比較平滑隨著截面個體的不同而改變，說明Gibbs迭代過程有效地模擬了模型中各參數的邊緣后驗分布.根據抽樣結果，利用MCMC算法可以進行動態面板數據模型的貝葉斯估計，表3給出了各截面參數α，β，γ的后驗均值、標準差、MC誤差、2.5％分位數和97.5％分位數的貝葉斯估計值.

圖7 參數γ的GR統計量收斂性診斷圖

圖8 參數α的后驗密度曲線

圖9 參數β的后驗密度曲線

圖10 參數γ的后驗密度圖

通過對表3的分析，可以得出如下結論：

1）參數α表示自發消費水平，α4的后驗估計值為161.5，可以發現湖南的自發消費水平在六省中最高；

2）參數γ表示前期消費對當期消費的影響程度，即消費慣性，從參數估計值上看，γ4的后驗估計值為0.5342，六省中消費慣性最高的為湖南，說明對于湖南來說，前期消費水平對當期消費水平有較大影響，應當重視消費慣性對居民消費水平的影響；

3）參數β表示可支配收入對當期消費的影響程度，從參數估計值上看，β2的后驗均值為0.398 1，說明可支配收入對消費影響最高者也是廣東，即由于廣東有較高的人均可支配收入，對當期消費可形成較大解釋力度.

另外，從參數估計的其他性質來看，所構建的隨機系數動態面板數據模型較好地擬合了不同截面個體的城鎮居民消費模式，其中參數估計的標準差和MC誤差都比較小，說明貝葉斯隨機系數動態面板數據模型是有效的.

5 結束語

基于貝葉斯理論構建了含外生變量的隨機系數動態面板數據模型，考慮經濟運行實際將經濟變量初始值假定為服從平穩分布；為保證待估計參數的平穩性，假定自回歸模型參數服從Logit正態分布.利用層次先驗方法基于Gibbs抽樣算法實現了對隨機系數動態面板數據模型的參數估計，解決了貝葉斯方法在應用中遇到的數值計算問題，模擬了各參數的MCMC迭代軌跡，參數的貝葉斯估計和參數的邊緣后驗分布.各參數的MCMC迭代軌跡是收斂的，說明Gibbs抽樣方法很好地模擬了參數的邊緣后驗分布，且各參數模擬的標準差和MC誤差都比較小，且參數的邊緣后驗密度呈鐘形，由此可見貝葉斯隨機系數動態面板數據模型的有效性.如何利用貝葉斯方法分析高階動態面板數據模型以及在經濟金融中的應用有待進一步研究.

參考文獻

［1］ J Heckman，ELeamer. Handbook of Econometrics ［M］. Amsterdam: Elsevier， 2001， 3229-3296.

［2］ J Imbs，H Mumtaz，MRavn， et al. PPP strikes back:aggregation and the Real Exchange Rate ［J］. The Quarterly Journal of Economics， 2005， 120(1): 1-43.

［3］ X Hu，S Ng. Estimating Covariance Structures of Dynamic Heterogeneous Panels［R］.Colambia:Columbia University Working Paper.

［4］ N Haque，MPesaran，S Sharma. Neglected Heterogeneity and Dynamics in Cross-country Saving Regressions［R］. IMF Working Paper， 1999.

［5］ J Yu，R Jong，L Lee.Quasi-maximum likelihood estimators for spatial dynamic panel data with fixed effects when both n and T are large ［J］. Journal of Econometrics， 2008， 146(1): 118-134.

［6］ M Pesaran，R Smith. Estimating long-run relationship from dynamic heterogeneous panels ［J］. Journal of Econometrics， 1995， 68(1): 79-113.

［7］ C Hsiao，M Pesaran，Atahmiscioglu. bayes estimation of short-run coefficients in dynamic panel data models［C］//C Hsiao，KLahiri，LLee， et al. analysis of panels and limited dependent variables:a volume in honour of G S Maddala. Cambridge: Cambridge University Press， 1999: 268-296.

［8］ B Park，R Sickles，LSimar. Semiparametric efficient estimation of dynamic panel data models ［J］. Journal of Econometrics， 2007， 136(1): 281-301.

［9］ 朱慧明，李素芳. 基于VAR模型的中國居民消費水平貝葉斯單位跟檢驗 ［J］. 財經理論與實踐， 2008， 29(154): 97-101.

［10］汪朝暉， 劉萬榮. 利用面板數據模型分析湖南城鎮居民的收入與消費之間的關系 ［J］. 經濟數學， 2007， 24(1): 54-58.

［11］王立平， 王建. 中國產業結構變遷對區域經濟增長分析——基于空間動態面板數據模型 ［J］. 統計與信息論壇， 2010， 25(7): 92-97.

［12］A Ciarreta，A Zarraga. Economic growth-electricity consumption causality in 12 European countries: A dynamic panel data approach ［J］.Energy Policy， 2010， 38(7): 3790-3796.

［13］V Sarafidis，TYamagata，D Robertson. A test of cross section dependence for a linear dynamic panel model with regressors ［J］.Journal of Econometrics， 2009， 148(2): 149-161.

［14］C Hsiao.Analysis of panel data ［M］. Cambridge: Cambridge University Press. 2003:69-109.

［15］B Nandram，JPetruccelli. A bayesian analysis of autoregressive time series panel data ［J］. Journal of Business Economic Statistics， 1997， 15(3): 328-334.

注：本文中所涉及到的圖表、注解、公式等內容請以PDF格式閱讀原文