基于SAS的回歸分析異方差的檢驗(yàn)與消除

王純杰，　鄭　茜，　李　群，　林珊屹

(長春工業(yè)大學(xué) 基礎(chǔ)科學(xué)學(xué)院， 吉林 長春　130012)

?

基于SAS的回歸分析異方差的檢驗(yàn)與消除

王純杰，鄭茜*，李群，林珊屹

(長春工業(yè)大學(xué) 基礎(chǔ)科學(xué)學(xué)院， 吉林 長春130012)

摘要：利用SAS軟件的REG過程步實(shí)現(xiàn)了回歸分析中異方差的White檢驗(yàn)，利用加權(quán)最小二乘法消除了異方差的影響。通過與殘差圖法和等級相關(guān)系數(shù)法的結(jié)果比較，該檢驗(yàn)效果更精確。利用2013年全國31個(gè)省市自治區(qū)保險(xiǎn)賠付和保費(fèi)數(shù)據(jù)進(jìn)行實(shí)證分析。

關(guān)鍵詞：回歸分析； 異方差檢驗(yàn)； 加權(quán)最小二乘估計(jì)

0引言

根據(jù)回歸分析中最小二乘估計(jì)的結(jié)果，可以得到自變量與因變量之間的關(guān)系，但是模型中可能存在整體模型或者個(gè)別變量不顯著的情況，由多種原因?qū)е拢绠惓Ｖ蹬c強(qiáng)影響點(diǎn)，多重共線性等。異方差是常見原因之一，消除異方差可以使模型顯著。異方差問題分為兩大類，即異方差的診斷和消除。關(guān)于異方差診斷，常見的方法有殘差圖檢驗(yàn)法和斯皮爾曼等級相關(guān)系數(shù)檢驗(yàn)法。這兩種方法可以判定是否存在異方差，但是它們有一定的局限性，可見探索其他的方法很有意義。

1回歸方程和普通最小二乘估計(jì)[1]

用y表示因變量，X=(X1,X2,…,Xp)T是對y有影響的p個(gè)自變量，回歸模型一般表述為:

(1)

(2)

即

(3)

回歸模型的基本假設(shè)如下：

1)E(ε)=0，var(ε)=σ2I;

2)E(X′ε)=0;

3)rank(X)=p+1。

于是有E(y)=Xβ,cov(εi,εj)=σ2，線性回歸模型簡記為(y,Xβ,σ2In)[1]。

定義離差平方和

(4)

回歸系數(shù)β的估計(jì)為

回歸方程的估計(jì)式為:

(5)

根據(jù)模型假設(shè)，回歸模型要檢驗(yàn)3個(gè)方面：

1)回歸模型的顯著性；

2)每個(gè)自變量的顯著性；

3)殘差是否為“純隨機(jī)”。

2 異方差的檢驗(yàn)

2.1異方差產(chǎn)生的原因和影響

經(jīng)典回歸假設(shè)中一個(gè)重要的假設(shè)是var(εi)=σ2，即擾動項(xiàng)的同方差性。當(dāng)模型中隨機(jī)擾動項(xiàng)的方差不是常數(shù)時(shí)，稱為異方差或方差非齊性。記為：

(6)

隨機(jī)誤差項(xiàng)包含眾多因素對解釋變量的影響，如果它們隨著解釋變量觀測值的變化而變化，則會對解釋變量產(chǎn)生影響。異方差現(xiàn)象存在時(shí)，參數(shù)估計(jì)值雖然是無偏的，但不是有效的，即不再具有最小方差性。以一個(gè)簡單線性回歸模型為例，模型如下:

(7)

對回歸系數(shù)β1進(jìn)行普通最小二乘(OLS)估計(jì)，可得:

(8)

(9)

(10)

(11)

則

(12)

則

2.2異方差的檢驗(yàn)

關(guān)于異方差性的檢驗(yàn)[2]，統(tǒng)計(jì)前輩們進(jìn)行了大量的研究，提出的診斷方法很多，但是沒有一個(gè)公認(rèn)的最好的方法。下面先介紹一下殘差圖分析法和等級相關(guān)系數(shù)法，這是兩種比較常見的方法[2]。

