多元相关与回归分析及R使用

Author:龙箬
Data Science and Big Data Technology
Change the world with data!
CSDN@weixin_43975035
很多人喜欢把心事扔进河里,就变成了石头

多元相关与回归分析及R使用

1.变量间的关系分析

简单相关分析的R计算
> x1=c(171,175,159,155,152,158,154,164,168,166,159,164) #身高
> x2=c(57,64,41,38,35,44,41,51,57,49,47,46) #体重
> plot(x1,x2) #作散点图

在这里插入图片描述

> lxy<-function(x,y){n=length(x);sum(x*y)-sum(x)*sum(y)/n} #离均差乘积和函数
> lxy(x1,x1) #x1的离均差平方和
[1] 556.9167
> lxy(x2,x2) #x2的离均差平方和
[1] 813
> lxy(x1,x2) #x1与x2得离均差乘积和
[1] 645.5
> r=lxy(x1,x2)/sqrt(lxy(x1,x1)*lxy(x2,x2)) #显示用离均差乘积和计算得相关系数
> r
[1] 0.9593031

也可以使用相关系数函数cor() 计算相关系数

> cor(x1,x2) #计算相关系数
[1] 0.9593031

计算相关系数r的t值:

> tr=r/sqrt((1-r^2)/(n-2)) #相关系数假设检验t统计量
> tr
[1] 10.74298

也可以使用相关系数检验函数cor.test() 计算相关系数r的t值

> cor.test(x1,x2) #相关系数假设检验

	Pearson's product-moment correlation

data:  x1 and x2
t = 10.743, df = 10, p-value = 8.21e-07
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.8574875 0.9888163
sample estimates:
      cor 
0.9593031
简单回归分析的R计算

以上述数据介绍建立直线回归方程的步骤

> x=x1
> y=x2
> b=lxy(x,y)/lxy(x,x) #线性回归方程斜率
> a=mean(y)-b*mean(x) #线性回归方程截距
> c(a=a,b=b) #显示线性回归方程估计值
         a          b 
-140.36436    1.15906
> plot(x,y) #作散点图
> lines(x,a+b*x) #添估计方程线

在这里插入图片描述

方差分析

> SST=lxy(y,y) #因变量的离均差平方和
> SSR=b*lxy(x,y) #回归平方和
> SSE=SST-SSR # 误差平方和
> MSR=SSR/1 #回归均方
> MSE=SSE/(n-2) #误差均方
> F=MSR/MSE #F统计量
> c(SST=SST,SSR=SSR,SSE=SSE,MSR=MSR,MSE=MSE,F=F) #显示结果
       SST        SSR        SSE        MSR        MSE          F 
813.000000 748.173425  64.826575 748.173425   6.482657 115.411531

t检验

> sy.x=sqrt(MSE) #估计标准差
> sb=sy.x/sqrt(lxy(x,x)) #离均差平方和
> t=b/sb #t统计量
> ta=qt(1-0.05/2,n-2) #t分位数
> c(sy.x=sy.x,sb=sb,t=t,ta=ta) #显示结果
      sy.x         sb          t         ta 
 2.5461063  0.1078901 10.7429759  2.2281389

也可以使用 线性回归拟合函数lm() 来拟合线性模型
例如如下用lm()函数拟合1978-2008年税收与财政收入的数据
(1)读入数据:

> yx=read.table("clipboard",header=T)
> yx
            y        x
1978  11.3262   5.1928
1979  11.4638   5.3782
1980  11.5993   5.7170
1981  11.7579   6.2989
1982  12.1233   7.0002
1983  18.6695   7.5559
1984  16.4286   9.4735
1985  20.0482  20.4079
1986  21.2201  20.9073
1987  21.9935  21.4036
1988  23.5724  23.9047
1989  26.6490  27.2740
1990  29.3710  28.2187
1991  31.4948  29.9017
1992  34.8337  32.9691
1993  43.4895  42.5530
1994  52.1810  51.2688
1995  62.4220  60.3804
1996  74.0799  69.0982
1997  86.5114  82.3404
1998  98.7595  92.6280
1999 114.4408 106.8258
2000 133.9523 125.8151
2001 163.8604 153.0138
2002 189.0364 176.3645
2003 217.1525 200.1731
2004 263.9647 241.6568
2005 316.4929 287.7854
2006 387.6020 348.0435
2007 513.2178 456.2197
2008 613.3035 542.1962

(2)拟合模型:

> fm=lm(y~x,data=yx) #一元线性回归模型
> fm

Call:
lm(formula = y ~ x, data = yx)

Coefficients:
(Intercept)            x  
     -1.197        1.116 

(3)做回归曲线:

> plot(y~x,data=yx) #作散点图

