趙煥麗, 何幼樺

(上海大學理學院, 上海 200444)

Cox[1]提出的logit 模型是二分類問題最常用的解決方案.而對于多分類或多等級問題,McCullagh[2]擴展二分類logit 模型為比例優(yōu)勢模型.在比例優(yōu)勢模型中, 設響應變量z為具有有限個等級1<2< ··· < K的順序變量, 在自變量x ∈Rk×1處, 記pj(x) =P(z=j|x),j=1,2,··· ,K, 其累積概率

由此, 在x處,z所屬等級落入{1,2,··· ,j}或{j+1,j+2,··· ,K}的概率分別為γj(x)和1?γj(x).將{1,2,··· ,j}和{j+1,j+2,··· ,K}作為兩個分類, 由logit 模型[1]得到針對順序變量的比例優(yōu)勢模型

式中:αj關于等級j單調遞增;β ∈Rk×1為待估計參數.式(2)是一個參數回歸模型.

當響應變量與解釋變量之間的關系不能用有限個參數描述時, 一般引入非參數模型.當響應變量與部分解釋變量滿足線性關系, 而與其他解釋變量的關系不能用有限個參數描述時, 則建立半參數回歸模型是合適的.熊笛等[3]用一個連續(xù)非線性函數替換式(2)中的線性部分, 建立了一種半參數順序回歸模型

并構造了參數部分αj的最小二乘估計和非參數部分f(·)的局部線性估計量.

1986 年, Engle等[4]提出了一般形式的半參數模型

式中: (x,u)為解釋變量;β和f(·)為待估計參數和待估計函數.

考慮到半參數模型適用范圍廣泛, 本工作在模型(4)的基礎上, 建立了更一般形式的半參數順序回歸模型

式中:γj(x,u) =P(z≤j|x,u)為給定解釋變量x,u時z所屬等級不超過j的概率.模型(5)是模型(2)與模型(3)的推廣.

對于半參數模型的估計, 研究人員給出了多種經典估計方法[5-7].關于貝葉斯方法,Koop等[8]對局部線性正態(tài)回歸模型的半參數推斷發(fā)展了一種新的貝葉斯方法, 但該方法只討論了取共軛先驗下的正態(tài)線性模型, Koop等[9]使用了類似的方法.李琪琪[10]導出了半參數回歸模型中參數的貝葉斯最小風險線性無偏估計, 并討論了相對于最小二乘加權估計(least squares weighted estimation, LSWE)的優(yōu)良性.Dimitrakopoulos[11]在誤差項服從Dirichlet 過程時發(fā)展了一種估計含隨機波動的時變參數回歸模型的貝葉斯半參數方法, 并將其應用在通貨膨脹問題的分析中.Chow等[12]在非參數部分服從Dirichlet 過程時利用數值方法計算出貝葉斯估計值.Kim等[13]構造了正態(tài)先驗下參數與非參數部分的一種貝葉斯估計.

本工作針對建立的半參數順序回歸模型(5), 以非參數部分隨機過程的有限維分布作為先驗, 構造了參數部分與非參數部分的貝葉斯估計量.基于Kim等[13]提出的方法, 在正態(tài)情形下推導出估計量的解析表達式, 最后通過仿真模擬與實證分析評價貝葉斯估計量的表現.

1 貝葉斯估計

根據樣本{(xi,ui,zi),i=1,2,··· ,n}構造模型(5)中α,β,σ2,f(·)的貝葉斯估計, 從而在新狀態(tài){(xi,ui),i=n+1,n+2,··· ,n+t}下估計zi的值.下面分3 步構造α,β,σ2,f(·)的貝葉斯估計量: ①構造pj(x,u)的貝葉斯估計; ②根據式(1)計算γj(x,u); ③結合先驗分布推導α,β,σ2,f(·)的后驗分布, 從而構造貝葉斯估計量.

本工作僅就正態(tài)情形給出推導過程, 但估計方法適用于一般分布.在正態(tài)情形下, 可寫出估計量的解析表達式; 在非正態(tài)情形下, 可利用蒙特卡羅-馬爾科夫(Monte Carol-Markov chain, MCMC)等方法進行數值計算.

