崔 樂, 吳 迪, 成麗波
(長春理工大學(xué) 理學(xué)院, 長春 130022)
近年來,異常值檢測一直是統(tǒng)計診斷中一個比較活躍的研究課題。在實際生活中,由于種種原因不能夠和其他數(shù)據(jù)一起用于多元統(tǒng)計的線性回歸模型中的值,稱為異常值。換句話說,參數(shù)估計、模型構(gòu)建及預(yù)測都會受到異常值的影響,因此,這篇文章主要研究的是基于逐步回歸的穩(wěn)健估計和異常值的檢測。
比較常見的統(tǒng)計量有殘差、學(xué)生化殘差、Cook距離以及W-K統(tǒng)計量等。殘差和學(xué)生化殘差[1-3]都可以用來檢測異常值,但是殘差沒有考慮到異方差性,而學(xué)生化殘差考慮到了這一點,因此用學(xué)生化殘差檢測異常值比用普通殘差更有效。而Cook距離和W-K統(tǒng)計量[4-6]僅可以判斷出數(shù)據(jù)是否有強影響點,但是異常點和強影響點之間沒有必然的關(guān)系,即異常點不一定是強影響點,強影響點也不一定是異常點。
本文基于逐步回歸模型,從殘差平方和的角度出發(fā),研究出檢測異常值的方法,并估計異常值的大小。通過與傳統(tǒng)方法的比較,驗證該方法的有效性。
在多元線性回歸模型中,自變量的選擇實質(zhì)上就是模型的選擇?,F(xiàn)設(shè)一切可供選擇的變量是t個,它們組成的回歸模型稱為全模型(記m=t+1),在獲得n組觀測數(shù)據(jù)后,有模型
其中:Y是n×1的觀測值,β是m×1未知參數(shù)向量,X是n×m結(jié)構(gòu)矩陣,并假定X的秩為m。
β=(βp,βq)′,X=(Xp?Xq)
下面從不同的角度給出自變量選擇的準則[7]。
準則1 平均殘差平方和達到最小
準則2CP統(tǒng)計量達到最小
該準則是由馬斯洛于1964年從預(yù)測的角度提出:
準則3 AIC準則
該準則由日本統(tǒng)計學(xué)家赤池弘次于1974年根據(jù)極大似然估計原理提出:
AIC=nln(SSEP)+2P
其中:SSEm是m個自變量x1,x2,…,xm所對應(yīng)的殘差平方和;SSEP是p個自變量x1,x2,…,xp所對應(yīng)的殘差平方和。
快速選擇變量的算法有很多,例如:向前法、向后法、逐步回歸法[7]等,其中逐步回歸法是應(yīng)用最廣的一種方法。它的具體做法是先將變量一個一個的引入,當引入到第3個變量之后的每一步,首先對已引入的變量進行剔除。這樣,自變量將不斷的引入、剔除、再引入、再剔除……直到自變量不能被剔除,同時也無法引入自變量為止。
經(jīng)過逐步回歸后生成最優(yōu)多元線性回歸模型如下:
yi=β0+β1xi1+β2xi2+…+βpxip+εi,i=1,2,…,n
(1)
其中:p為解釋變量的數(shù)目;βj(j=0,1,…,p)為回歸系數(shù);εi為隨機誤差,且εi~N(0,σ2),i=1,2,…,n;β=(β0,β1,…βp)T;σ為未知參數(shù)。
均值漂移模型[8-9]是在第i個數(shù)據(jù)點上增加一個漂移項δ,即在這個數(shù)據(jù)點yi處的均值發(fā)生了非隨機漂移,若δ顯著不等于零,則yi處的均值發(fā)生了漂移,說明此點為異常點。
由于事先不知道在線性模型中出現(xiàn)異常值,因此,可以先假定模型中沒有異常值,其線性模型的矩陣形式為
Y=Xβ+ε
(2)
其中
SSE=YT(I-H(X))Y
殘差平方和(SSE)會受到異常值大小的影響,即隨著異常值的增加,殘差平方和(SSE)也增加。
基于均值漂移模型,假設(shè)模型中的第i個觀測值為異常值,即異常值的大小為δi,其余的皆為正常數(shù)據(jù),計算此時的殘差平方和(SSEi)為
(3)
(4)
(5)
M估計是基于最小二乘估計發(fā)展起來的一種抗差估計(Robust Estimation)方法[10-12],是由huber于1964年最先提出來的,也稱為廣義最大似然估計。M估計已經(jīng)成為最經(jīng)典的一種穩(wěn)健估計方法。
M估計的估計方程寫成矩陣形式是這樣的:
XTWXβ=XTWY
(6)
迭代公式如下:
(7)
其中:W是以ωi,i=1,2,…,n為對角線的權(quán)矩陣;X是解釋變量矩陣,X=(x1,x2,…,xn);Y是因變量向量,Y=(y1,y2,…yn)T。
假設(shè)有一多元線性回歸模型y=3+5x1-4x2+4x3-3x4+x5-2x6+2x7-x8+3x9+6x10+ε,現(xiàn)利用計算機模擬產(chǎn)生100個數(shù)據(jù)樣本。下面考慮3種方案進行實驗。
方案1 將大小為5.5的異常值加入到第49個樣本觀測值中。
方案2 將大小為6,-5的異常值分別加入到第31和69個樣本觀測值中。
方案3 將大小為-6.5, 5.5, 4.5, -3的異常值分別加入到第48~52個樣本觀測值中。
1) 異常值檢測
對于3種方案的D統(tǒng)計量如圖1所示。圖1是基于M估計的D統(tǒng)計量,由圖1(a)可以看出,第50個觀測值對應(yīng)的D統(tǒng)計量遠大于3,判定為異常值。同樣由圖1(b)可以看出,第32個和第70個觀測值為異常值;由圖1(c)或者看出,從第49~53個觀測值均被檢測出異常值。由此可以看出該檢測異常值的統(tǒng)計量可以很好地檢測出異常值,并且由于方案3是5個連續(xù)的異常值,因此說明此方法對異常值的遮蔽現(xiàn)象有一定的作用,能夠有效地檢測出連續(xù)幾個異常值。
圖1 方案1~3的基于M估計的D統(tǒng)計量Fig.1 D Statistics Based on M Estimation for Schemes 1~3
圖2是方案1~3中每個樣本觀測值的學(xué)生化殘差示意圖。按照傳統(tǒng)的方法,將觀測值的學(xué)生化殘差的絕對值大于3的認為是異常值。從圖2(a)中可以看出,無法準確地檢測出方案1中異常值,圖2(b)中,可以很好地檢測出方案2中的異常值,從圖2(c)中,學(xué)生化殘差只能檢測出第49, 50, 52, 53個觀測值為異常值。由此可以說明,用新構(gòu)建的基于穩(wěn)健估計的D統(tǒng)計量比傳統(tǒng)的學(xué)生化殘差檢測異常值更加有效。
圖2 方案1~3的學(xué)生化殘差Fig.2 Student Residual in Scheme 1~3
2) 異常值大小的估計
通過公式(4)來計算異常值的大小,表1是3種不同方案的異常值大小的評估值。從表中可以看出,對于方案1~2,本文構(gòu)建的統(tǒng)計量都可以很好地檢測出異常值的大小。對于方案3異常值大小估計精確度不如方案1~2,但仍然還是比較準確的。
表1 基于方案1~3的異常值大小估計值Tab.1 Outlier size estimates based on Schemes 1~3
在房地產(chǎn)行業(yè)中,影響房屋價值的因素有很多。本實例中,根據(jù)房屋價值的影響因素及人們的偏好,記錄的影響房屋價值的指標包括面積、戶型、當前樓層、總樓層、朝向、裝修情況。其中,戶型為幾室?guī)讖d,可以分開作為單獨變量;當前樓層分為底層、中層和高層;朝向分為北、東北、東、東南、南、西南、南北、西、西北、東西;裝修情況分為毛坯、普通裝修、精裝修。
房價和房屋面積需要以10為底取對數(shù),目的是降低房價和面積的數(shù)量級,提高模型估計精度。戶型是一個連續(xù)的數(shù)值變量,可以不用修改。其他的字符型變量則需要進行量化。樓層包括3種,根據(jù)計量經(jīng)濟學(xué)知識,本方案將采用2個0-1變量(構(gòu)成一個二維行向量)度量該3種類型。例如,(0,0)代表低層,(1,0)代表中層,(0,1)代表高層。其他分類指標同樣需要進行類似處理。本文選取河北省邯鄲市某小區(qū)的樣本數(shù)據(jù)。
表2 河北省邯鄲市某小區(qū)的樣本數(shù)據(jù)Tab.2 Sample data of a residential district in Handan city, Hebei province
該實例的線性回歸模型為
經(jīng)過逐步回歸之后所得到的“最優(yōu)”回歸模型為
圖3 基于M估計的D統(tǒng)計量Fig.3 D StatisticsbBased on M estimation
1) 異常值檢測
從圖3可以看出,第24,38,42,43個樣本的對應(yīng)的D統(tǒng)計量大于3,可判斷為異常值,其中雖然第42和43的樣本連續(xù)異常,但仍然能夠被很好地檢測出來,再一次證明此方法對異常值的遮蔽現(xiàn)象有一定的作用,能夠有效地將連續(xù)幾個異常值檢測出來。
2) 異常值大小的估計
同樣,由公式(4)估計出異常值的大小,估計結(jié)果如表3所示。
表3 異常值大小的估計值Tab.3 Estimates of the size of outliers
異常值檢測是當前數(shù)據(jù)分析研究中的一個熱點問題。通過大量的模擬實驗和實例分析,得到這樣的結(jié)論:基于M估計的D統(tǒng)計量可以很好地檢測出異常值,尤其是對異常值的遮蔽現(xiàn)象有一定的作用----能夠有效地檢測出連續(xù)的幾個異常值。