在这里插入图片描述

> abline(fm) #添加回归线

在这里插入图片描述

(4)回归方程的假设检验
①模型的方差分析(ANOVA)

> anova(fm) #模型方差分析
Analysis of Variance Table

Response: y
          Df Sum Sq Mean Sq F value    Pr(>F)    
x          1 712077  712077   27428 < 2.2e-16 ***
Residuals 29    753      26                      
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

②回归系数的t检验

> summary(fm) #回归系数t检验

Call:
lm(formula = y ~ x, data = yx)

Residuals:
   Min     1Q Median     3Q    Max 
-6.630 -3.692 -1.535  5.338 11.432 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -1.19656    1.16125   -1.03    0.311    
x            1.11623    0.00674  165.61   <2e-16 ***
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 5.095 on 29 degrees of freedom
Multiple R-squared:  0.9989,	Adjusted R-squared:  0.9989 
F-statistic: 2.743e+04 on 1 and 29 DF,  p-value: < 2.2e-16

2.多元线性回归分析

多元线性回归模型建立

读入数据

> yx=read.table("clipboard",header=T) #加载数据
> yx
            y       x1       x2       x3     x4
1978  11.3262   36.241   5.1928    3.550 406.82
1979  11.4638   40.382   5.3782    4.120 415.92
1980  11.5993   45.178   5.7170    5.700 429.03
1981  11.7579   48.603   6.2989    8.904 441.65
1982  12.1233   53.018   7.0002   12.801 456.74
1983  18.6695   59.574   7.5559   15.903 467.07
1984  16.4286   72.067   9.4735   18.202 484.33
1985  20.0482   89.891  20.4079   20.667 501.12
1986  21.2201  102.014  20.9073   26.019 515.46
1987  21.9935  119.545  21.4036   32.202 530.60
1988  23.5724  149.223  23.9047   41.600 546.30
1989  26.6490  169.178  27.2740   49.802 557.07
1990  29.3710  185.984  28.2187   55.601 653.23
1991  31.4948  216.625  29.9017   72.258 660.91
1992  34.8337  266.519  32.9691   91.196 667.82
1993  43.4895  345.605  42.5530  112.710 674.68
1994  52.1810  466.700  51.2688  203.819 681.35
1995  62.4220  574.949  60.3804  234.999 688.55
1996  74.0799  668.505  69.0982  241.338 697.65
1997  86.5114  731.427  82.3404  269.672 708.00
1998  98.7595  769.672  92.6280  268.577 720.87
1999 114.4408  805.794 106.8258  298.963 727.91
2000 133.9523  882.281 125.8151  392.742 739.92
2001 163.8604  943.464 153.0138  421.933 744.32
2002 189.0364 1203.327 176.3645  513.782 753.60
2003 217.1525 1358.228 200.1731  704.835 760.75
2004 263.9647 1598.783 241.6568  955.391 768.23
2005 316.4929 1832.174 287.7854 1169.218 778.77
2006 387.6020 2119.235 348.0435 1409.714 782.44
2007 513.2178 2495.299 456.2197 1667.402 786.45
2008 613.3035 3006.700 542.1962 1778.898 790.48

建立多元线性回归模型

> fm=lm(y~x1+x2+x3+x4,data=yx) #显示多元线性回归模型
> fm

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

Coefficients:
(Intercept)           x1           x2           x3  
 23.5321088   -0.0033866    1.1641150    0.0002919  
         x4  
 -0.0437416  

标准化偏回归系数

> library(mvstats)
> coef.sd(fm) #标准化偏回归系数结果
$coef.sd
           x1            x2            x3            x4 
-0.0174513678  1.0423522972  0.0009628564 -0.0371053994 
多元线性回归模型检验
> summary(fm) #多元线性回归系数t检验

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

Residuals:
    Min      1Q  Median      3Q     Max 
-5.0229 -2.1354  0.3297  1.2639  6.9690 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) 23.5321088  4.5990714   5.117 2.47e-05 ***
x1          -0.0033866  0.0080749  -0.419    0.678    
x2           1.1641150  0.0404889  28.751  < 2e-16 ***
x3           0.0002919  0.0085527   0.034    0.973    
x4          -0.0437416  0.0092638  -4.722 7.00e-05 ***
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.79 on 26 degrees of freedom
Multiple R-squared:  0.9997,	Adjusted R-squared:  0.9997 
F-statistic: 2.289e+04 on 4 and 26 DF,  p-value: < 2.2e-16

3.多元线性相关分析

①矩阵相关分析
> cor(yx) #多元数据相关系数矩阵
           y        x1        x2        x3        x4
