笔记&代码 | 统计学——基于R（第四版）第十章多元线性回归

R语言不能阻挡我对spss的爱

大数据界Olu

9594人浏览 · 2021-12-30 17:28:36

大数据界Olu · 2021-12-30 17:28:36 发布

在这里插入图片描述

多元线性回归模型及其参数估计

多元线性回归建模的步骤

确定所关注的因变量𝑦和影响因变量的𝑘个自变量
假定因变量𝑦与𝑘个自变量之间为线性关系，并建立线性关系模型
对模型进行估计和检验
判别模型中是否存在多重共线性，如果存在，进行处理
利用回归方程进行预测
对回归模型进行诊断

回归模型与回归方程

参数的最小二乘估计

25家餐馆的调查数据，建立多元线性回归模型，并解释各回归系数的含义

#回归模型的拟合
> model1<-lm(y~x1+x2+x3+x4+x5,data=example10_1)
> summary(model1)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4 + x5, data = example10_1)
Residuals:
     Min       1Q   Median       3Q      Max 
-16.7204  -6.0600   0.7152   3.2144  21.4805 
Coefficients:
              Estimate Std. Error t value Pr(>|t|)   
(Intercept)  4.2604768 10.4679833   0.407  0.68856   
x1           0.1273254  0.0959790   1.327  0.20037   
x2           0.1605660  0.0556834   2.884  0.00952 **
x3           0.0007636  0.0013556   0.563  0.57982   
x4          -0.3331990  0.3986248  -0.836  0.41362   
x5          -0.5746462  0.3087506  -1.861  0.07826 . 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 10.65 on 19 degrees of freedom
Multiple R-squared:  0.8518,	Adjusted R-squared:  0.8128 
F-statistic: 21.84 on 5 and 19 DF,  p-value: 2.835e-07

多元线性回归方程：
y=4.2604768
+0.1273254x1
+0.1605660x2
+0.0007636x3
-0.3331990x4
-0.5746462*x5

#计算回归系数的置信区间
> confint(model1,level=0.95)
                    2.5 %       97.5 %
(Intercept) -17.649264072 26.170217667
x1           -0.073561002  0.328211809
x2            0.044019355  0.277112598
x3           -0.002073719  0.003600932
x4           -1.167530271  0.501132297
x5           -1.220868586  0.071576251

对于x1来说，其含义是：
x2345不变的条件下，周边居民每变动1万人，日平均营业额变动在 -0.073561002~0.328211809之间

#输出方差分析表
> anova(model1)
Analysis of Variance Table
Response: y
          Df  Sum Sq Mean Sq F value    Pr(>F)    
x1         1 10508.9 10508.9 92.7389 9.625e-09 ***
x2         1  1347.1  1347.1 11.8878  0.002696 ** 
x3         1    85.4    85.4  0.7539  0.396074    
x4         1    40.5    40.5  0.3573  0.557082    
x5         1   392.5   392.5  3.4641  0.078262 .  
Residuals 19  2153.0   113.3                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

拟合优度和显著性检验

模型的拟合优度

多重决定系数R^2=SSR/SST
估计标准误Se

模型的显著性检验

线性关系检验
整体显著性检验
回归系数检验
每个回归系数分别进行t检验

对例题模型的线性关系和回归系数分别进行显著性检验（a=0.05）

线性关系的显著性检验：假设H0=β12345=0
第一个代码块中可得F=21.84，p<0.05，表示关系显著
回归系数的显著性检验：假设H0：任意参数β=0
第一个代码块中可得t12345，只有p2<0.05，所以只有x2影响显著

模型诊断

绘制残差图

plot(model1)

在这里插入图片描述
点2、4、16具有较大残差

在这里插入图片描述

表示残差的正态性假定存在问题，点2、4影响较大，意味着上述多元线性回归模型有可能遗漏了一些重要的回归项

多重共线性及其处理

指回归模型中>=2个自变量彼此相关

多重共线性及其识别

多重共线性所产生的问题

可能会使回归的结果造成混乱，甚至会把分析引入歧途
可能对参数估计值的正负号产生影响，特别是各回归系数的正负号有可能同预期的正负号相反

多重共线性的识别

模型中各对自变量之间显著相关

> library(psych)
> corr.test(example10_1[3:7],use="complete")

Call:corr.test(x = example10_1[3:7], use = "complete")
Correlation matrix 
      x1    x2    x3    x4    x5
x1  1.00  0.74  0.88 -0.62 -0.28
x2  0.74  1.00  0.55 -0.54 -0.32
x3  0.88  0.55  1.00 -0.52 -0.29
x4 -0.62 -0.54 -0.52  1.00  0.10
x5 -0.28 -0.32 -0.29  0.10  1.00
Sample Size 
[1] 25
Probability values (Entries above the diagonal are adjusted for multiple tests.) 
     x1   x2   x3   x4   x5
