半?yún)?shù)回歸模型的廣義交互改進

半?yún)?shù)回歸模型的廣義交互改進

在醫(yī)學(xué)科學(xué)研究數(shù)據(jù)的統(tǒng)計分析過程中，參數(shù)回歸模型的一些假設(shè)是不完全滿足的。例如，反應(yīng)變量與解釋變量之間的具體依存關(guān)系不明確，反應(yīng)變量的分布難以確定。此時,參數(shù)回歸模型難以進行擬合處理,而非參數(shù)回歸模型則能進行有效的分析。簡單的非參數(shù)回歸模型研究的是反應(yīng)變量Y與單一解釋變量t的依存關(guān)系,它能夠解決醫(yī)學(xué)與衛(wèi)生研究工作中的許多重要問題,但是,在實際工作中,有許多事物或現(xiàn)象受多個變量的影響,因此,需要研究多個變量間的相互關(guān)系。經(jīng)典統(tǒng)計模型在研究受多個解釋變量影響的依存關(guān)系時常常采用多重回歸,而多重回歸的更一般模型即為線性模型:yi=x′iβ+εi,為了放寬該線性模型中的某一個解釋變量的線性假定,使模型在假定方面具有較強的適應(yīng)性,本文對半?yún)?shù)回歸模型進行了研究?；貧w系數(shù)向量的估計假定對每一個觀察值yi,有p+1個解釋變量,其中p維向量xi和數(shù)量變量t,如果反應(yīng)變量y線性相關(guān)于解釋變量x,則有以下模型yi=x′iβ+g(ti)+εi(1)yi=x′iβ+g(ti)+εi(1)其中β為未知的p維回歸系數(shù)向量,g(t)為未知的光滑函數(shù)(如光滑樣條),x為線性變量,t為樣條變量,ε與(x,t)相互獨立,且E(ε)=0,V(ε)=σ2(未知),顯然,xi不含常數(shù)1,常數(shù)項可以包含在g(t)中,則以上模型被稱為半?yún)?shù)回歸模型(semiparametricregressionmodel)。半?yún)?shù)回歸模型可通過懲罰最小二乘方法進行求解,β和g(t)的估計使得以下加權(quán)懲罰平方和最小Sw(β,g)=n∑i=1wi{yi-x′iβ-g(ti)}2+α∫g″(t)2dt(2)Sw(β,g)=∑i=1nwi{yi?x′iβ?g(ti)}2+α∫g′′(t)2dt(2)其中光滑參數(shù)α>0,wi>0,不加權(quán)時,可令wi=1。令Y=(y1,…,yn)′,W=diag(w1,…,wn),X為n×p階矩陣,其第i行為x′i,為了考慮相持情況,假定t1,t2,…,tn可由s1,s2,…,sq來表示,表示它們之間關(guān)系的矩陣叫關(guān)聯(lián)矩陣(incidencematrix),用N來表示,N為n×q階矩陣,其元素為Nij,當(dāng)ti=sj時,Nij=1,否則,Nij=0。假定點ti不全相同,則q≥2。令αj=g(sj),j=1,2,…,q。則待估計向量g為(α1,α2,…,αq)′。同理,假定s1<s2<…<sq,而且αj=g(sj),則可以定義兩個矩陣Q和R,只不過要用s1,s2,…,sq來代替t1,t2,…,tn。令K=QR-1Qt,則∫g″(s)2ds=g′Kg。若用矩陣符號來表示Sw(β,g),則Sw(β,g)=(Y-Xβ-Ng)′W(Y-Xβ-Ng)+αg′Kg(3)當(dāng)β和g為以下分塊矩陣方程的解時,上式取最小值。[X′WXX′WΝΝ′WXΝ′WΝ+αΚ](βg)=[X′Ν′]WY(4)方程(4)是一個(p+q)元方程組,直接解方程組不方便,也很不實際,實際工作中,一般將方程(4)化為以下形式X′WXβ=X′W(Y-Ng)(5)(N′WN+αK)g=N′W(Y-Xβ)(6)求解時可采用不需迭代的直接法(directmethod)進行求解。由(6)可得:Ng=S(Y-Xβ)(7)其中S=N(N′WN+αK)-1N′W,(Ng)i=g(ti)。將(7)代入(5),化簡得X′W(I-S)Xβ=X′W(I-S)Y(8)這是廣義最小二乘正規(guī)方程組,用來估計β,加權(quán)矩陣為非對角陣W(I-S),解得β后,就可通過(7)求得g和Ng,因此,可得到光滑曲線g(t)。對于回歸系數(shù)向量β的估計值,可進行假設(shè)檢驗,β=(β1,β2,…,βp)′。檢驗假設(shè)為H0:βi=0,i=1,2,…,p備擇假設(shè)為H1:βi≠0,α=0.05檢驗統(tǒng)計量為t=?βi√Cii?σ2(9)其中Cii表示(X′W(I-S)X)-1的對角線上第i個元素,?σ2=n∑i=1(yi-?yi)2tr{Ι-A}?A為帽子陣。A=S+(I-S)X{X′W(I-S)X}-1X′W(I-S)(10)當(dāng)H0成立時,t～tυ,υ=tr{I-A}。在半?yún)?shù)回歸模型中,對于光滑參數(shù)的自動選擇需要計算廣義交互有效GCV(generalizedcross-validation)得分函數(shù)。GCV得分函數(shù)為GCV(α)=n∑i=1wi(yi-?yi)2(1-n-1trA)2(11)其中trA=trS+tr[{X′W(I-S)X}-1X′W(I-S)2X]。另外,半?yún)?shù)模型的誤差自由度EDF=tr{I-A}=n-trA,均方差MSE=n∑i=1(yi-?yi)2tr{Ι-A},殘差平方和SSE=n∑i=1(yi-?yi)2,令ˉy=1nn∑i=1yi,則擬合優(yōu)度R2=1-SSEn∑i=1(yi-ˉy)2。本文利用6.11版SAS軟件的IML模塊進行編程來實現(xiàn)以上分析過程?；貧w模型的建立為說明半?yún)?shù)模型的擬合效果,本文用SAS程序進行模擬抽樣實驗,取p=2,n=60,t由1變化到60,x1～N(12.66,2.572),x2～N(6.7,1.872),誤差項ε相互獨立且服從分布N(0.52),y=3.4x1-5.2x2+0.1(t-30)2+30.2+ε,則用SAS模擬抽樣程序可得到一個樣本模擬數(shù)據(jù)(表1)。如果假定y與x1,x2存在線性依存關(guān)系,對該數(shù)據(jù)人為地進行參數(shù)線性模型擬合,則可以得到回歸方程:?y=49.0545+0.1282t+4.4925x1-6.0078x2,雖然該回歸方程有意義(P≈0.0005),但擬合效果差,SSE=45494.6052,R2=0.2692,誤差均方為812.4037,從下面的圖1可知,殘差與t之間存在二次曲線趨勢,即殘差中仍然蘊含有用的回歸信息。如果采用半?yún)?shù)回歸模型進行擬合,則計算得到的α值為148.75,x1和x2的回歸系數(shù)分別為3.7976和-5.2356,標準誤分別為0.2385和0.2958,檢驗結(jié)果均有顯著意義(P<0.01),SSE=980.6252,MSE=19.2357,R2=0.9842,模型擬合的殘差情況見圖2,由上述計算結(jié)果和圖2可以看出,半?yún)?shù)模型的擬合效果得到大大提高,并且正確地反映了y與t的關(guān)系?；緮?shù)學(xué)模型的基本思想半?yún)?shù)回歸模型可看作是參數(shù)線性模型和非參數(shù)回歸模型的混合模型,半?yún)?shù)回歸模型較參數(shù)線性模型有較強的適應(yīng)性。由于實際工作中經(jīng)常會遇到某個變量有影響,但表現(xiàn)為未知函數(shù)的情況,因此,半?yún)?shù)回歸模型是線性模型的一個擴展,它放寬了線性模型中的某一個解釋變量的線性假定,使模型適應(yīng)數(shù)據(jù)變化的能力更強。實際應(yīng)用半?yún)?shù)回歸模型時,反應(yīng)變量線性相關(guān)于線性變量應(yīng)以專業(yè)理論知識或以往經(jīng)驗為依據(jù),樣條變量t的處理不同于其他線性變量,它是采用非參數(shù)的形式進行處理。方程(4)是一個(p+q)元方程組,多元方程組的解法很多,但直接解方程組不方便,也很不實際,實際工作中,也可采用backfitting方法求解方程組,backfitting是一個迭代求解的方法,它在上述二個方程(5)和(6)之間交替迭代求解,直至收斂為止。該方法的收斂速度取決于α的大小以及有關(guān)矩陣的特征值的大小,雖然該矩陣特征值的絕對值都小于1,最終也會收斂,但實際應(yīng)用中,經(jīng)常發(fā)生最大特征值很接近于1,從