R语言回归分析相关流程

数据预处理

中心化标准化

data <- as.data.frame(scale(data, center = T, scale = T))    # 标准化
data <- as.data.frame(scale(data, center = T, scale = F))    # 中心化

变量相关性检验

使用Hmisc

library(Hmisc)
rcorr(as.matrix(data), type = 'pearson')
R <- rcorr(as.matrix(data), type = 'pearson')$r  # 相关系数阵
P <- rcorr(as.matrix(data), type = 'pearson')$P  # 相关性检验的P值

使用psych

library(psych)
corr.test(data,use="complete")

热力图

library(Hmisc)
R <- rcorr(as.matrix(data), type = 'pearson')$r # 相关阵
library(corrplot)
corrplot(R, method="color", addCoef.col='grey', tl.col='black') # 热力图

使用ggplot绘制热图

library(Hmisc)
R <- rcorr(as.matrix(data), type = 'pearson')$r # 相关阵
library(reshape2)  
melted_r <- melt(R)  # 将相关系数矩阵拉直
library(ggplot2)
ggplot(data = melted_r, aes(x=Var1,y=Var2,fill = value)) + # ggplot2 基本图层
  geom_raster() + # 绘制热力图
  scale_fill_gradient2(low = rgb(10,8,32,max = 255),mid = rgb(165,24,90,max = 255),
                       high = rgb(249,228,209,max = 255),limit = c(0,1),
                       name="Pearson\nCorrelation") + 
  theme( # 调整图表样式
    axis.title = element_blank(), # 删除ggplot2的坐标标签
    panel.background = element_blank(), # 删除ggplot2的灰色背景
    panel.grid.major = element_blank(), # 删除ggplot2的面板网格
    panel.border = element_blank(),
    axis.ticks = element_blank() # 删除ggplot2的轴刻度
  ) + 
  geom_text(aes(Var2,Var1,label = round(value,2)),
            color = '#228B22',size = 2) # 加上相关系数数值

偏相关系数

在多要素所构成的系统中，当研究某一个要素对另一个要素的影响或相关程度时，把其他要素的影响视作常数（保持不变），即暂时不考虑其他要素影响，单独研究两个要素之间的相互关系的密切程度，所得数值结果为偏相关系数

计算偏相关系数并绘图

library(corpcor)
cor2pcor(R) # 计算得到偏相关阵
library(corrplot)
corrplot(cor2pcor(R), method="color", addCoef.col='grey', tl.col='black') # 热力图

偏相关系数图

建立模型

最小二乘回归

model <- lm(y~x1+x2+x3+x4+x5, data)     # 最小二乘回归
#model <- lm(y~0+x1+x2+x3+x4+x5, data)   # 无截距项最小二乘回归
fitted(model)                           # 拟合值
model$coefficients                      # 回归系数
model$residuals                         # 残差

计算杠杆值

hii <- hatvalues(mo)                                # 杠杆值

查看杠杆值(按照X大小顺序排序)

计算模型残差

SSE <- sum(mo$residuals^2)                          # 残差平方和
sqrt(SSE/(length(data$x)-2))

回归模型检验

显著性检验

summary(model)$r.squared         # R-square
summary(model)$adj.r.squared     # Adjusted R-square
sum((result - mean(result))^2)   # SSR
sum((data$y - result)^2)         # SSE
sum((data$y - mean(data$y))^2)   # SST, and SST = SSR + SSE
plot(fitted(model), model$residuals, pch = 19, xlab = "fitted values",ylab = "residuals") # 残差图

模型分析

点预测

mo <- lm(y~x,data)                                                  # 最小二乘回归
none_interval <- predict.lm(mo, data, interval="none")              # 计算点预测值

区间预测

confidence_interval <- predict.lm(mo, data, interval="confidence")  # 计算置信区间预测值
prediction_interval <- predict.lm(mo, data, interval="prediction")  # 计算预测区间预测值

绘制区间图

# 对齐数据顺序
data_order <- data[order(data$x),]
confidence_interval_order <- confidence_interval[order(data$x),]
prediction_interval_order <- prediction_interval[order(data$x),]

开始绘图

plot(data$x, data$y, xlab="x", ylab="y", pch = 16)  # 绘制轮廓
polygon(c(data_order$x,rev(data_order$x)),c(prediction_interval_order[,2],rev(prediction_interval_order[,3])),col = rgb(221,234,243,max = 255),border = NA)
polygon(c(data_order$x,rev(data_order$x)),c(confidence_interval_order[,2],rev(confidence_interval_order[,3])),col = "gray",border = NA)  # 其中polygon为绘制多边形的函数，rev为排序函数
abline(mo, col="red", lwd=2)  
points(data$x, data$y, pch = 16)

置信区间和预测区间

预测相关

新值预测

new <- data.frame(x = 3.5)  # 生成新的数据
predict.lm(mo, new, interval="prediction")  # 点预测+预测上下界
predict.lm(mo, new, interval="confidence")  # 点预测+置信上下界

回归分析