基于模態(tài)回歸的半?yún)?shù)部分線性模型的穩(wěn)健估計(jì)

高佳佳，何曉霞

(武漢科技大學(xué)理學(xué)院，湖北 武漢，430065)

半?yún)?shù)回歸模型結(jié)合了線性回歸和非線性回歸模型，既包含參數(shù)分量又包含非參數(shù)分量。參數(shù)分量用于對(duì)確定性因素進(jìn)行分析，而非參數(shù)部分能夠?qū)﹄S機(jī)干擾因素進(jìn)行刻畫。與線性模型相比，半?yún)?shù)回歸模型更具靈活性，能更好地解釋每一個(gè)變量的效應(yīng)，因此，其在理論研究和實(shí)際應(yīng)用中都有重要意義。

針對(duì)具體問題，通常假設(shè)數(shù)據(jù)來自于某一參數(shù)模型或非參數(shù)模型，然而許多實(shí)際問題并沒有那么簡(jiǎn)單。例如，影響考察對(duì)象(指標(biāo)Y)的因素(解釋變量)可分為兩部分，即(X1,X2,…,Xp)及T，根據(jù)經(jīng)驗(yàn)或歷史資料可以認(rèn)為因素 (X1,X2,…,Xp) 是主要的，Y與(X1,X2,…,Xp)線性相關(guān)，而T則是某種干擾因素(或看作為協(xié)變量)，它同Y的關(guān)系是完全未知的，而且沒有理由將其納入誤差項(xiàng)，如果用非參數(shù)回歸加以處理，則會(huì)丟失太多的信息，若采用線性回歸方法，一般擬合情況很差，這種情況下就可采用半?yún)?shù)回歸模型。

為了解決異常值的問題，研究人員開始考慮穩(wěn)健估計(jì)方法。Zhu等[7]提出針對(duì)大維度協(xié)變量的穩(wěn)健估計(jì)方法。Yao等[8]基于局部模態(tài)回歸建立了一種用于非參數(shù)模型的估計(jì)方法，能根據(jù)觀測(cè)數(shù)據(jù)自動(dòng)調(diào)整參數(shù)。該估計(jì)方法不僅在數(shù)據(jù)集包含異常值或者誤差分布重尾的時(shí)候具有穩(wěn)健性，還能滿足數(shù)據(jù)集沒有異常值或者誤差分布為正態(tài)分布時(shí)的漸進(jìn)最小方差性。Zhao等[9]基于模態(tài)回歸研究了半?yún)?shù)變系數(shù)部分線性模型。本文擬運(yùn)用模態(tài)回歸來研究半?yún)?shù)部分線性模型中參數(shù)和非參數(shù)部分的估計(jì)問題，探討估計(jì)量的大樣本性質(zhì)。

1 模型描述

半?yún)?shù)部分線性模型的一般形式可以表示為

(1)

式中：Yi為響應(yīng)變量；Xi=(xi1,xi2,…,xip)T；β=(β1,β2,…,βp)T為待估未知參數(shù)；(Xi,Ti)是獨(dú)立同分布的隨機(jī)設(shè)計(jì)或固定非隨機(jī)設(shè)計(jì)點(diǎn)列；εi是獨(dú)立同分布的隨機(jī)誤差項(xiàng)；g(·)是定義在R上的未知光滑函數(shù)。

2 穩(wěn)健估計(jì)方法

2.1 模態(tài)回歸

作為衡量指標(biāo)，均值、中值和模是誤差分布的3個(gè)重要數(shù)值特征，其中，中值和模在處理異常值上具有同等的穩(wěn)健性。另外，模態(tài)回歸對(duì)于大多數(shù)的數(shù)據(jù)可以提供有意義的點(diǎn)預(yù)測(cè)，并且對(duì)于相同長(zhǎng)度的區(qū)間估計(jì)，當(dāng)誤差分布不規(guī)則時(shí)，模態(tài)回歸相比于其他方法能預(yù)測(cè)更大的范圍，同時(shí)預(yù)測(cè)結(jié)果更有意義。

