R语言 CHAR 10

多元线性回归

  1. 确定因变量y和影响因变量的k个自变量
  2. 假定为线性关系,建立线性关系模型
  3. 对模型进行估计与检验
  4. 判别是否有多重共线性,如果存在,进行处理
  5. 利用回归方程进行预测
  6. 对回归模型进行诊断

回归模型拟合

mode<-lm(y~x1+x2+x3+x4+x5,data=table)
summary(mode)

计算回归系数的置信区间

confint(mode,level=0.95)

输出方差分析表

anova(mode)

绘制残差图诊断模型

plot(mode)

对模型中各对自变量之间的相关系数进行检验

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

容忍度与膨胀因子

library(car)
vif(mode) #膨胀因子
1/vif(mode) #容忍度

逐步回归

mode1<-step(mode)

标准化回归系数

library(lm.beta)
mode.beta<-lm.beta(mode)
summary(mode.beta)

anova(要求模型相互嵌套)

model1<-lm(y~x1+x2+x3+x4+x5,data=table)
model2<-lm(y~x1+x2+x5,data=table)
anova(model1,model2)

AIC(不要求嵌套)

model1<-lm(y~x1+x2+x3+x4+x5,data=table)
model2<-lm(y~x1+x2+x5,data=table)
AIC(model1,model2)


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

10.1.1 回归模型与回归方程

  • 多元线性回归模型

请添加图片描述

e为误差项
对e的三个基本假定

  • 正态性
  • 方差齐性
  • 独立性
  • 多元线性回归方程
    请添加图片描述
  • 估计的多元线性回归方程
    请添加图片描述

10.1.2 参数的最小二乘估计

请添加图片描述
请添加图片描述

请添加图片描述
请添加图片描述

  • 回归模型拟合
> table<-read.csv("/Users/zhourui/Documents/example10_1.csv")
> mode<-lm(y~x1+x2+x3+x4+x5,data=table)
> summary(mode)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4 + x5, data = table)

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 = 0.1273254 * x1 + 0.1605660 * x2+0.0007636 * x3 -0.3331990 * x4 -0.5746462 * x5+4.2604768

x2对y值影响较大

  • 计算回归系数的置信区间