2.2.1殘差圖分析法

殘差圖分析法非常簡單、直觀，它是一種非正式的檢驗(yàn)方法。如果在檢驗(yàn)前不能判斷是否存在異方差，可以采用該法。做法是先不考慮存在異方差的假定下構(gòu)造回歸模型,然后以殘差ei為縱坐標(biāo)，以其他合適的變量為橫坐標(biāo)，畫散點(diǎn)圖。常用的橫坐標(biāo)有3種：

1)以解釋變量xi(i=1,2,…,p)為橫坐標(biāo)；

3)以觀測值或序號為橫坐標(biāo)。

如果回歸模型適合樣本數(shù)據(jù)，那么殘差ei應(yīng)該反映εi所假定的性質(zhì)，即殘差圖上的n個(gè)點(diǎn)應(yīng)該是隨機(jī)分布和沒有任何規(guī)律的。如果回歸模型存在異方差，殘差圖上點(diǎn)的分布會呈現(xiàn)一定的趨勢，如殘差ei的值隨y值的增大而增大(或減小)，具有明顯的規(guī)律，由此可以判斷模型的隨機(jī)誤差項(xiàng)εi的方差是非齊性的，存在異方差。

2.2.2等級相關(guān)系數(shù)法

等級相關(guān)系數(shù)法又稱斯皮爾曼(Spearman)檢驗(yàn)，是一種應(yīng)用比較廣泛的方法。進(jìn)行等級相關(guān)系數(shù)檢驗(yàn)有3個(gè)步驟：

1)求出y關(guān)于x的普通最小二乘估計(jì)，求出εi的估計(jì)值，即ei的值。

2)取ei的絕對值，即|ei|，把xi和|ei|按遞增(或遞減)的次序排列，分成等級，按下式計(jì)算等級相關(guān)系數(shù):

(13)

式中:n----樣本容量;

di----對應(yīng)xi和|ei|的等級的差數(shù)。

3)做等級相關(guān)系數(shù)的顯著性檢驗(yàn)。在n>8的情況下，用下式對樣本等級相關(guān)系數(shù)rs做t檢驗(yàn)，檢驗(yàn)統(tǒng)計(jì)量為:

(14)

2.2.3White檢驗(yàn)

(15)

(16)式中：n----觀測值的個(gè)數(shù);

p----回歸參數(shù)的個(gè)數(shù)，包括常數(shù)項(xiàng)。

(17)

(18)

研究發(fā)現(xiàn)，當(dāng)樣本量小于250時(shí)，HC3的性質(zhì)最優(yōu)。

在SAS中,可以通過在model語句后面加入HCCMETHOD=0(,1,2,3)選項(xiàng)來選擇異方差性相容協(xié)方差矩陣的估計(jì)方法，默認(rèn)選項(xiàng)是HC0。

model語句后面的ACOV選項(xiàng)的作用是在參數(shù)估計(jì)的結(jié)果中顯示異方差性相容協(xié)方差矩陣，而且增加了異方差性標(biāo)準(zhǔn)差，就是White標(biāo)準(zhǔn)差。如果在model語句后面選擇了HCC選項(xiàng)或者WHITE選項(xiàng)，但是沒有選擇ACOV選項(xiàng)，結(jié)果中會顯示異方差性標(biāo)準(zhǔn)差，但是不會顯示異方差性相容協(xié)方差矩陣。SPEC選項(xiàng)是對模型做一個(gè)檢驗(yàn)，原假設(shè)是回歸估計(jì)的誤差項(xiàng),是獨(dú)立同方差的。備擇假設(shè)是誤差項(xiàng)具有異方差?？梢灾付╩odel語句后面的SPEC，ACOV，HCC，White選項(xiàng)，來選擇檢驗(yàn)的類型，檢驗(yàn)給出的一致協(xié)方差矩陣都是漸近的。ACOV選項(xiàng)和SPEC選項(xiàng)可以用在model語句或者print語句中。

