线性回归是一种预测建模技术，旨在根据一个或多个输入预测变量预测结果变量的值。目标是在预测变量和响应变量之间建立线性关系（数学公式），因此当预测变量值已知时，我们可以使用它来估计响应的值。

介绍

对于此分析，我们将默认使用R附带的“cars”数据集。'cars'是一个标准的内置数据集，可以方便地以简单易懂的方式演示线性回归。只需在R控制台中输入“cars”即可访问此数据集，它包含50个观察（行）和2个变量（列） - Dist和Speed。让我们在这里看前六个观察结果。

head(cars)

## speed dist
## 1 4 2
## 2 4 10
## 3 7 4
## 4 7 22
## 5 8 16
## 6 9 10

图形分析

本节的目的是建立一个简单的回归模型，我们可以通过建立与速度到距离的线性关系来练习。但在进入建模之前，让我们尝试以图形方式理解这些变量。

通常，对于每个独立变量（预测变量），绘制以下图表以显示后续行为：

  • 散点图：可视化预测变量和响应之间的线性关系
  • 箱形图：在变量中发现任何异常值观测值。在预测器中使用异常值会极大地影响预测，因为它们很容易影响最佳拟合线的方向/斜率。
  • 密度图：查看预测变量的分布。理想地，接近正态分布（钟形曲线），而不是向左或向右倾斜是优选的。

让我们看看如何制作每一个。

散点图

散点图可以帮助可视化依赖（响应）变量和独立（预测变量）变量之间的任何线性关系。理想情况下，如果您有多个预测变量，则会针对响应绘制每个预测变量的散点图，以及最佳线，如下所示。

scatter.smooth(x=cars$speed, y=cars$dist, main="Dist ~ Speed")

散点图与上面的平滑线一起表明'dist'和'speed'变量之间的线性增加关系。这是一件好事，因为线性回归的一个基本假设是响应和预测变量之间的关系是线性和加性的。

BoxPlot - 检查异常值

默认情况下，位于1.5 *'四分位间距（'1.5 * IQR）之外的任何数据点都被视为异常值，其中，IQR计算公式为该变量的第25百分位数和第75百分位数值之间的距离。

par(mfrow=c(1, 2)) # 将画分为两个

boxplot(cars$speed, main="Speed", sub=paste("Outlier rows: ", paste(boxplot.stats(cars$speed)$out, collapse=" "))) #画速度的箱线图

boxplot(cars$dist, main="Distance", sub=paste("Outlier rows: ", paste(boxplot.stats(cars$dist)$out, collapse=" "))) # 画距离的箱线图

密度图 - 检查响应变量是否接近正常

library(e1071)
par(mfrow=c(1, 2))
# 画速度的密度图
plot(density(cars$speed), main="Density Plot: Speed", ylab="Frequency", sub=paste("Skewness:", round(e1071::skewness(cars$speed), 2)))
polygon(density(cars$speed), col="red")
# 距离的密度图
plot(density(cars$dist), main="Density Plot: Distance", ylab="Frequency", sub=paste("Skewness:", round(e1071::skewness(cars$dist), 2)))
polygon(density(cars$dist), col="red")

关联

相关性是一种统计量度，表明两个变量之间的线性相关性水平，即成对出现 - 就像我们在速度和dist中所拥有的那样。相关可以取-1到+1之间的值。

如果我们观察速度增加的每个实例，距离也随之增加，那么它们之间就存在很高的正相关性，因此它们之间的相关性将更接近1.反之关系则相反，其中例如，变量之间的相关性将接近-1。任何接近0的值都表明变量之间的关系较弱。

值得注意的是，相关性是一种相对的衡量标准。在构建回归模型时，由于“低”相关性，您无法拒绝变量。低相关性可能表明响应变量中仍然存在无法解释的方差，在这种情况下，我们应该寻找新的变量。

cor(cars$speed, cars$dist) # 计算速度和距离的相关性

0.8068949＃相当大！

建模

现在我们已经在散点图中以图形方式看到了线性关系，并通过计算相关性，让我们看一下构建线性模型的语法。用于构建线性模型的函​​数是lm（）。lm（）函数接受两个主要参数，即：1公式。2数据，数据通常是data.frame，公式是类公式的对象。但最常见的惯例是直接写出公式代替参数，如下所示。

linearMod <- lm(dist ~ speed, data=cars)

## Call:
## lm(formula = dist ~ speed, data = cars)
## 
## Coefficients:
## (Intercept)        speed  
##     -17.579        3.932

现在我们已经构建了线性模型，我们还以Distance（dist）的数学公式的形式建立了预测器和响应之间的关系，作为速度的函数。

