利用swiss数据集进行多元线性回归研究(先查看各变量间的散点图)

使用瑞士数据集的多元线性回归研究

#先看变量之间的散点图

pairs(swiss, panel = panel.smooth, main = “swiss data”,

col = 3 + (swiss$Catholic > 50))

#使用所有变量创建多元线性回归

a=lm(生育力 ~ . , 数据 = 瑞士)

总结（一）

##

##调用：

##lm(公式 = 生育能力 ~ ., 数据 = 瑞士)

##

## 残差：

##Min1Q Median3QMax

## -15.274 -5.262 0.503 4.120 15.321

##

## 系数：

##估计标准。误差t值Pr(>|t|)

##（拦截）66.9152 10.7060 6.25 1.9e-07 ***

## 农业-0.17210.0703 -2.45 0.0187 *

##考试-0.25800.2539-1.02 0.3155

## 教育-0.87090.1830-4.76 2.4e-05 ***

## 天主教0.10410.0353 2.95 0.0052 **

## 婴儿。死亡率 1.07700.3817 2.82 0.0073 **

##—

## 意义。代码：0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1

##

## 残差标准误差：7.17 on 41 自由度

## 多个 R 平方：0.707，调整后的 R 平方：0.671

## F-statistic: 19.8 on 5 and 41 DF, p-value: 5.59e-10

#从结果来看，Education变量的p值没有一颗星，说明对模型极其不显着。

# R中提供了add1 drop1函数来增加或减少线性模型的变量

drop1(a)

##单词删除

##

##型号：

##生育能力~农业+考试+教育+天主教+

##Infant.Mortality

##Df Sq RSS AIC 之和

##2105 191

##农业1308 2413 195

##考试153 2158 190

## 教育11163 3268 209

## 天主教1448 2553 198

## 婴儿。死亡率 1409 2514 197

#从结果来看，去掉Education这个变量后，AIC最小，所以下一步可以去掉这个变量进行建模。

b=update(a,.~.-教育)

总结(b)

##

##调用：

##lm(公式=生育能力~农业+考试+天主教+

##婴儿。死亡率，数据 = 瑞士）

##

## 残差：

##Min1Q Median3QMax

## -23.919 -3.553 -0.649 6.596 14.177

##

## 系数：

##估计标准。误差t值Pr(>|t|)

##（拦截）59.6027 13.0425 4.57 4.2e-05 ***

##农业-0.04760.0803-0.59 0.55669

##考试-0.96800.2528-3.83 0.00042***

## 天主教0.02610.0384 0.68 0.50055

## 婴儿。死亡率 1.39600.4626 3.02 0.00431 **

##—

## 意义。代码：0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1

##

## 剩余标准误差：8.42 自由度上的 82

## 多个 R 平方：0.545，调整后的 R 平方：0.501

## F-statistic: 12.6 on 4 and 42 DF, p-value: 8.27e-07

#从下面的结果看，有两个变量不显着，R平方只有0.53，模型效果极差。需要进一步研究。

#好在R有一个step函数用r语言做多元线性回归用r语言做多元线性回归，可以根据AIC最小原理自动过滤模型的变量

b=step(a,direction=”backward”)

## 开始：AIC=190.7

##生育能力~农业+考试+教育+天主教+

##Infant.Mortality

##

##Df Sq RSS AIC 之和

## – 考试153 2158 190

##2105 191

## – 农业1308 2413 195

## – 婴儿。死亡率 1409 2514 197

## – 天主教1448 2553 198

## – 教育11163 3268 209

##

## 步骤：AIC=189.9

## 生育能力 ~ 农业 + 教育 + 天主教 + 婴儿。死亡率

##

##Df Sq RSS AIC 之和

## 2158 190

## – 农业1264 2422 193

## – 婴儿。死亡率 1410 2568 196

## – 天主教1957 3115 205

## – 教育12250 4408 221

总结(b)

##

##调用：

##lm(公式 = 生育能力 ~ 农业 + 教育 + 天主教 +

##婴儿。死亡率，数据 = 瑞士）

##

## 残差：

##Min1Q Median3QMax

## -14.676 -6.052 0.751 3.166 16.142

##

## 系数：

##估计标准。误差t值Pr(>|t|)

##（拦截）62.10139.6049 6.47 8.5e-08 ***

## 农业-0.15460.0682 -2.27 0.0286 *

## 教育-0.98030.1481-6.62 5.1e-08 ***

## 天主教0.12470.0289 4.31 9.5e-05 ***

## 婴儿。死亡率 1.07840.3819 2.82 0.0072 **

##—

## 意义。代码：0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1

##

## 残差标准误差：7.17 on 42 自由度

## 多重 R 平方：0.699，调整后的 R 平方：0.671

## F-statistic: 24.4 on 4 and 42 DF, p-value: 1.72e-10

接下来，建模变量和模型的回归诊断研究

首先，对自变量进行正态性检验

shapiro.test(swiss$Agriculture)

##

## Shapiro-Wilk 正态性检验

##

## 数据：瑞士$农业

## W = 0.9664，p 值 = 0.193

shapiro.test(swiss$Examination)

##

## Shapiro-Wilk 正态性检验

##

## 数据：swiss$Examination

## W = 0.9696，p 值 = 0.2563

shapiro.test(swiss$Education)

##

## Shapiro-Wilk 正态性检验

##

## 数据：swiss$Education

## W = 0.7482，p 值 = 1.312e-07

shapiro.test(swiss$Catholic)

##

## Shapiro-Wilk 正态性检验

##

## 数据：swiss$Catholic

## W = 0.7463，p 值 = 1.205e-07

shapiro.test(swiss$Infant.Mortality)

##

## Shapiro-Wilk 正态性检验

##

## 数据：swiss$Infant。死亡率

## W = 0.9776，p 值 = 0.4978

从各个变量的正态性检验结果来看，变量Education和Catholic的p值均小于0.05，所以这两个变量的数据不符合正态分布.

现在，模型的残差也进行了正态性检验（回归模型的残差也应该符合正态分布）

b.res

shapiro.test(b.res)

##

## Shapiro-Wilk 正态性检验

##

## 数据：b.res

## W = 0.9766，p 值 = 0.459

从结果看，p值为0.459，模型残差符合正态分布

接下来，绘制回归值和残差的残差图（应该符合均匀分布，即无论回归值如何，残差都具有相同的分布）

par(mfrow=c(1,2))

#绘制残差图

绘图(b.res~predict(b))

#绘制标准残差

绘图(rstandard(b)~predict(b))

par(mfrow=c(1,1))

从残差图看，效果不明显。

其实可以直接画残差图

par(mfrow=c(2,2))

情节(b)

par(mfrow=c(1,1))

结束。