當(xi,ui)處的觀察值較少時,pij=pj(xi,ui)的極大似然估計值可能會取不合理的0 或1,因此首先對pij進行貝葉斯估計.取pi= (pi1,pi2,··· ,piK)的先驗為Dirichlet 分布, 響應變量服從多項分布, 即

則后驗

因此,pij的貝葉斯估計[14]為

在(xi,ui)下, 模型(5)滿足因此, 在不計一個常數的情況下,{αj,1 ≤j

當n=2,t=1,K=3,k=1 時, 模型(5)的樣本模型為

T=為n階單位矩陣的每一行重復K?1 行,f(u)=(f(u1),f(u2),··· ,f(un+t))′,ε=(ε11,··· ,ε1,K?1,ε21,··· ,εn,K?1)′, 則

為對t個新樣本預測其分類, 需估計參數β?和函數f(·).但f(·)是一個隨機過程, 有無窮多個參數, 因此無法直接估計函數f(·), 但可估計f(·)在un+1,un+2,··· ,un+t處的函數值下面對正態(tài)情形構造參數部分與非參數部分的貝葉斯估計量.

在擾動項ε服從正態(tài)分布的情況下, 有

為了構造參數和非參數部分的貝葉斯估計量, 設參數部分β?的先驗服從多維正態(tài)分布,函數f(·)服從高斯過程(Gaussian process, GP), 則f(·)在u處的函數值服從高斯過程的有限維分布N(·,·).σ的先驗分布取Jeffreys 先驗, 即

以下推導過程中將f(u),f0(u)分別簡記為f,f0.

定理1記的后驗分布為正態(tài)-逆伽馬分布.

證明 根據貝葉斯理論, 在先驗(15), (16)和(17)滿足的情況下, (β?,f,σ)的后驗分布為

式中:

式(18)的推導主要用到矩陣二次型的化簡.

注意到, 對任意v ?=0∈RK?1+k+n+t,

式(23)中的“>”號由v1,v2不同時為0 保證.因此, 矩陣A正定, 式(19)中的A?1存在.

由式(18)可知,θ,σ|Y,x,u服從正態(tài)-逆伽馬分布(normal-inverse-gamma distribution,NIG), 即分別對σ2,θ積分, 可求得θ,σ2的邊際分布.

推論1后驗分布

在二次損失下,θ的貝葉斯估計為后驗期望估計即可得^θ?,f.

推論2后驗分布σ2|Y,x,u服從逆伽馬分布(inverse-gamma distribution, IGa), 即

在二次損失下,σ2的貝葉斯估計為后驗期望估計

對新的樣本(x,u), 根據

預測

2 數值模擬

觀察在不同先驗分布、樣本點處的觀察次數下^β?,f,^σ2的表現.設回歸模型樣本量n= 70, 預測t= 30 步, 研究不同先驗分布, 不同觀察次數下的表現.先驗設置中, 有

表1 N =1 000 時的參數估計Table 1 Estimation of parameters when N =1 000

(2) 為探究先驗分布對估計量效果的影響, 以f(u)為例, 取先驗均值f0(u)=16u(0.5?u),重復觀察次數m=3(隨機種子為123), 由推論1 估計得到的與真實的f(u) = sin(2πu)的比較如圖1(b)所示.

圖1 非參數部分f 與估計的比較Fig.1 Comparing of nonparametric f and estimation

圖1(a)表明, 在f(u)的先驗分布的均值f0(u) = 0 時,會對先驗分布有一定修正, 估計值趨向于真實的f(u).每個樣本點處的重復觀察次數m越多,的表現越好.表1 中與真實的α2,α3差異不大,與真實的β,σ2有一定差異, 但會隨著重復觀察次數m的增大越來越接近真實值.從,的1 000 次模擬的方差來看, 每個樣本點處的重復觀察次數m越多,,,,的估計值越穩(wěn)定.

對比圖1(b)所示的f0(u)=16u(0.5?u)與f(u)=sin(2πu)的曲線,在區(qū)間[0,0.7]內f0(u)接近f(u),估計得到的接近真實值;在區(qū)間[0.7,1]內f0(u)偏離f(u)很大,端點處f0(1)=?8,估計得到的很大程度上修正了先驗分布與真實值之間的偏離,端點處接近?2.3.因此在先驗分布接近真實值時,的估計效果較理想, 當先驗分布偏離真實值時,也有不錯的表現.同理, 可探究先驗分布的選取對的影響.

3 實證分析

對收入等級預測問題建立一個半參數順序回歸模型, 以食品、衣著、居住、家庭設備及用品、交通通信、文教娛樂、醫(yī)療保健和其他共8 項主要生活性消費支出[15]的占比作為解釋變量, 根據消費結構預測收入等級.在很多經濟調查中, 真實收入的收集非常困難, 本工作根據消費習慣來預測收入等級, 有助于基于收入進行數據分析.

根據國家統(tǒng)計局2002—2012 年人均消費支出占比的統(tǒng)計數據(http://data.stats.gov.cn/easyquery.htm?cn=C01), 對不同收入等級居民的消費結構進行實證分析.8 項消費支出占比用每項的消費支出/總消費支出計算得到.農村居民家庭的收入等級采取國家統(tǒng)計局的五等份劃分法, 分為低收入、中等偏下收入、中等收入、中等偏上收入、高收入家庭(分別用1, 2, 3, 4,5 表示); 城鎮(zhèn)居民家庭的收入等級在國家統(tǒng)計局的劃分法基礎上稍作改動, 將低收入與較低收入家庭劃分為低收入家庭, 較高收入與最高收入家庭劃分為高收入家庭, 這樣城鎮(zhèn)居民家庭由原來的7 個收入等級改為5 個收入等級.

選取2002—2011 年的數據(共120 組)作為訓練樣本, 用半參數順序回歸模型(5)擬合8 項人均消費支出占比與家庭收入等級之間的關系.令u為8 項消費支出占比,

則由模型

從表2 的實驗結果可以看到, 對收入五等級問題的12 組樣本預測中, 預測準確率為58.33%.預測錯誤的5 組樣本中, 有4 組樣本的預測等級與實際等級只相差1 個等級.

表2 半參數順序回歸模型的外推收入等級Table 2 Extrapolation income level of semiparanetric ordinal regression mode

已有研究表明, 在諸多影響消費需求的因素中, 收入水平始終是影響消費需求的最重要因素.而本工作實證分析的結果表明, 以居民家庭的消費結構為解釋變量可以較準確地預測收入等級, 因此消費結構反過來也反映了家庭的收入情況.以8 項生活性消費支出考慮家庭的收入等級的方式, 比采用恩格爾系數更為全面.綜上可知, 在小樣本情況下, 貝葉斯估計由于利用了先驗信息, 往往有更好的估計效果.

4 結束語

本工作在正態(tài)情形下構造了半參數回歸模型中參數與非參數部分的貝葉斯估計, 多次模擬結果表明, 在先驗均值均取0 時仍然有不錯的估計效果, 在先驗分布接近真實值時, 估計效果會更理想.

相比模型(3), 本工作建立的模型是基于比例優(yōu)勢模型和半參數模型(4)所形成的更一般的半參數順序回歸模型.模型(3)只考慮了被解釋變量與非參數部分變量u的關系, 無法反映參數部分變量x的影響, 而本模型建立了解釋變量與被解釋變量的半參數關系, 適用范圍更廣泛.例如, 當x為外生變量時, 本模型可反映不同情況下u對被解釋變量的影響.在根據消費結構預測家庭收入等級的實例中, 本工作將線性部分取為啞變量, 使模型對城鎮(zhèn)居民和農村居民的收入等級預測問題都適用, 同時線性部分的參數反映了相同收入等級的城鎮(zhèn)居民家庭與農村居民家庭的消費結構差異.以2002—2011 年間的數據作為樣本, 對2012 年的數據作預測,即在外推情形下, 對五等級預測問題的預測準確率達到58.33%, 預測錯誤的樣本中, 預測等級與實際等級大多只相差1 個等級.