背景
线性回归是一种最基本的预测方法。线性回归根据数值型或类别型输入(称为自变量)建模一个数值型的期望值(称作因变量)。
数据
数据使用的为PUMS数据小样本。这里使用的输入变量包括年龄(AGEP)、性别(SEX)、员工类别(COW)以及受教育程度(SCHL),输出变量是个人收入(PINCP)。
实现过程
加载数据
load("psub.RData")
dtrain <- subset(psub, ORIGRANDGROUP >= 500)
dtest <- subset(psub, ORIGRANDGROUP < 500)
构建线性回归模型
model <- lm(log(PINCP, base = 10) ~ AGEP + SEX + COW + SCHL, data = dtrain)
dtest$predLogPINCP <- predict(model, newdata = dtest)
dtrain$predLogPINCP <- predict(model, newdata = dtrain)
预测
在上一步中,我们已经建模成功并且做了预测,在这个我们通过可视化来观察模型效果。
# 将收入的对数视为被预测收入对数的函数,绘图
ggplot(data = dtest, aes(x = predLogPINCP, y = log(PINCP, base = 10))) +
geom_point(alpha = 0.2, color = "black") +
geom_smooth(aes(x = predLogPINCP, y = log(PINCP, base = 10)), color = "black") +
geom_line(aes(x = log(PINCP, base = 10), y = log(PINCP, base = 10)), color = "blue", linetype = 2) +
scale_x_continuous(limits = c(4, 5)) +
scale_y_continuous(limits = c(3.5, 5.5))
# 将残差收入(预测收入与实际收入之差)视为被预测收入对数的函数
ggplot(data = dtest, aes(x = predLogPINCP, y = predLogPINCP - log(PINCP, base = 10))) +
geom_point(alpha = 0.2, color = "black") +
geom_smooth(aes(x = predLogPINCP, y = predLogPINCP - log(PINCP, base = 10)), color = "black")
上面的两幅图都可以看出来,数据点离散的很厉害,这可能是模型拟合的质量不高。
除了上图,还可以是使用两个指标:R-平方、均方根误差(RMSE)。
# 计算R-平方
rsq <- function(y, f) {
1 - sum((y - f)^2) / sum((y - mean(y))^2)
}
> rsq(log(dtrain$PINCP, base = 10), predict(model, newdata = dtrain))
[1] 0.3382568
> rsq(log(dtest$PINCP, base = 10), predict(model, newdata = dtest))
[1] 0.2605496
R-平方应该相当大(最大可以取1.0),并且在测试集和训练集上的R-平方应该相似。如果测试集的值明显降低,则模型是过拟合的,在上面显示的值,两个都很低(高的为0.7~1.0),说明这个模型是低质量的,但不是过拟合的。
# 计算均方根误差
rmse <- function(y, f) {
sqrt(mean((y - f)^2))
}
> rmse(log(dtrain$PINCP, base = 10), predict(model, newdata = dtrain))
[1] 0.2651856
> rmse(log(dtest$PINCP, base = 10), predict(model, newdata = dtest))
[1] 0.2752171
均方根误差越小越好,一种方法是引入更有用的解释性变量。
> summary(model)
Call:
lm(formula = log(PINCP, base = 10) ~ AGEP + SEX + COW + SCHL,
data = dtrain)
Residuals:
Min 1Q Median 3Q Max
-1.29220 -0.14153 0.02458 0.17632 0.62532
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.973283 0.059343 66.954 < 2e-16 ***
AGEP 0.011717 0.001352 8.666 < 2e-16 ***
SEXF -0.093133 0.023405 -3.979 7.80e-05 ***
COWFederal government employee 0.059727 0.060927 0.980 0.327343
COWLocal government employee -0.034417 0.048030 -0.717 0.473928
COWPrivate not-for-profit employee -0.133524 0.039223 -3.404 0.000709 ***
COWSelf-employed incorporated -0.072427 0.068093 -1.064 0.287928
COWSelf-employed not incorporated -0.060609 0.069244 -0.875 0.381779
COWState government employee -0.081284 0.057796 -1.406 0.160146
SCHLAssociate's degree 0.221432 0.052094 4.251 2.49e-05 ***
SCHLBachelor's degree 0.393834 0.043249 9.106 < 2e-16 ***
SCHLDoctorate degree 0.306554 0.160127 1.914 0.056058 .
SCHLGED or alternative credential 0.137768 0.078192 1.762 0.078612 .
SCHLMaster's degree 0.477284 0.050895 9.378 < 2e-16 ***
SCHLProfessional degree 0.678774 0.087321 7.773 3.52e-14 ***
SCHLRegular high school diploma 0.101746 0.042628 2.387 0.017316 *
SCHLsome college credit, no degree 0.163152 0.042729 3.818 0.000149 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.2691 on 578 degrees of freedom
Multiple R-squared: 0.3383, Adjusted R-squared: 0.3199
F-statistic: 18.47 on 16 and 578 DF, p-value: < 2.2e-16
模型概要主要包括模型调用概要、残差概要、系数和模型质量概要。
线性回归总结
- 线性回归是针对数值量首选的统计建模方法;
- 总是应该首先尝试线性回归,只有更复杂方法确实优于线性回归模型时才使用更复杂的方法;
- 当涉及的数据集中存在大量的变量、或者类别型变量有大量的类别值时,线性回归将遇到麻烦;
- 可以通过加入新变量或转换某些变量来提高线性回归的处理能力(如我们使用的y的log()变换,但变换y时要慎重,因为它会改变其误差模型);
- 使用线性回归,就要考虑残差,要寻找与误差相关的变量并设法加入到模型中,消除系统性的建模误差;
- 即便存在相关的变量,线性回归也能比较好地进行预测,但相关变量的存在会降低所给出建议的质量;
- 过大的系数量级、系数估计方面过大的标准误差以及系数上的错误符号都将预示输入变量具有相关性;
- 线性回归包中有某些可用的最佳内嵌诊断工具,但最有效的安全检查方式仍然是在测试数据上重新检查所构建的模型。
参考文献
[1][数据科学 理论、方法与R语言实践]