> confint(mode,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
  • 输出方差分析表
> anova(mode)
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

x1、x2对y值影响较大

10.2 拟合优度和显著性检验

10.2.1 模型的拟合优度

多重决定系数评估

1. 多重决定系数

R方

请添加图片描述

Multiple R-squared:  0.8518
调整后的R方

考虑样本量和模型中自变量的个数
请添加图片描述

Adjusted R-squared:  0.8128 

2. 估计标准误

请添加图片描述

Residual standard error: 10.65 on 19 degrees of freedom

10.2.2 模型的显著性检验

  • 在一元线性回归中

只有一个自变量,F检验(线性关系检验)t检验(回归系数检验)是等价的

  • 多元线性回归

F检验:检验因变量同多个自变量的整体线性关系是否显著(有一个就显著)
t检验:对每个回归系数分别检验,以判断每个自变量对因变量的影响是否显著

1.线性关系检验

请添加图片描述

2.回归系数检验

请添加图片描述
请添加图片描述
请添加图片描述

3.模型诊断

  • 绘制残差图诊断模型
> par(mfrow=c(1,2),mai=c(0.8,0.8,0.4,0.1),cex=0.8,cex.main=0.7)
> plot(mode,which=1:2)

请添加图片描述

10.3 多重共线性及其处理

10.3.1 多重共线性及其识别

  • 多重共线性:当回归模型中两个或两个以上的自变量彼此相关

1. 多重共线性所产生的问题

  • 变量高度相关时,可能会造成回归结果的混乱
  • 对参数估计值的正负号产生影响

2. 多重共线性的识别和处理

  • 对模型中各对自变量之间的相关系数进行检验
> library(psych)
> corr.test(table[3:7],use = "complete")
Call:corr.test(x = table[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对x1 0.47
x5对x2 0.46
x5对x3 0.47
x5对x4 0.65
均大于0.05(P值)
x5与其他自变量间关系不显著

use = “complete”(一对多,形成矩阵)
use = “pair”(一对一)

  • 考察各回归系数的显著性

当模型的F检验显著
多个回归系数的t检验却不显著

(同10_1)

  • 分析回归系数的正负号

回归系数正负号与预期相反

  • 用容忍度和方差膨胀因子来识别

容忍度小于0.1,存在严重共线性
膨胀因子(VIF)越大,共线性越严重,VIF大于10,存在严重共线性

> library(car)
> vif(mode) #膨胀因子
      x1       x2       x3       x4       x5 
8.233159 2.629940 5.184365 1.702361 1.174053 
> 1/vif(mode) #容忍度
       x1        x2        x3        x4        x5 
0.1214601 0.3802368 0.1928877 0.5874195 0.8517500 

回归模型的共线性并不严重

10.3.2 变量选择与逐步回归

在建立模型前有选择的确定进入模型的自变量,可以避免多重共线性的问题

1. 向前选择

  • 引入F最大(P最小)
  • 加入的变量一定会被保留
    请添加图片描述

2.向后剔除

  • 剔除F最小(P最大)
  • 剔除的变量不会再进入
    请添加图片描述
    请添加图片描述

3. 逐步回归

    • 上述两种方法相结合
    • 变量有可能被踢出,也有可能重新加入
      请添加图片描述
  • AIC

请添加图片描述

  • 逐步回归变量选择
> mode1<-step(mode)
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

在最后一步中去掉三个变量中的任何一个都会增大AIC,所以运行停止

重新拟合后进行分析:

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

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

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
  • 逐步回归的方差分析表
> anova(mode12)
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value    Pr(>F)    
x1         1 10508.9 10508.9 97.3392 2.452e-09 ***
x2         1  1347.1  1347.1 12.4775  0.001976 ** 
x5         1   404.3   404.3  3.7447  0.066558 .  
Residuals 21  2267.2   108.0                      
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

y = -1.68928 + 0.19022 * x1 + 0.15763 * x2 -0.56979 * x5

  • 模型诊断
> par(mfrow=c(1,2),mai=c(0.8,0.8,0.4,0.1),cex=0.8,cex.main=0.7)
> plot(mode12,which=1:2)

请添加图片描述

10.4 相对重要性和模型比较

哪些自变量对预测更重要

10.4.1 自变量的相对重要性

评估自变量相对重要性:比较标准回归系数
请添加图片描述
请添加图片描述

标准化回归系数越大,该自变量对因变量预测越重要

  • 标准化回归系数
> library(lm.beta)
> mode.beta<-lm.beta(mode)
> summary(mode.beta)

Call:
lm(formula = y ~ x1 + x2 + x3 + x4 + x5, data = table)

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

请添加图片描述

10.4.2 模型比较

简约模型:最佳模型
逐步回归方法虽然剔除不必要自变量,但是无法确认是否更好拟合,需要进一步比较
目的:找简约模型

  • 嵌套模型

请添加图片描述

  • 检验
    请添加图片描述
    请添加图片描述
    请添加图片描述
    请添加图片描述
  • anova(要求模型相互嵌套)
> model1<-lm(y~x1+x2+x3+x4+x5,data=table)
> model2<-lm(y~x1+x2+x5,data=table)
> 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

RSS 残差平方和
Sum of Sq 交互效应残差和

P=0.6121,不拒绝H0,没有证据显示两个模型有显著差异
由于简约原则选择model2

  • AIC(不要求嵌套)

> model1<-lm(y~x1+x2+x3+x4+x5,data=table)
> model2<-lm(y~x1+x2+x5,data=table)
> AIC(model1,model2)
       df      AIC
model1  7 196.3408
model2  5 193.6325

AIC 小说明模型用较少参数

逐步回归模型更好

10.5 利用回归方程进行预测

请添加图片描述
请添加图片描述

> pre<-predict(model2)
> res<-residuals(model2)
> zre<-rstandard(model2)
> x<-table[,c(3,4,7)]
> con_int<-predict(model2,x,interval = "confidence",level=0.95)
> pre_int<-predict(model2,x,interval = "prediction",level=0.95)
> musummary<-data.frame(营业额=table$y,点预测值=pre,残差=res,标准化残差=zre,置信下限=con_int[,2],置信上限=con_int[,3],预测下限=pre_int[,2],预测上限=pre_int[,3])
> round(musummary,3)
> round(musummary,3)
   营业额 点预测值    残差 标准化残差 置信下限 置信上限 预测下限 预测上限
1    53.2   52.189   1.011      0.102   45.457   58.921   29.557   74.822
2    18.5   -4.501  23.001      2.359  -11.969    2.967  -27.363   18.361
3    11.3   21.963 -10.663     -1.072   15.685   28.240   -0.539   44.464
4    84.7   65.114  19.586      2.766   49.300   80.928   38.337   91.891
5     7.3    6.128   1.172      0.122   -2.283   14.540  -17.059   29.316
6    17.9   22.408  -4.508     -0.466   14.581   30.235   -0.574   45.390
7     2.5    1.278   1.222      0.125   -6.113    8.670  -21.559   24.116
8    27.3   34.719  -7.419     -0.747   28.400   41.038   12.206   57.232
9     5.9   10.628  -4.728     -0.472    4.904   16.351  -11.726   32.981
10   23.9   37.843 -13.943     -1.390   32.193   43.493   15.509   60.178
11   69.4   62.852   6.548      0.709   52.939   72.766   39.079   86.626
12   20.6   18.296   2.304      0.229   13.036   23.556   -3.943   40.535
13    1.9   -5.510   7.410      0.771  -13.694    2.674  -28.616   17.596
14    3.0   14.956 -11.956     -1.202    8.687   21.224   -7.543   37.455
15    7.3    8.860  -1.560     -0.169   -1.050   18.770  -14.913   32.632
16   46.2   29.831  16.369      1.699   21.745   37.917    6.759   52.903
17   78.8   78.016   0.784      0.112   62.105   93.927   51.182  104.850
18   11.1   13.167  -2.067     -0.209    6.420   19.915   -9.470   35.805
19    8.6   15.847  -7.247     -0.815    4.670   27.024   -8.481   40.175
20   48.9   50.727  -1.827     -0.184   44.205   57.250   28.156   73.299
21   22.1   27.461  -5.361     -0.536   21.671   33.251    5.090   49.831
22   11.1    0.233  10.867      1.354  -13.486   13.953  -25.362   25.829
23    8.6   22.627 -14.027     -1.395   17.181   28.073    0.344   44.911
24   48.9   48.961  -0.061     -0.006   42.722   55.200   26.470   71.452
25   22.1   27.005  -4.905     -0.490   21.220   32.790    4.636   49.374
  • 求x为具体值时,日均营业额的点预测值、置信区间和预测区间(新值预测)
> x0<-data.frame(x1=50,x2=100,x5=10)
> predict(model2,newdata = x0)
       1 
17.88685 
> predict(model2,x0,interval = "confidence",level=0.95)
       fit      lwr      upr
1 17.88685 10.98784 24.78585
> predict(model2,x0,interval = "prediction",level=0.95)
       fit       lwr      upr
1 17.88685 -4.795935 40.56963

10.6 哑变量回归

哑变量 == 分类类型变量

10.6.1 在模型中引入哑变量

请添加图片描述

10.6.2 含有一个哑变量的回归

请添加图片描述
请添加图片描述

> table<-read.csv("/Users/zhourui/Documents/example10_72.csv")
> model<-lm(y~X1,data=table)
> summary(model)

Call:
lm(formula = y ~ X1, data = table)

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    
X1           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
  • 二元回归
    请添加图片描述
> model<-lm(y~X1+x2,data=table)
> summary(model)

Call:
lm(formula = y ~ X1 + x2, data = table)

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 .  
X1           0.28641    0.04145   6.909 6.15e-07 ***
x2方便      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
  • 方差分析表
> anova(model)
Analysis of Variance Table

Response: y
          Df  Sum Sq Mean Sq F value    Pr(>F)    
X1         1 10085.9 10085.9  68.062 3.529e-08 ***
x2         1  1181.5  1181.5   7.973  0.009889 ** 
Residuals 22  3260.1   148.2                      
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

请添加图片描述
请添加图片描述
请添加图片描述

  • 分析
    请添加图片描述

请添加图片描述
请添加图片描述
请添加图片描述

  • 0
    点赞
  • 6
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值