(2)

2.2 估計(jì)方法

假設(shè)(Xi,Ti)i=1,…,n是模型(1)的獨(dú)立同分布樣本。由于g(·)是未知的非參數(shù)函數(shù)，在文獻(xiàn)[8]中采用局部多項(xiàng)式來近似g(·)。對(duì)于半?yún)?shù)模型，局部多項(xiàng)式估計(jì)有兩個(gè)缺點(diǎn)：①β是一個(gè)全局參數(shù)，為了得到它的最優(yōu)相合估計(jì)，需要采用兩步估計(jì)法；②局部多項(xiàng)式估計(jì)的計(jì)算量非常大，尤其是在高維模型中。

(3)

3 估計(jì)量的性質(zhì)

(C1) 指標(biāo)變量T具有有界支撐Ω，其密度函數(shù)fT(·)為正，并且有連續(xù)的二階導(dǎo)數(shù)。不失一般性，這里假設(shè)Ω=[0,1]。

(C2) 函數(shù)g(·)是區(qū)間[0,1]上r階連續(xù)可微的函數(shù)，其中r>2。

(C4) 令t1,…,tK為[0,1]區(qū)間的內(nèi)部節(jié)點(diǎn)，t0=0,tK+1=1，ξi=ti-ti-1，ξ=max{ξi}，存在常數(shù)C0使得

(C5)F(x,t,h)關(guān)于(x,t)連續(xù)。

(C6) 對(duì)于任意的h>0，有F(x,t,h)<0。

(4)

Ξ(β,θ)=

針對(duì)I1進(jìn)行泰勒展開，有

結(jié)合條件(C4)、(C7)和‖U(ti)‖=O(K-r)，可以得到

=Op(nK-r‖v‖)

(5)

因此I1=Op(nδK-(r+1)‖v‖)=Op(nδ2K-1‖v‖)。

對(duì)于I2，可以證明

I2=E[F(X,T,h)]Op(nK-1δ2‖v‖2)

(6)

因此，若選擇足夠大的常數(shù)C，則I2可通過‖v‖=C控制I1。

同樣可以證明

I3=Op(nK-1δ3‖v‖3)

(7)

當(dāng)n→的時(shí)候，有δ→0，因此δ‖v‖→0，從而有I3=op(I2)。故I2通過‖v‖=C控制I3。

因?yàn)镮1、I2、I3均可以通過‖v‖=C控制，并且由條件(C6)知F(x,t,h)<0，所以通過選擇足夠大的常數(shù)C，有Ξ(β,θ)<0，從而有Q(β,θ)0，式(4)成立。

因此存在一個(gè)局部最大化，使得

(8)

以上為定理1中結(jié)論(I)的證明，接下來證明結(jié)論(II)。

由于

4 帶寬選擇方法和估計(jì)算法

4.1 帶寬選擇

不失一般性，基于式(3)，假設(shè)誤差變量與Xi、Ti是獨(dú)立的，并且類似于文獻(xiàn)[12]中給出的最小二乘B樣條估計(jì)(LSB)的漸進(jìn)方差，這里給出BSMR估計(jì)量的漸進(jìn)方差的比率：

(9)

可通過下式來選擇h:

hopt=argminhr(h)=argminhG(h)F-2(h)

(10)

由式(10)可以知道，hopt僅由ε的條件誤差分布來決定。

需要指出的是，根據(jù)r(h)的表達(dá)式：當(dāng)h>0的時(shí)候，如果誤差服從標(biāo)準(zhǔn)正態(tài)分布，infhr(h)=1；如果不考慮誤差分布，infhr(h)≤1。因此，BSMR方法要優(yōu)于或至少不劣于LSB方法，尤其是當(dāng)誤差分布有重尾或者大方差的時(shí)候，BSMR的性能要比LSB的性能好很多。