3異方差的消除

如果模型存在異方差，就違背了回歸模型的基本假定，不能使用普通最小二乘法估計(jì)參數(shù)了。如果對原有的模型進(jìn)行變換，使其滿足同方差性的假設(shè)，再進(jìn)行參數(shù)估計(jì)，可以得到理想的模型。最常用的方法是加權(quán)最小二乘法[4]。

3.1加權(quán)最小二乘估計(jì)

假設(shè)觀測數(shù)據(jù)y與參數(shù)β1,β2，…，βp之間服從如下線性關(guān)系:

(19)

其中,ε是擾動項(xiàng)，進(jìn)行n次觀測可以得到n個(gè)線性方程:

(20)

式(20)可以用矩陣表達(dá)，即:

(21)

(22)

(23)

令式(23)為0，解得:

(24)

式(24)就是最小二乘估計(jì)的表達(dá)式，它與觀測數(shù)據(jù)y和系數(shù)陣X有關(guān)[5]。

(25)

(26)

解得:

(27)

3.2加權(quán)最小二乘估計(jì)權(quán)數(shù)的確定

加權(quán)最小二乘估計(jì)誤差的方差推導(dǎo)過程如下：

(28)

我們的目標(biāo)是求出使式(28)達(dá)到最小的w, 即加權(quán)系數(shù)，由此引入矩陣型許瓦茲不等式。

定理1設(shè)A,B分別為m×n和n×l矩陣，且AAT是可逆矩陣，則有：

(29)

當(dāng)且僅當(dāng)存在一個(gè)m×l維矩陣P，使得B=ATP時(shí)式(29)取等。

證明：

BTB-2BTAT(AAT)-1AB+BTAT(AAT)-1AAT(AAT)-1AB=

(30)

所以BTB≥(AB)T(AAT)-1AB得證。

當(dāng)且僅當(dāng)B=AT(AAT)-1AB時(shí),式(30)右端為零矩陣。令P=(AAT)-1(AB)，有B=ATP；反之,如果B=ATP，則有AT(AAT)-1AB=AT(AAT)-1AATP=ATP=B，由此證明了式(28)中等號成立的充要條件是B=ATP，證畢。

下面利用這個(gè)定理求加權(quán)系數(shù)。因?yàn)関ar(ε)=σε是對稱的正定矩陣,令σε=PTP,其中P為n×n的可逆矩陣。取l=n,令A(yù)=XTP-1,B=PwX(XTwX)-1,可以看出AB=1,所以對任意的加權(quán)矩陣可以得到如下不等式：

BTB≥

(31)

比較式(31)的兩邊，有

(32)

此時(shí)加權(quán)最小二乘估計(jì)的誤差方差達(dá)到最小,即：

(33)

此時(shí)的估計(jì)為：

(34)

3.3加權(quán)最小二乘估計(jì)的步驟

1)將回歸方程的各項(xiàng)除以比例因子Zi，得到新的回歸方程，同時(shí)滿足誤差項(xiàng)的同方差性；

2)對新回歸方程進(jìn)行最小二乘估計(jì)，求得回歸參數(shù)。

4實(shí)證分析

SAS系統(tǒng)(Statistics Analysis System)的重要組成部分和核心功能是統(tǒng)計(jì)分析，最新版本是9.4版。 SAS可以提供多個(gè)統(tǒng)計(jì)分析過程，每個(gè)過程均含有豐富的選項(xiàng)，滿足客戶的不同需求[6]。

若想對某一實(shí)際現(xiàn)象作回歸分析，首先要檢驗(yàn)數(shù)據(jù)是否存在異方差現(xiàn)象，然后選擇合適的方法進(jìn)行消除，最后觀察模型擬合數(shù)據(jù)的程度。2013年全國31個(gè)省市自治區(qū)保險(xiǎn)保費(fèi)收入和賠付支出情況見表1。

表1　2013年全國31個(gè)省市自治區(qū)保險(xiǎn)保費(fèi)收入和賠付支出情況