x1 0.00 0.00 0.00 0.01 0.47
x2 0.00 0.00 0.03 0.03 0.46
x3 0.00 0.00 0.00 0.04 0.47
x4 0.00 0.01 0.01 0.00 0.65
x5 0.18 0.12 0.16 0.65 0.00

 To see confidence intervals of the correlations, print with the short=FALSE option

可得只有x5与其它4个自变量之间关系不显著

当模型的线性关系检验(F检验)显著时，几乎所有回归系数的t检验却不显著
回归系数的正负号与预期的相反
用容忍度与方差扩大因子VIF识别

> library(carData)
> library(car)
> vif(model1)
      x1       x2       x3       x4       x5 
8.233159 2.629940 5.184365 1.702361 1.174053 
> 1/vif(model1)
       x1        x2        x3        x4        x5 
0.1214601 0.3802368 0.1928877 0.5874195 0.8517500

可得容忍度均>0.1，vif均<10，说明共线性不严重

变量选择与逐步回归

向前选择
向后剔除
逐步回归

例题使用逐步回归法建立回归模型

#先进行变量选择
> model2<-step(model1)

Start:  AIC=123.39
y ~ x1 + x2 + x3 + x4 + x5

       Df Sum of Sq    RSS    AIC
- x3    1     35.96 2189.0 121.81
- x4    1     79.17 2232.2 122.30
<none>              2153.0 123.39
- x1    1    199.42 2352.4 123.61
- x5    1    392.54 2545.6 125.58
- x2    1    942.22 3095.2 130.47

Step:  AIC=121.81
y ~ x1 + x2 + x4 + x5

       Df Sum of Sq    RSS    AIC
- x4    1     78.22 2267.2 120.69
<none>              2189.0 121.81
- x5    1    445.69 2634.7 124.44
- x2    1    925.88 3114.9 128.63
- x1    1   1133.27 3322.3 130.24

Step:  AIC=120.69
y ~ x1 + x2 + x5

       Df Sum of Sq    RSS    AIC
<none>              2267.2 120.69
- x5    1    404.28 2671.5 122.79
- x2    1   1050.90 3318.1 128.21
- x1    1   1661.83 3929.0 132.43

选择x125

#拟合逐步回归模型
> model2<-lm(y~x1+x2+x5,data=example10_1)
> summary(model2)

Call:
lm(formula = y ~ x1 + x2 + x5, data = example10_1)

Residuals:
    Min      1Q  Median      3Q     Max 
-14.027  -5.361  -1.560   2.304  23.001 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1.68928    6.25242  -0.270  0.78966    
x1           0.19022    0.04848   3.923  0.00078 ***
x2           0.15763    0.05052   3.120  0.00518 ** 
x5          -0.56979    0.29445  -1.935  0.06656 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 10.39 on 21 degrees of freedom
Multiple R-squared:  0.8439,	Adjusted R-squared:  0.8216 
F-statistic: 37.85 on 3 and 21 DF,  p-value: 1.187e-08

得到估计方程：
y=-1.68928
+0.19022x1
+0.15763x2
-0.56979x5

最后进行诊断，判断模型是否满足各种假定，绘制模型诊断图

plot(model2)

在这里插入图片描述
看出残差具有某种曲线关系，意味着可能需要在模型中加二次项

显示残差的正态性假定存在问题

相对重要性和模型比较

自变量的相对重要性

标准化回归方程

例题计算标准化回归系数，分析各自变量对预测日均营业额的相对重要性

> library(lm.beta)
> model1.beta<-lm.beta(model1)
> summary(model1.beta)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4 + x5, data = example10_1)
Residuals:
     Min       1Q   Median       3Q      Max 
-16.7204  -6.0600   0.7152   3.2144  21.4805 
Coefficients:
              Estimate Standardized Std. Error t value Pr(>|t|)   
(Intercept)  4.2604768    0.0000000 10.4679833   0.407  0.68856   
x1           0.1273254    0.3361822  0.0959790   1.327  0.20037   
x2           0.1605660    0.4130034  0.0556834   2.884  0.00952 **
x3           0.0007636    0.1132753  0.0013556   0.563  0.57982   
x4          -0.3331990   -0.0963203  0.3986248  -0.836  0.41362   
x5          -0.5746462   -0.1781104  0.3087506  -1.861  0.07826 . 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 10.65 on 19 degrees of freedom
Multiple R-squared:  0.8518,	Adjusted R-squared:  0.8128 
F-statistic: 21.84 on 5 and 19 DF,  p-value: 2.835e-07

标准化回归系数是Standardized下的数字，按绝对值比大小，最大的是最重要的变量

模型比较

嵌套模型分为完全模型、简化模型