对于上面的输出，你可以注意到'Coefficients'部分有两个部分：截距：-17.579，斜率：3.932

这些也称为β系数。换一种说法，

dist = Intercept + (β∗speed)
=> dist = −17.579 + 3.932∗speed

线性回归诊断

现在建立线性模型，我们有一个公式，如果已知相应的速度，我们可以用它来预测dist值。这足以真正使用这个模型吗？还不行！
在使用回归模型之前，您必须确保它具有统计意义。你如何确保这一点？让我们首先打印linearMod的摘要统计信息。

summary(linearMod)

## Call:
## lm(formula = dist ~ speed, data = cars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -29.069  -9.525  -2.272   9.215  43.201 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -17.5791     6.7584  -2.601   0.0123 *  
## speed         3.9324     0.4155   9.464 1.49e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.38 on 48 degrees of freedom
## Multiple R-squared:  0.6511, Adjusted R-squared:  0.6438 
## F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12

p值：检查统计显着性

这个摘要统计数据告诉我们很多事情。其中之一是模型p值（底部最后一行）和各个预测变量的p值（“系数”下的极右列）。p值非常重要，因为我们可以认为只有当这两个p值小于预定的统计显着性水平（理想情况下为0.05）时，线性模型才具有统计学意义。这可以通过行尾的重要星形进行视觉解释。变量的p值旁边的恒星越多，变量越重要。

在线性回归中，空假设是与变量相关的系数等于零。的备择假设是，系数不等于零（即存在有问题的自变量和因变量之间的关系）。

我们可以解释像这样的t值。较大的t值表示系数不太可能纯粹偶然地等于零。

Pr（> | t |）是当零假设（系数等于零或没有关系）为真时，得到t值高于或高于观测值的概率。因此，如果Pr（> | t |）为低，则系数是显着的（与零显着不同）。如果Pr（> | t |）为高，则系数不显着。

这对我们意味着，当p值小于显着性水平时，我们可以安全地拒绝预测因子的系数β为零的零假设。在我们的例子中，线性模型，这两个p值都远低于0.05阈值，因此我们可以得出结论，我们的模型确实具有统计显着性。

在我们继续使用它来预测（或估计）因变量之前，模型具有统计显着性非常重要，否则，来自该模型的预测值的置信度降低并且可能被解释为偶然事件。 ANCE

如何计算t统计量和p值？

当模型系数和标准误差已知时，计算t统计量和p值的公式如下：

t−Statistic = β−coefficient / Std.Error
modelSummary <- summary(linearMod) # 把模型转为object
modelCoeffs <- modelSummary$coefficients #模型系数
beta.estimate <- modelCoeffs["speed", "Estimate"] # 得到最好的速度预估
std.error <- modelCoeffs["speed", "Std. Error"] #计算速度的标准误差
t_value <- beta.estimate/std.error #计算t检验
p_value <- 2*pt(-abs(t_value), df=nrow(cars)-ncol(cars)) # 计算p值

"t Statistic:  9.46398999029837"
"p Value:  1.4898364962951e-12"
f_statistic <- linearMod$fstatistic[1] # F检验
f <- summary(linearMod)$fstatistic 
model_p <- pf(f[1], f[2], f[3], lower=FALSE)

"Model F Statistic:  "
"Model p-Value:  1.48983649629509e-12"

R-Squared和调整后的R-Squared

数据中的实际信息是它包含的总变化，还记得吗？R-Squared告诉我们的是该模型解释的因变量（响应）变量的变化比例。此外，我们不一定丢弃基于低R-Squared值的模型。在决定模型的功效时，查看AIC和验证样本的预测准确性是一种更好的做法。

现在这是关于R-Squared。调整后的R-Squared怎么样？在向模型添加项时，新模型的R-Squared值将始终大于其子集的R-Squared值。这是因为，由于原始模型中的所有变量也存在，它们对解释因变量的贡献仍然保留在超集中，因此，无论我们添加什么新变量，都只能增强（如果不是很明显）已经解释的内容。以下是调整后的R-Squared值的帮助。Adj R-Squared惩罚模型中术语数（读取预测变量）的总值。因此，在比较嵌套模型时，最好在R平方上查看adj-R平方值。

AIC和BIC

Akaike的信息标准AIC（Akaike，1974）和贝叶斯信息标准BIC（Schwarz，1978）是估计统计模型的拟合优度的度量，并且还可以用于模型选择。两个标准都取决于估计模型的似然函数L的最大值。
AIC定义为：

AIC =（ - 2）•ln（L）+ 2•k
其中k是模型参数的数量，BIC定义为：
BIC =（ - 2）•ln（L）+ k• ln（n）
其中n是样本大小。
对于模型比较，优选具有最低AIC和BIC分数的模型。