文中利用的數(shù)據(jù)是2013年全國31個(gè)省市自治區(qū)保險(xiǎn)賠付支出和保費(fèi)收入的數(shù)據(jù)，該數(shù)據(jù)是橫截面數(shù)據(jù)(來自中國統(tǒng)計(jì)年鑒，2014)。設(shè)保費(fèi)收入為自變量x，賠付支出為因變量y，單位是億元，由此建立回歸模型，分析自變量對因變量產(chǎn)生的影響。

在進(jìn)行回歸分析之前，需要先進(jìn)行相關(guān)分析，利用SAS軟件的CORR過程步可以實(shí)現(xiàn)。因變量支出與自變量收入的相關(guān)系數(shù)是0.989 75，在顯著性水平為0.05下，通過檢驗(yàn)，說明自變量與因變量之間有很大的相關(guān)性, 可以進(jìn)行回歸分析[7]。 SAS程序如下：

data expend2013;

input expend income;

cards;

318.17994.44

102.00276.80

…

24.0472.70

106.59273.49

;;

proc corr data=expend2013;

run;

proc reg data=expend2013;

model expend=income;

output out=expend2013_r r=residual;

run;

quit;

普通回歸分析結(jié)果見表2。

表2　方差分析表

因?yàn)閿?shù)據(jù)的單位小，每個(gè)變量的值很大，所以導(dǎo)致總平方和、離差平方和、回歸平方和都很大，F(xiàn)值為1 393.24，P值<0.000 1，表明回歸方程高度顯著，說明自變量對因變量有高度顯著影響。通過復(fù)決定系數(shù)R2=0.979 6，也可以判定回歸方程高度顯著。在顯著性水平為0.05以下，自變量的t檢驗(yàn)通過了，說明它對因變量有顯著影響，常數(shù)項(xiàng)也通過了檢驗(yàn)，下面開始檢驗(yàn)異方差[8]，SAS程序如下：

proc gplot data=expend2013_r;

plot residual*income=1;

symbol1 c=blue v=circle i=none;

label income="保費(fèi)" residual="殘差";

run;

殘差關(guān)于自變量收入的散點(diǎn)圖如圖1所示。

圖1殘差圖

可以看出殘差不是隨機(jī)分布的，有向外擴(kuò)張的趨勢，初步判定存在異方差現(xiàn)象。下面進(jìn)行等級相關(guān)系數(shù)檢驗(yàn)，程序如下：

data r2013_abs;

set expend2013_r;

abs_r=abs(residual);

run;

proc corr data=r2013_abs spearman;

var abs_r income;

label income="保費(fèi)" abs_r="殘差絕對值";

run;

等級相關(guān)系數(shù)檢驗(yàn)結(jié)果見表3。

表3通過SAS軟件得到自變量收入與殘差絕對值|ei|的相關(guān)系數(shù)是0.572 18，P值為0.000 8，自變量收入與殘差絕對值|ei|顯著相關(guān)，即存在異方差現(xiàn)象，與殘差圖的結(jié)果吻合。接下來進(jìn)行white檢驗(yàn)，程序如下：

proc reg data=expend2013;

model expend=income / white;

run;

quit;

表3　等級相關(guān)系數(shù)檢驗(yàn)結(jié)果

white檢驗(yàn)結(jié)果見表4。

表4　white檢驗(yàn)結(jié)果

由表4可以判斷存在異方差現(xiàn)象(P值<0.000 1)。

下面開始進(jìn)行加權(quán)最小二乘的處理，分別采用殘差絕對值和殘差平方項(xiàng)作為權(quán)重，經(jīng)比較殘差平方項(xiàng)作權(quán)重的效果更好。利用殘差平方項(xiàng)加權(quán)的SAS程序如下：

procregdata=resid;

modelsqresid=income;

outputout=v_weightsp=v_hat;

run;

quit;

datav_weights;

setv_weights;

v_weight=/v_hat;

labelv_weight="weightsusingsquaredresiduals";