例题用多元线性回归和逐步回归法建立的两个模型的回归结果进行比较（a=0.05）
5个变量的看作完全模型
3个变量的看作简化模型

#上述两个模型的查看
> model1
Call:
lm(formula = y ~ x1 + x2 + x3 + x4 + x5, data = example10_1)
Coefficients:
(Intercept)           x1           x2           x3           x4  
  4.2604768    0.1273254    0.1605660    0.0007636   -0.3331990  
         x5  
 -0.5746462  
 
> model2
Call:
lm(formula = y ~ x1 + x2 + x5, data = example10_1)
Coefficients:
(Intercept)           x1           x2           x5  
    -1.6893       0.1902       0.1576      -0.5698

利用模型比较法，如果没有差异说明逐步回归模型预测效果不差，根据简约原则应选择逐步回归法

> anova(model1,model2)

Analysis of Variance Table
Model 1: y ~ x1 + x2 + x3 + x4 + x5
Model 2: y ~ x1 + x2 + x5
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1     19 2153.0                           
2     21 2267.2 -2   -114.17 0.5038 0.6121

p>0.05不拒绝原假设，没有证据显示两个模型有显著差异
使用anova（）函数要求模型必须嵌套

还可以用AIC准则

> AIC(model2,model1)
       df      AIC
model2  5 193.6325
model1  7 196.3408

结果显示model2的AIC更小，意味着逐步模型回归比包含全部5个变量好

利用回归方程进行预测

例题用逐步回归模型得到的回归方程求日均营业额95%的置信区间和预测区间

> x<-example10_1[,c(3,4,7)] 
> pre<-predict(model2)
> res<-residuals(model2)
> zre<-rstandard(model2)
> con_int<-predict(model2,x,interval="confidence",level=0.95)
> pre_int<-predict(model2,x,interval="prediction",level=0.95)
> mysummary<-data.frame(营业额=example10_1$y, 点预测值=pre, 残差=res, 标准化残差=zre,置信下限=con_int[,2],置信上限=con_int[,3],预测下限=pre_int[,2],预测上限=pre_int[,3])
> round(mysummary,3)

在这里插入图片描述

> x0<-data.frame(x1=50,x2=100,x5=10)
> predict(model2,newdata=x0)
       1 
17.88685

逐步回归模型的回归方程求x1=50,x2=100,x5=100时日均营业额的点预测值、置信区间、预测区间

> predict(model2,data.frame(x1=50,x2=100,x5=10),interval="confidence",level=0.95)
       fit      lwr      upr
1 17.88685 10.98784 24.78585
> predict(model2,data.frame(x1=50,x2=100,x5=10),interval="prediction",level=0.95)
       fit       lwr      upr
1 17.88685 -4.795935 40.56963

哑变量回归

需要将文字用代码表示的类别自变量

在模型中引入哑变量

含有一个哑变量的回归

沿用例题。假定在分析影响日均营业额的因素中，再考虑“交通方便程度”变量，并设其取值为“方便”和“不方便”。为便于理解，原来的5个自变量我们只保留用餐平均支出一个数值自变量。假定调查得到的数据表

建立两个模型
日均营业额与用餐平均支出的一元回归模型

> model_s<-lm(日均营业额~用餐平均支出,data=example10_7)
> summary(model_s)

Call:
lm(formula = 日均营业额 ~ 用餐平均支出, data = example10_7)
Residuals:
     Min       1Q   Median       3Q      Max 
-19.7604 -10.7832   0.7195   4.3343  28.9301 
Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -5.75023    5.25068  -1.095    0.285    
用餐平均支出  0.32394    0.04482   7.227 2.34e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 13.9 on 23 degrees of freedom
Multiple R-squared:  0.6943,	Adjusted R-squared:  0.681 
F-statistic: 52.23 on 1 and 23 DF,  p-value: 2.343e-07

p<0.05表明模型显著，回归系数0.32394
Multiple R-squared: 0.6943，表示误差69%

日均营业额与用餐平均支出和交通方便程度的二元回归模型

> model_dummy<-lm(日均营业额~用餐平均支出+交通方便程度,data=example10_7)
> summary(model_dummy)

Call:
lm(formula = 日均营业额 ~ 用餐平均支出 + 交通方便程度, data = example10_7)

Residuals:
    Min      1Q  Median      3Q     Max 
-19.443 -11.579  -1.256   8.607  23.456 
Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
(Intercept)      -8.45413    4.69817  -1.799  0.08568 .  
用餐平均支出      0.28641    0.04145   6.909 6.15e-07 ***
交通方便程度方便 14.62088    5.17802   2.824  0.00989 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 12.17 on 22 degrees of freedom
Multiple R-squared:  0.7756,	Adjusted R-squared:  0.7552 
F-statistic: 38.02 on 2 and 22 DF,  p-value: 7.269e-08