R语言普通最小二乘回归分析

回归分析都是统计学的核心。它其实是一个广义的概念，通指那些用一个或多个预测变量（也称自变量或解释变量）来预测响应变量（也称因变量、效标变量或结果变量）的方法。普通最小二乘（OLS）回归 包括简单线性回归、多项式回归和多元线性回归。

OLS回归的数据必须满足以下统计假设:

(1)正态性:对于固定的自变量值，因变量值成正态分布。
(2)独立性: Yi值之间相互独立。
(3)线性:因变量与自变量之间为线性相关。
(4)同方差性:因变量的方差不随自变量的水平不同而变化。

如果违背了以上假设，你的统计显著性检验结果和所得的置信区间很可能就不精确。

1. 模型的拟合

### 1. 模型的拟合
## 简单线性回归
# ctl <- c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
# trt <- c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
# group <- gl(2, 10, 20, labels = c("Ctl","Trt"))
# # gl:Generate Factor Levels
# weight <- c(ctl, trt)
# lm.D9 <- lm(weight ~ group)
# lm.D90 <- lm(weight ~ group - 1) # omitting intercept

fit <- lm(weight ~ height,data = women)
# fit <- lm(weight ~ height - 1,data = women) # omitting intercept

## 多项式回归
fit2 <- lm(weight ~ height + I(height^2),data = women)

## 多元线性回归
fit3 <- lm(mpg ~ hp + wt + cyl, data = mtcars)

## 有显著交互项的多元线性回归
fit4 <- lm(mpg ~ hp + wt + hp:wt, data = mtcars)

## R表达式中常用的符号
# ~ 分隔符号，左边为响应变量，右边为解释变量。例如，要通过x、z和w预测y，代码为y ~ x + z + w
# + 分隔预测变量
# : 表示预测变量的交互项。例如，要通过x、z及x与z的交互项预测y，代码为y ~ x + z + x:z
# * 表示所有可能交互项的简洁方式。代码y~ x * z * w可展开为y ~ x + z + w + x:z + x:w + z:w + 
#   x:z:w
# ^ 表示交互项达到某个次数。代码y ~ (x + z + w)^2可展开为y ~ x + z + w + x:z + x:w + z:w
# . 表示包含除因变量外的所有变量。例如，若一个数据框包含变量x、y、z和w，代码y ~ .可展开为y ~ x + 
#   z + w
# - 减号，表示从等式中移除某个变量。例如，y ~ (x + z + w)^2 – x:w可展开为y ~ x + z + w + 
#   x:z + z:w
# -1 删除截距项。例如，表达式y ~ x - 1拟合y在x上的回归，并强制直线通过原点
# I() 从算术的角度来解释括号中的元素。例如，y ~ x + (z + w)^2将展开为y ~ x + z + w + z:w。相反, 代码y 
# ~ x + I((z + w)^2)将展开为y ~ x + h, h是一个由z和w的平方和创建的新变量
# function 可以在表达式中用的数学函数。例如，log(y) ~ x + z + w表示通过x、z和w来预测log(y)

summary(fit) #模型的详细结果
# R平方值，表明模型可以解释响应变量方差的百分比，是实际和预测值之间的相关系数
# 残差标准误（Residual standard error）则可认为是模型用身高预测体重的平均误差。
# F统计量检验所有的预测变量预测响应变量是否都在某个几率水平之上

coefficients(fit) #模型参数
confint(fit) # 模型参数的置信区间（95%）
fitted(fit) # 模型预测值，实际值women$weight
residuals(fit) # 残差值，weight-fitted(fit)
predict(fit,women) # 预测,women必需为含有fit中自变量的dataframe

anova(fit) #拟合模型的方差分析表
vcov(fit)  # 模型参数的协方差矩阵

2. 模型可视化与模型评价

### 2. 模型可视化与模型评价
plot(women$height,women$weight,
     xlab="Height (in inches)",
     ylab="Weight (in pounds)")
# Add Straight Lines to a Plot
abline(fit) 
# Add Connected Line Segments to a Plot
lines(women$height,fitted(fit)) 

par(mfrow = c(2,2))
# 生成评价模型拟合情况的四幅图形
plot(fit)
# 残差图与拟合图（左上）：残差值与预测（拟合）值应该没有任何系统关联
# 正态Q-Q图（右上）：应该落在呈45度角的直线上
# 位置尺度图（左下）：水平线周围的点应该随机分布
# 残差与杠杆图（右下）：从图形可以鉴别出离群点、高杠杆值点和强影响点

3. car包评价回归模型

### 3. car包评价回归模型
library(car)
states<- as.data.frame(state.x77[,1:7])
fit <- lm(Murder ~ Population + Illiteracy + Income + Frost,
          data = states)
## 正态假设检验
qqPlot(fit,lables = rownames(states),id.method = 'identify',
       simulate = TRUE, main = "Q-Q Plot")

residuals(fit)
rstudent(fit)  # 学生化残差（studentized residual）

## 误差项独立性检测
durbinWatsonTest(fit) 
# p值不显著（p=0.268）说明无自相关性，误差项之间独立

## 线性假设检测
crPlots(fit) # 成分残差图
# 如果不是线性，重新建模，加交互项或变量变换

## 同方差性检测
ncvTest(fit)
# 若检验不显著，则说明误差同方差性存在，适合线性模型。
spreadLevelPlot(fit)
# 图中点在水平的最佳拟合曲线周围呈水平随机分布，则说明误差同方差性

## 预测变量多重共线性检测
# 多重共线性可用统计量VIF（Variance Inflation Factor，方差膨胀因子）进行检测
vif(fit)
sqrt(vif(fit)) > 2  # 为真，则预测变量存在多重共线性问题

## 离群点检测
# 离群点是指那些模型预测效果不佳的观测点。
# 标准化残差值大于2或者小于-2的点可能是离群点
outlierTest(fit)

4. 线性模型假设的综合检验

### 4. 线性模型假设的综合检验
library(gvlma)
gvmodel <- gvlma(fit)
# 数据满足OLS回归模型所有的统计假设(p>0.05)

## 改进措施
# 删除离群点 需谨慎
# 变量变换
# 增删变量
# 尝试其他方法

5.模型比较

### 5.模型比较
fit1 <- lm(weight ~ height,data = women)
fit2 <- lm(weight ~ height + I(height^2),data = women)
anova(fit1,fit2)
# anova()函数可以比较两个嵌套模型的拟合优度

AIC(fit1,fit2)
# AIC（Akaike Information Criterion，赤池信息准则）:
# 考虑了模型的统计拟合度以及用来拟合的参数数目
# AIC值越小的模型要优先选择，它说明模型用较少的参数获得了足够的拟合度。

6.变量选择

### 6.变量选择
##逐步回归法（stepwise method）
library(MASS)
fit <- lm(Murder ~ Population + Illiteracy + Income + Frost,
          data = states)
stepAIC(fit,direcetion = 'backward')
# 选择AIC值小的模型。虽然它可能会找到一个好的模型，
# 但是不能保证模型就是最佳模型，因为不是每一个可能的模型都被评价了。

##全子集回归（all-subsets regression）
library(leaps)
leaps <- regsubsets(Murder ~ Population + Illiteracy + Income + Frost,
           data = states)

plot(leaps,scale = 'adjr2')
# 纵坐标为：调整后的R平方，R平方含义是预测变量解释响应变量的程度
# 大部分情况中，全子集回归要优于逐步回归，因为考虑了更多模型。但是，当有大量预测变
# 量时，全子集回归会很慢。

参考：《R语言实战》，人民邮电出版社