在實(shí)際應(yīng)用中，若不知道誤差分布，則得不到G(h)和F(h)。通常用下式來估計(jì)G(h)和F(h)：

(11)

4.2 估計(jì)算法

步驟1(E-step)通過下式更新π(i/θ(l))：

π(i/θ(l))=

i=1,2,…,n

步驟2(M-step)更新γ：

=(WTZW)-1WTZY，

其中，W=(W1,W2,…,Wn)，Z=diag(π(1/γ(l)),…,π(n/γ(l)))，Y=(y1,y2,…,yn)。

由于該算法的收斂值可能會(huì)依賴于初始值，并且不能保證EM算法可以收斂到全局最優(yōu)解，因此需要對(duì)不同的初始值進(jìn)行計(jì)算,并從中選取局部最優(yōu)解。

另外，對(duì)于以上的估計(jì)過程，需要確定最優(yōu)的節(jié)點(diǎn)數(shù)K，本文選取最大化交叉驗(yàn)證函數(shù)的解作為最優(yōu)節(jié)點(diǎn)數(shù)，即

(12)

5 數(shù)值模擬和實(shí)例驗(yàn)證

5.1 蒙特卡洛模擬

對(duì)于非參數(shù)部分，使用均方誤的平方根(square root of average square errors, RASE)指標(biāo)來評(píng)價(jià)估計(jì)結(jié)果：

(13)

用文獻(xiàn)[11]中定義的廣義均方誤(generalized mean square error, GMSE)來評(píng)價(jià)參數(shù)部分的估計(jì)結(jié)果：

(14)

數(shù)值模擬結(jié)果如表1所示，為了檢驗(yàn)BSMR估計(jì)量的優(yōu)效性和穩(wěn)健性，表中還列出最小二乘B樣條估計(jì)(LSB)[12]的結(jié)果進(jìn)行對(duì)比分析。

從表1中可以看出，對(duì)于參數(shù)部分的估計(jì)，LSB和BSMR的估計(jì)誤差都很小，性能比較接近，而且隨著樣本量的增加，兩種估計(jì)的結(jié)果都會(huì)逐漸變好；對(duì)于非參數(shù)部分，當(dāng)誤差服從正態(tài)分布或者t分布時(shí)，BSMR估計(jì)誤差比LSB的小很多，當(dāng)誤差服從混合正態(tài)分布時(shí)，BSMR的估計(jì)結(jié)果也好于LSB，但是差距不是特別明顯?？偟膩砜矗珺SMR估計(jì)要優(yōu)于LSB估計(jì)。

5.2 實(shí)例分析

下面通過實(shí)例來驗(yàn)證BSMR估計(jì)方法的可行性。采用Nierenbrg等[16]收集的血漿中β胡蘿卜素的水平數(shù)據(jù)，該數(shù)據(jù)集有315個(gè)觀測(cè)值，本文研究血漿中β胡蘿卜素的水平與下列因素的關(guān)系：年齡，性別，吸煙狀況，克托萊指數(shù)，維生素服用情況，食物熱量(卡路里)，脂肪、膳食纖維、酒精飲料和膽固醇的攝入量。

應(yīng)用模型(1)進(jìn)行分析，其中T為年齡，協(xié)變量吸煙狀況和維生素服用情況是分類變量，重新將它們?cè)O(shè)置為虛擬變量；將以上的虛擬變量和離散變量(性別)以及酒精飲料攝入量作為參數(shù)部分的協(xié)變量。將年齡指標(biāo)T歸一化處理。

表2給出了模型的系數(shù)以及參數(shù)估計(jì)的標(biāo)準(zhǔn)差。由表2可見，BMSR的MAPE值比LSB的MAPE值小，即BSMR 的擬合效果優(yōu)于LSB，再次驗(yàn)證了本文方法的可行性及有效性。

表1 不同誤差分布下的模擬結(jié)果

表2 參數(shù)估計(jì)值及其MAPE