R语言与多元线性回归方程及各种检验

R语言与多元线性回归方程及各种检验

---

一、模型建立

数据放在最后，自取试建立y与x1,x2,x3,x4,x5的线性回归方程

代码如下：

a<-read.csv("eg2.1.csv",header=T);a
b<-lm(y~x1+x2+x3+x4+x5,data=a)

这里介绍一个包，可以把结果用表格显示：

install.packages("flextable")
library(flextable)#制作表格
as_flextable(b)

输出结果如下：

得到回归方程

$$y=-32+0.163x_1+0.228x_2+0.881x_3-0.05x_4+0.169x_5$$

二、多重共线性

（1）产生的背景：

1.经济变量间具有共同变化的趋势。例如，对于时间序列数据收入、消费、就业率等等，在经济上升时期均呈现增长的趋势，而在经济收缩期，又都呈现下降趋势。当这些变量同时作为解释变量进入模型时就会带来多重共线性问题。
2.模型中引入了滞后变量，变量X常常与滞后变量高度相关
3.利用截面数据建立模型也可能会出现多重共线性。利用截面数据建模时，许多变量变化与发展规模有关，会出现同步增长的趋势，如资本、劳动力、科技、能源等投入与产出的规模相关，这时容易产生多重共线性。
4.样本数据自身的原因。例如，抽样仅仅限于总体中解释变量取值的一个有限范围，使得变量间变异不大；或由于总体受限，多个解释变量的样本数据相关，这时都可能出现多重共线性。

（2）多重共线性的检验

1.简单相关系数法：

一般而言，如果每两个解释变量的简单的相关系数(0阶自相关系数)比较高，如大于0.8，则可认为存在严重的多重共线性。

ab<-a[,c(2:6)]#只保留解释变量
corr<-cor(ab);corr
corrplot(corr,method="number",type="full",mar=c(0,0,0,0) ,bg="black",tl.col = "blue")

其中corrplot函数详细见这篇文章：
相关系数热力图

由上图可知，x1和x3的相关系数为0.73，初步判断不存在多重共线性。

2.方差膨胀因子（vif）法

一般vif大于10即可认为存在严重的多重共线性：
通过建立x1~x2+x3+x4+x5的线性关系得到R方，那么x1的方差膨胀因子就是
$$VIF=\frac{1}{1-R^2}$$
其他的类似，r语言有vif这个函数可以直接求：

library(car)#检验多重共线性的包
vif(b)#b是线性回归的结果

可以看出：vif值都小于10，即可判断出没有多重共线性

3.矩阵$X^{T}X$的条件数k

条件数的定义为：
$$k\left(X^{T}X\right)=\Vert X^{T}X\Vert \cdot \Vert {\left(X^{T}X\right)}^{-1}\Vert =\frac{{\lambda }_{\max \left(X^{T}X\right)}}{\lambda \min \left(X^{T}X\right)}$$
若k<100,则认为多重共线性的程度很小，在100到1000之间则认为存在中度或较强的多重共线性，若大于1000，则认为存在严重的多重共线性。
在R软件中，用kappa（）计算矩阵的条件数，其使用方法为：

**kappa(z,exact=F)***

其中z是矩阵，exact=F即不精确计算条件数，近似计算

XX<-cor(ab)#先做出相关系数矩阵
kappa(XX)#对相关系数矩阵做
eigen(XX)#求特征值和特征向量

可以看出条件数的值为7，即不存在多重共线性。

（3）多重共线性的修正

本文虽然没有出现多重共线性的问题，但是如果出现了，该如何解决呢？
本文主要通过逐步回归来修正多重共线性

代码如下：

aaa<-step(b)
as_flextable(aaa)

得到最后的回归方程：

$$y=-18.43+0.249x_2+0.986x_3$$

二、异方差性的检验及修正

1.异方差性的实质

设模型为
$${Y}_i={\beta }_0+{\beta }_1{x}_{1i}+{\beta }_2{x}_{2i}+...+{\beta }_k{x}_{{}_{ki}}+u_i\left(i=1,2,...n\right)$$
如果其他假定均不变，但模型中随机误差项$u_i$的方差为
$$Var(u_i)={\sigma }_i^2(i=1,2,...n)$$
则称$u_i$具有异方差性（heteroscedasiticity）

2.异方差性的检验

1.图示检验法

一般把标准化残差的绝对值大于等于2的观测值认为是可疑点，而把标准化残差值绝对值大于等于3认为是异常值

rst <- rstandard(b)#标准化残差
fit<-predict(b)#预测值
library(ggplot2)#载入ggplot2包
mn<-data.frame(rst,fit)
ggplot(mn,aes(fit,rst))+geom_point()

通过绘制出散点图可以看出：

标准化残差值大致在区间（-2,3）内 ，存在可疑点。是否存在异方差性还需要进一步检测。

2.Goldfeld—Quandt检验

install.packages("lmtest")
library(lmtest)
as_flextable(gqtest(b))

可以看出p值为0.4，即接受原假设，即不存在异方差

3.White检验和H.glesjser检验

White检验的统计量：

$$w=n\ast R^{2}ChiSquare\left(k\right)$$
由于解释变量由5个，这个其实不太适合做white检验
可以参考这篇文章：
异方差与r语言实践

3.异方差的修正

加权最小二乘法：在r语言中代码非常简单
这里取权重为1/abs（e）

e<-resid(b)#e是标准化残差
b1<-lm(y~x1+x2+x3+x4+x5,weights=1/abs(e),data=a)
as_flextable(b1)

得到方程：
$$\frac{y}{\sqrt{e}}=\frac{-36}{\sqrt{e}}+0.187\frac{x_1}{\sqrt{e}}+0.205\frac{x_2}{\sqrt{e}}+0.883\frac{x_3}{\sqrt{e}}+0.004\frac{x_4}{\sqrt{e}}+0.158\frac{x_5}{\sqrt{e}}$$
其中e是残差

数据链接
链接：https://pan.baidu.com/s/1My6WTMXEZyIVfIM7A5oIbA?pwd=tjby
提取码：tjby
–来自百度网盘超级会员V4的分享