run;

procregdata=v_weights;

weightv_weight;

modelexpend=income/noint;

outputout=e2013_r_vr=residual;

labelincome="保費(fèi)"expend="賠付";

run;

quit;

datar2013_abs_v;

sete2000_r_v;

abs_r_v=abs(residual2);

run;

proccorrdata=r2013_abs_vspearman;

varabs_r_vincome;

run;

加權(quán)之后，需要去掉常數(shù)項(xiàng)，加權(quán)之后的等級相關(guān)系數(shù)的結(jié)果是自變量收入與殘差絕對值|ei|的相關(guān)系數(shù)降為0.243 16，P值為0.111 75，可以判定異方差現(xiàn)象消除了。

表5　加權(quán)后的的參數(shù)估計(jì)結(jié)果

圖2回歸擬合圖

5結(jié)語

通過SAS軟件可以實(shí)現(xiàn)殘差圖法、等級相關(guān)系數(shù)法、white檢驗(yàn)，還可以實(shí)現(xiàn)加權(quán)最小二乘法消除異方差。SAS軟件功能強(qiáng)大，處理異方差的方法的還有很多，SAS軟件在這個(gè)方面的應(yīng)用還有待開發(fā)。

參考文獻(xiàn)：

[1]何曉群.應(yīng)用回歸分析[M].北京:中國人民大學(xué)出版社,2011:96-116.

[2]龔秀芳.幾種異方差檢驗(yàn)方法的比較[J].菏澤師范專科學(xué)校學(xué)報(bào),2003(25):19-22.

[3]WhiteH.Aheteroskedasticity-consistentcovariancematrixestimatorandadirecttestforheteroskedastivity[J].Econometrica,1980(48):817-838.

[4]楊波.加權(quán)最小二乘估計(jì)中加權(quán)系數(shù)的確定[J].現(xiàn)代電子技術(shù),2002(12):45-46.

[5]王純杰.基于分位數(shù)回歸的長春市職工工資水平的分析[J].長春工業(yè)大學(xué)學(xué)報(bào):自然科學(xué)版,2010,31(4):367-373.

[6]董大鈞.SAS統(tǒng)計(jì)分析應(yīng)用[M].北京：電子工業(yè)出版社,2008.

[7]朱世武.SAS編程技術(shù)教程[M].北京：清華大學(xué)出版社,2013.

[8]InstituteS.SAS/STATuser’sguide[M]. [S.l.]:SASInstitute,1990.

Inspection and elimination of heteroscedasticity based on SAS

WANG Chunjie,ZHENG Xi*,LI Qun,LIN Shanyi

(School of Basic Sciences, Changchun University of Technology, Changchun 130012, China)

Abstract：Heteroscedasticity phenomenon is tested with white test of REG procedure in SAS, and the weighted least square estimation is used to eliminate heteroscedasticity. The test is compared with both the residual graph and Spearman correlation coefficient method to show more reasonable results. The data of insurance and premium in 2013 is applied to the example analysis.

Key words：regression analysis; heteroscedasticity; weighted least square estimate.

收稿日期：2016-01-21

基金項(xiàng)目：國家自然科學(xué)基金資助項(xiàng)目(11301031，11571051)； 吉林省教育廳“十三五”規(guī)劃項(xiàng)目(2016317)

作者簡介：王純杰(1978-)，女，漢族，遼寧遼陽人，長春工業(yè)大學(xué)副教授，博士，主要從事數(shù)理統(tǒng)計(jì)方向研究,E-mail:wangchunjie@ccut.edu.cn. *通訊作者：鄭茜(1989-)，女，漢族，吉林遼源人，長春工業(yè)大學(xué)碩士研究生，主要從事數(shù)理統(tǒng)計(jì)方向研究,E-mail:13069008450@163.com.

DOI:10.15923/j.cnki.cn22-1382/t.2016.3.01

中圖分類號：O 212.1

文獻(xiàn)標(biāo)志碼：A

文章編號：1674-1374(2016)03-0209-08