y  1.0000000 0.9871498 0.9994718 0.9912053 0.6956619
x1 0.9871498 1.0000000 0.9907018 0.9867664 0.7818066
x2 0.9994718 0.9907018 1.0000000 0.9917094 0.7154297
x3 0.9912053 0.9867664 0.9917094 1.0000000 0.7073820
x4 0.6956619 0.7818066 0.7154297 0.7073820 1.0000000

矩阵散点图函数pair()

> pairs(yx) # 多元数据散点图,yx为数值矩阵或数据框

在这里插入图片描述
相关矩阵检验函数pair()

> library(mvstats)
> corr.test(yx)  # 多元数据相关系数检验,yx为数值矩阵或数据框
corr test: 
         y     x1     x2    x3 x4
y    0.000  0.000  0.000 0.000  0
x1  33.267  0.000  0.000 0.000  0
x2 165.614 39.214  0.000 0.000  0
x3  40.336 32.772 41.560 0.000  0
x4   5.215  6.752  5.514 5.389  0
lower is t value,upper is p value
②复相关分析
> R2=summary(fm)$r.sq #显示多元线性回归模型决定系数
> R2
[1] 0.9997
> R=sqrt(R2) # 显示多元数据复相关系数
> R
[1] 0.9999

4.回归变量的选择方法

变量选择准则
> library(leaps) #加载leaps包
> varsel=regsubsets(y~x1+x2+x3+x4,data=yx) # 多元数据线性回归变量选择模型
> result=summary(varsel) #变量选择方法结果
> data.frame(result$outmat,RSS=result$rss,R2=result$rsq) #RSS和决定系数准则结果展示
         x1 x2 x3 x4   RSS     R2
1  ( 1 )     *       752.9 0.9989
2  ( 1 )     *     * 203.9 0.9997
3  ( 1 )  *  *     * 202.3 0.9997
4  ( 1 )  *  *  *  * 202.3 0.9997
> data.frame(result$outmat,adjR2=result$adjr2,Cp=result$cp,BIC=result$bic) #调整决定系数,Cp和BCI准则结果展示
         x1 x2 x3 x4  adjR2     Cp    BIC
1  ( 1 )     *       0.9989 69.745 -205.6
2  ( 1 )     *     * 0.9997  1.199 -242.6
3  ( 1 )  *  *     * 0.9997  3.001 -239.4
4  ( 1 )  *  *  *  * 0.9997  5.000 -236.0
逐步回归分析

①向前引入法

> fm=lm(y~x1+x2+x3+x4,data=yx) #多元数据线性回归模型
> fm.step=step(fm,direction = "forward") #向前引入法变量选择结果
Start:  AIC=68.15
y ~ x1 + x2 + x3 + x4

②向后剔除法

> fm.step=step(fm,direction = "backward") #向后剔除法变量选择结果
Start:  AIC=68.15
y ~ x1 + x2 + x3 + x4

       Df Sum of Sq  RSS   AIC
- x3    1         0  202  66.2
- x1    1         1  204  66.4
<none>               202  68.2
- x4    1       174  376  85.4
- x2    1      6433 6635 174.4

Step:  AIC=66.16
y ~ x1 + x2 + x4

       Df Sum of Sq  RSS   AIC
- x1    1         2  204  64.4
<none>               202  66.2
- x4    1       197  400  85.3
- x2    1      7382 7585 176.5

Step:  AIC=64.39
y ~ x2 + x4

       Df Sum of Sq    RSS   AIC
<none>                 204  64.4
- x4    1       549    753 102.9
- x2    1    367655 367859 294.8

③逐步筛选法

> fm.step=step(fm,direction = "both") #逐步筛选法变量选择结果
Start:  AIC=68.15
y ~ x1 + x2 + x3 + x4

       Df Sum of Sq  RSS   AIC
- x3    1         0  202  66.2
- x1    1         1  204  66.4
<none>               202  68.2
- x4    1       174  376  85.4
- x2    1      6433 6635 174.4

Step:  AIC=66.16
y ~ x1 + x2 + x4

       Df Sum of Sq  RSS   AIC
- x1    1         2  204  64.4
<none>               202  66.2
+ x3    1         0  202  68.2
- x4    1       197  400  85.3
- x2    1      7382 7585 176.5

Step:  AIC=64.39
y ~ x2 + x4

       Df Sum of Sq    RSS   AIC
<none>                 204  64.4
+ x1    1         2    202  66.2
+ x3    1         0    204  66.4
- x4    1       549    753 102.9
- x2    1    367655 367859 294.8

参考致谢:
王斌会.多元统计分析及R语言建模(第四版)

如有侵权,请联系侵删
需要本实验源数据及代码的小伙伴请联系QQ:2225872659

  • 2
    点赞
  • 44
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值