R语言检验经典线性回归方程具有一阶自相关并消除自相关的方法
2016-03-24 10:46阅读：
若Cov(ui , uj) = E(uiuj) = 0,
i≠j不成立，即线性回归模型扰动项的方差—协方差矩阵的非主对角线元素不全为0，则称为扰动项自相关，或序列相关。
自相关的后果与异方差类似：系数估计量的无偏性、一致性和渐近正态性仍成立。
(1) 在扰动项自相关的情况下，不再具有最小方差的性质, 即不是BLUE。
(2) RSS/(n-K-1)不是
的一致估计，OLS估计量的标准误差不再是标准差的一致估计量，
t检验统计量即使在大样本下也不再服从N(0,1)，F检验也不再服从F分布，LM检验不服从卡方分布。无法再信赖参数的置信区间或假设检验结果。
德宾—沃森检验法(Durbin－Watson d
test)
要检验的假设：
检验统计量：
R语言代码：
yhat<-lm(y~x)
yhat
plot(y~x);abline(y1)
ahat=summary(y1)$coef
![[1]]( "[1]")
;bhat=summary(y1)$coef
![[2]]( "[2]")
et=summary(y1)
$resi
E=rep(0,n-1)
for(t in 2:n){E[t-1]=et ![[t]]( "[t]")
-et[t-1]}
DW=sum(E^2)/sum(et^2)
DW
消除自相关的方法：
科克伦—奥克特法(Cochrane—Orcutt)
它是一个迭代过程，步骤如下：
①估计原模型，计算OLS残差et 。
②et对et-1回归，即估计et=ρet-1
εt，得到ρ的估计值phat
③用
phat
产生
Ytp=Yt-phat*Yt-1 ；Xtp=Xt-phat*Xt-1
然后估计
Ytp=ahatp bhatp*Xtp
εt ，得估计值
ahat 和 bhat ，以及残差et。
④重复第②步、第③ 步。
此过程不断更新phat ，ahat
和
bhat ，直至收敛。
R语言代码：
k=50;p=1
Allr=rep(0,k);Alla=rep(0,k);Allb=rep(0,k)
while(p<=k){
e=summary(y1)$resi
#获取et
e1=rep(0,n);e1
![[1]]( "[1]")
=0
for (i in 2:n){e1
![[i]]( "[i]")
=e[i-1]}
#获取et-1
E1=lm(e~e1)
phat=summary(E1)
}
y1=lm(Yp~Xp)
Allp
![[p]]( "[p]")
=phat;Alla
![[p]]( "[p]")
=summary(lm(Yp~Xp))'>$coef
![[2]]( "[2]")
Yp=rep(0,n);Yp
![[1]]( "[1]")
=y
![[1]]( "[1]")
for (i in 1:(n-1)){Yp[i 1]=y[i 1]-phat*y
![[i]]( "[i]")
}
Xp=rep(0,n);Xp
![[1]]( "[1]")
=x
![[1]]( "[1]")
for (i in 1:(n-1)){Xp[i 1]=x[i 1]-phat*x
![[i]]( "[i]")
}
y1=lm(Yp~Xp)
Allp
![[p]]( "[p]")
=phat;Alla
![[p]]( "[p]")
=summary(lm(Yp~Xp))$coef
![[1]]( "[1]")
;Allb
![[p]]( "[p]")
=summary(lm(Yp~Xp))$coef
![[2]]( "[2]")
print(summary(y1)$coef)
print(c(估计次数=p,p估计=phat,a估计=summary(lm(Yp~Xp))$coef
![[1]]( "[1]")
,b估计=summary(lm(Yp~Xp))$coef
![[2]]( "[2]")
))
p=p 1
}
#绘图看是否收敛
par(mfrow=c(2,2))
plot(Allp)
plot(Alla)
plot(Allb)
如上图所示，当循环40次以后，估计值就开始收敛到一个恒定的值上。用此时的估计值做系数就消除了自相关性。