AIC(linearMod) # AIC => 419.1569
BIC(linearMod) # BIC => 424.8929

如何评价整体的模型成果？

选择模型时要考虑的最常见指标是：

统计 标准
R平方 越高越好
Adj R-Squared 越高越好
AIC 越低越好
BIC 越低越好
Mallows cp 应该接近模型中预测变量的数量
MAPE（平均绝对百分比误差） 越低越好
MSE（均方误差） 越低越好
Min_Max准确度 - >平均值（min（实际，预测）/ max（实际，预测）） 越高越好

使用线性模型进行预测

到目前为止，我们已经了解了如何使用整个数据集构建线性回归模型。如果我们以这种方式构建它，就无法告诉模型将如何预测新数据。因此，首选做法是将数据集拆分为80:20，然后在80％样本上构建模型，20%样本通过模型预测结果。通过计算精度测量（如min_max精度）和误差率（MAPE或MSE），我们可以找出模型的预测精度。现在来开始实操

步骤1：从原始数据创建训练（开发）和测试（验证）数据样本。

set.seed(100) # 设定随机数种子
trainingRowIndex <- sample(1:nrow(cars), 0.8*nrow(cars)) #80%训练集的行
trainingData <- cars[trainingRowIndex, ] # 将训练集提出来
testData <- cars[-trainingRowIndex, ] # 这个是测试集

步骤2：在训练数据上构建模型并使用其来预测

lmMod <- lm(dist ~ speed, data=trainingData) #构建模型
distPred <- predict(lmMod, testData) # 预测结果

第3步：诊断模型

summary (lmMod)

## Call:
## lm(formula = dist ~ speed, data = trainingData)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -23.350 -10.771  -2.137   9.255  42.231 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -22.657      7.999  -2.833  0.00735 ** 
## speed          4.316      0.487   8.863 8.73e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 15.84 on 38 degrees of freedom
## Multiple R-squared:  0.674,  Adjusted R-squared:  0.6654 
## F-statistic: 78.56 on 1 and 38 DF,  p-value: 8.734e-11


# 计算AIC值
AIC (lmMod) 338.4489

从模型总结来看，模型p值和预测变量的p值小于显着性水平，因此我们知道我们有一个统计上显着的模型。

第4步：计算预测准确度和错误率

实际值和预测值之间的简单相关可以用作精度测量的一种形式。较高的相关性意味着实际值和预测值由相似的方向移动，即当实际值增加时，预测值也增加，反之亦然。 现在让我们计算最小最大准确度和MAPE值

actuals_preds <- data.frame(cbind(actuals=testData$dist, predicteds=distPred)) # 构建预测数据框
correlation_accuracy <- cor(actuals_preds) # 82.7%



MinMaxAccuracy = Average(Min(Actuals,Predicteds) / Max(Actuals,Predicteds))

MeanAbsolutePercentageError(MAPE) = Average(abs(Predicteds−Actuals) / Actuals))
min_max_accuracy <- mean (apply(actuals_preds, 1, min) / apply(actuals_preds, 1, max)) # => 58.42%, 最小最大准确度
mape <- mean(abs((actuals_preds$predicteds - actuals_preds$actuals))/actuals_preds$actuals) # => 48.38%,MAPE值

 

对样本进行交叉验证（k-fold交叉验证）

假设模型在20％的测试数据上预测令人满意，那我们是否可以认为i模型始终表现出色？最佳做法是尽可能严格地测试模型的性能。一种方法是将模型“建立”在不同的训练数据子集上并在剩余数据上进行预测时。具体的做法就是将数据随机拆分为'k'份互斥的样本部分。每次我们将其中一个部分作为测试数据，在剩余的（k-1部分）数据上建立模型并计算预测的均方误差。

通过这样做，我们可以检查模型的预测精度是否对于任何一个特定样本都没有太大的变化，并且最佳拟合线是否与斜率和水平相比变化不大，它们应该是平行的并且尽可能彼此接近。

library(DAAG)
par(new=TRUE); par(bg="aliceblue"); par(mar=c(4,4,4,4)) # set margins
# 检查虚线是否平行，并且对于任何一种特定颜色，小符号和大符号不会过分散。
cvResults <- suppressWarnings(CVlm(df=cars, form.lm=dist ~ speed, m=5, dots=FALSE, seed=29, legend.pos="topleft", printit=FALSE, main="Small symbols are predicted values while bigger ones are actuals.")); # 进行交叉验证，5折

# 计算均方误差
attr(cvResults, 'ms') 251.2783