乐趣区

关于算法:R语言用线性回归模型预测空气质量臭氧数据

原文链接:http://tecdat.cn/?p=11387

只管线性模型是最简略的机器学习技术之一,但它们依然是进行预测的弱小工具。这尤其是因为线性模型特地容易解释这一事实。在这里,我将探讨应用空气质量数据集的一般最小二乘回归示例解释线性模型时最重要的方面。

空气质量数据集

空气质量数据集蕴含以下四个空气质量指标的 154 次测量:

  • 臭氧:均匀臭氧程度,以十亿分之一为单位
  • Solar.R:太阳辐射 
  • 风:均匀风速,每小时英里
  • 温度:每日最高温度,以华氏度为单位

咱们将通过删除所有NA 并排除  Month 和Day 列来清理数据集,抉择局部预测变量。

data(airquality)
ozone <- subset(na.omit(airquality), 
        select = c("Ozone", "Solar.R", "Wind", "Temp"))

数据摸索和筹备

预测工作如下:依据太阳辐射,风速和温度,咱们能够预测臭氧程度吗?要查看线性模型的假如是否适宜手头的数据,咱们将计算变量之间的相关性:

# 散点图矩阵
plot(ozone)

 

# 成对变量的相关性
cors <- cor(ozone)
print(cors)
##              Ozone    Solar.R       Wind       Temp
## Ozone    1.0000000  0.3483417 -0.6124966  0.6985414
## Solar.R  0.3483417  1.0000000 -0.1271835  0.2940876
## Wind    -0.6124966 -0.1271835  1.0000000 -0.4971897
## Temp     0.6985414  0.2940876 -0.4971897  1.0000000
# 哪些变量高度相干,排除自相干

print(cor.names)
## [1] "Wind+Ozone: -0.61" "Temp+Ozone: 0.7"   "Ozone+Wind: -0.61"
## [4] "Ozone+Temp: 0.7"

因为臭氧参加两个线性相互作用,即:

  • 臭氧与温度呈正相干
  • 臭氧与风速负相关

这表明应该有可能应用其余特色来造成预测臭氧程度的线性模型。

分为训练和测试集

咱们将抽取 70%的样本进行训练,并抽取 30%的样本进行测试:

set.seed(123)
N.train <- ceiling(0.7 * nrow(ozone))
N.test <- nrow(ozone) - N.train
trainset <- sample(seq_len(nrow(ozone)), N.train)
testset <- setdiff(seq_len(nrow(ozone)), trainset)

钻研线性模型

为了阐明解释线性模型的最重要方面,咱们将通过以下形式训练训练数据的一般最小二乘模型:

为了解释模型,咱们应用以下 summary 函数:

model.summary <- summary(model)
print(model.summary)
## 
## Call:

## Residuals:
##     Min      1Q  Median      3Q     Max 
## -36.135 -12.670  -2.221   9.420  65.914 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -65.76604   22.52381  -2.920 0.004638 ** 
## Solar.R       0.05309    0.02305   2.303 0.024099 *  
## Temp          1.56320    0.25530   6.123 4.03e-08 ***
## Wind         -2.61904    0.68921  -3.800 0.000295 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 18.17 on 74 degrees of freedom
## Multiple R-squared:  0.5924, Adjusted R-squared:  0.5759 
## F-statistic: 35.85 on 3 and 74 DF,  p-value: 2.039e-14

残差

咱们取得的第一条信息是残差。

残差值表明,该模型通常预测的臭氧值略高于观测值。然而,最大值很大,表明某些离群值预测太低了。查看数字可能有点形象,因而让咱们依据察看值绘制模型的预测:

res <- residuals(model)

# 向图中增加残差
segments(obs, pred, obs, pred + res)

系数

当初咱们理解了残差,让咱们看一下系数。咱们能够应用该coefficients 函数来获取模型的拟合系数:

##  (Intercept)      Solar.R         Temp         Wind 
## -65.76603538   0.05308965   1.56320267  -2.61904128

请留神,模型的截距值非常低。这是在所有独立值均为零的状况下模型预测的值。低系数  Solar.R 示意太阳辐射对预测臭氧程度没有重要作用,这难能可贵,因为在咱们的探索性剖析中,它与臭氧程度没有很大的相关性。系数 Temp示意温度高时臭氧程度高(因为臭氧会更快造成)。系数  Wind 通知咱们风速快时臭氧程度会升高(因为臭氧会被吹走)。

与系数关联的其余值提供无关预计的统计确定性的信息。

##                 Estimate  Std. Error   t value     Pr(>|t|)
## (Intercept) -65.76603538 22.52380940 -2.919845 4.638426e-03
## Solar.R       0.05308965  0.02305379  2.302860 2.409936e-02
## Temp          1.56320267  0.25530453  6.122894 4.034064e-08
## Wind         -2.61904128  0.68920661 -3.800081 2.946349e-04
  • Std. Error 是系数预计的标准误差
  • t value 以标准误差示意系数的值
  • Pr(>|t|) 是 t 测验的 p 值,示意测验统计量的重要性

标准误差

系数的标准误差定义为特色方差的标准偏差:

在 R 中,能够通过以下形式计算模型预计的标准误差:

##              (Intercept)       Solar.R         Temp          Wind
## (Intercept) 507.32198977  2.893612e-02 -5.345957524 -9.940961e+00
## Solar.R       0.02893612  5.314773e-04 -0.001667748  7.495211e-05
## Temp         -5.34595752 -1.667748e-03  0.065180401  6.715467e-02
## Wind         -9.94096142  7.495211e-05  0.067154670  4.750058e-01
## (Intercept)     Solar.R        Temp        Wind 
## 22.52380940  0.02305379  0.25530453  0.68920661

当初,你可能想晓得这些值的 vcov。它定义为矩阵的方差 - 协方差矩阵,该矩阵按误差的方差标准化:

##              (Intercept)       Solar.R         Temp          Wind
## (Intercept) 507.32198977  2.893612e-02 -5.345957524 -9.940961e+00
## Solar.R       0.02893612  5.314773e-04 -0.001667748  7.495211e-05
## Temp         -5.34595752 -1.667748e-03  0.065180401  6.715467e-02
## Wind         -9.94096142  7.495211e-05  0.067154670  4.750058e-01

用于标准化的方差 - 协方差矩阵的方差是误差的预计方差,其定义为

 cov.unscaled 参数是方差 - 协方差矩阵:

# 通过 'model.matrix' 将截距作为特色
X <- model.matrix(model) # design matrix
# 求解 X ^ T%*%X = I 求 X ^ T * X 的逆
unscaled.var.matrix <- solve(crossprod(X), diag(4))
print(paste("Is this the same?", isTRUE(all.equal(unscaled.var.matrix, model.summary$cov.unscaled, check.attributes = FALSE))))
## [1] "Is this the same? TRUE"

t 值

t 值定义为

在 R 中 

## (Intercept)     Solar.R        Temp        Wind 
##   -2.919845    2.302860    6.122894   -3.800081

p 值

在所有系数 βi=0 的假如下计算 p 值。t 值遵循 t 散布

model.df <- df.residual(model)

自由度。线性模型的自由度定义为

其中 n 是样本数,p 是特色数(包含 inctercept)。p 值示意取得的系数预计纯正是偶尔地与零不同的可能性。因而,低 p 值表明变量与后果之间存在显着关联。

进一步统计

summary 函数提供以下附加统计信息:R 方,调整后的 R 方和 F 统计。

残余标准偏差

顾名思义,_残余标准偏差_是模型的均匀 RSS(MSE)的平方根:

## [1] 18.16979

_残余标准偏差_仅示意模型的均匀精度。在这种状况下,该值非常低,表明该模型具备良好的拟合度。

R 方

R 方示意确定系数。它定义为估计值与察看到的后果之间的相关性的平方:

## [1] 0.5924073

与 [-1,1] 中的相关性相同,R 平方在[0,1] 中。

调整后的 R 方

调整后的 R 方值会依据模型的复杂性来调整 R 方:

其中 n 是察看数,p 是特色数。因而,调整后的 R 方能够像这样计算:

n <- length(trainset) # 样本数

print(r.squared.adj)
## [1] 0.5758832

如果 R 平方和调整后的 R 方之间存在相当大的差别,则表明能够思考缩小特色空间。

F 统计

F 统计量定义为已解释方差与无法解释方差的比率。为了进行回归,F 统计量始终批示两个模型之间的差别,其中模型 1(p1)由模型 2(p2)的特色子集定义:

F 统计量形容模型 2 的预测性能(就 RSS 而言)优于模型 1 的水平。报告的默认 F 统计量是指训练后的模型与仅截距模型之间的差别:

## 
## Call:
## 
## Coefficients:
## (Intercept)  
##       36.76

 

因而,测试的零假如是惟一的截距 - 模型的拟合和指定的模型是相等的。如果能够回绝原假如,则意味着指定模型比原模型具备更好的拟合度。

让咱们通过手工计算得出这个想法:

rss <- function(model) {return(sum(model$residuals^2))
}

# 比拟仅截距模型和具备 3 个模型
f.statistic(null.model, model)
## [1] 35.85126

在这种状况下,F 统计量具备较大的值,这表明咱们训练的模型显著优于仅截距模型。

置信区间

置信区间是解释线性模型的有用工具。默认状况下,confint 计算 95%置信区间(±1.96σ^±1.96σ^):

ci <- confint(model) 
##                (Intercept)                    Solar.R 
## "95% CI: [-110.65,-20.89]"       "95% CI: [0.01,0.1]" 
##                       Temp                       Wind 
##      "95% CI: [1.05,2.07]"    "95% CI: [-3.99,-1.25]"

这些值表明模型对截距的预计不确定。这可能表明须要更多数据能力取得更好的拟合度。

检索估计值的置信度和预测区间

通过提供自interval 变量,能够将线性模型的预测转换为区间。这些区间给出了对预测值的置信度。区间有两种类型:置信区间和预测区间。让咱们将模型利用于测试集,应用不同的参数作为  interval 参数,以查看两种区间类型之间的差别:

# 计算预测的置信区间(CI)preds.ci <- predict(model, newdata = ozone[testset,], interval = "confidence") 
##      fit                 lwr                 upr                Method
## [1,] "-4.42397219873667" "-13.4448767773931" "4.59693237991976" "CI"  
## [2,] "-22.0446990408131" "-61.0004555440645" "16.9110574624383" "PI"

置信区间是窄区间,而预测 区间是宽区间。它们的值基于level 参数指定的提供的重要性程度(默认值:0.95)。

它们的定义略有不同。给定新的观测值 x,CI 和 PI 定义如下

其中 tα/ 2,dftα/ 2,df 是 df = 2 自由度且显着性程度为 α 的 t 值,σerr 是残差的标准误差,σ2xσx2 是独立特色的方差,x(x)示意特色的平均值。


最受欢迎的见解

1.R 语言多元 Logistic 逻辑回归 利用案例

2. 面板平滑转移回归 (PSTR) 剖析案例实现

3.matlab 中的偏最小二乘回归(PLSR)和主成分回归(PCR)

4.R 语言泊松 Poisson 回归模型剖析案例

5.R 语言回归中的 Hosmer-Lemeshow 拟合优度测验

6.r 语言中对 LASSO 回归,Ridge 岭回归和 Elastic Net 模型实现

7. 在 R 语言中实现 Logistic 逻辑回归

8.python 用线性回归预测股票价格

9.R 语言如何在生存剖析与 Cox 回归中计算 IDI,NRI 指标

退出移动版