r语言做线性回归模型（用R语言做数据分析）

在实际问题中，影响因变量Y的因素很多，人们可以从中挑选若干个变量建立回归方程，这便涉及变量选择的问题。

一般来讲，如果在一个回归方程中忽略了对Y有显著影响的自变量，那么所建立的方程必与实际有较大的偏离，但变量选得过多，使用就不方便，特别是方程中含有对Y影响不大的变量时。因此适当地选择变量以建立一个“最优”的回归方程是十分重要的。

什么是“最优”回归方程呢？对于这个问题有许多不同的准则，在不同的准则下，“最优”回归方程也可能不同，这讲的“最优”是指从可供选择的所有变量中选出对Y有显著影响的变量建立方程，且在方程中不含有对Y无显著影响的变量。

有许多方法准则可以获得“最优”回归方程，如“一切子集回归法”、"前进法"、“后退法”、“逐步回归法”等，其中“逐步回归法”由于计算程序简便，因而使用较为普遍。

R语言提供了较为方便的“逐步回归”计算函数step(),它是以AIC信息统计量为准则，通过选择最小的AIC信息统计量，来达到删除或增加变量的目的。step()函数的调用格式为：

其中object是回归模型，scope是确定逐步搜索的区域，scale是用于AIC统计量，direction确定逐步搜索的方向，缺省值为“both”，是“一切子集回归法”，“backward”是“后退法”，“forward”是“前进法”。下面通过一个简单的例子，介绍如何使用R语言来完成逐步回归的过程，从而达到选择“最优”方程的目的。

例子：某种水泥在凝固时放出的热量Y与水泥中四种化学成分X1、X2、X3、X4有关，现测得13组数据，求从中选出主要的变量，建立Y关于它们的线性回归方程。

首先作多元线性回归方程

> X1=c( 7, 1, 11, 11, 7, 11, 3, 1, 2, 21, 1, 11, 10);

> X2=c(26, 29, 56, 31, 52, 55, 71, 31, 54, 47, 40, 66, 68);

> X3=c( 6, 15, 8, 8, 6, 9, 17, 22, 18, 4, 23, 9, 8);

> X4=c(60, 52, 20, 47, 33, 22, 6, 44, 22, 26, 34, 12, 12);

> Y =c(78.5, 74.3, 104.3, 87.6, 95.9, 109.2, 102.7, 72.5,93.1,115.9,83.8,113.3,109.4);

> cement<-data.frame(X1,X2,X3,X4,Y)

> lm.sol<-lm(Y~X1 X2 X3 X4,data=cement)

> summary(lm.sol)

Call:

lm(formula = Y ~ X1 X2 X3 X4, data = cement)

Residuals:

Min 1Q Median 3Q Max

-3.1750 -1.6709 0.2508 1.3783 3.9254

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 62.4054 70.0710 0.891 0.3991

X1 1.5511 0.7448 2.083 0.0708 .

X2 0.5102 0.7238 0.705 0.5009

X3 0.1019 0.7547 0.135 0.8959

X4 -0.1441 0.7091 -0.203 0.8441

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.446 on 8 degrees of freedom

Multiple R-squared: 0.9824, Adjusted R-squared: 0.9736

F-statistic: 111.5 on 4 and 8 DF, p-value: 4.756e-07

从上述计算中可以看到，如果选择全部变量作回归方程，效果是不好的，因为回归方程的系数没有一项通过检验。接下来使用step()作逐步回归：

> lm.step<-step(lm.sol)

Start: AIC=26.94

Y ~ X1 X2 X3 X4

Df Sum of Sq RSS AIC

- X3 1 0.1091 47.973 24.974

- X4 1 0.2470 48.111 25.011

- X2 1 2.9725 50.836 25.728

<none> 47.864 26.944

- X1 1 25.9509 73.815 30.576

Step: AIC=24.97

Y ~ X1 X2 X4

Df Sum of Sq RSS AIC

<none> 47.97 24.974

- X4 1 9.93 57.90 25.420

- X2 1 26.79 74.76 28.742

- X1 1 820.91 868.88 60.629

从程序运行结果可以看到，用全部变量作回归方程时，AIC值为26.94。如果去掉变量X3，得到的回归方程的AIC值为24.974，如果去掉变量X4，得到的回归方程的AIC值为25.011。后面类推，由于去掉变量X3可以使AIC达到最小，因此R语言自动去掉变量X3，进行下一轮计算。

> summary(lm.step)

Call:

lm(formula = Y ~ X1 X2 X4, data = cement)

Residuals:

Min 1Q Median 3Q Max

-3.0919 -1.8016 0.2562 1.2818 3.8982

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 71.6483 14.1424 5.066 0.000675 ***

X1 1.4519 0.1170 12.410 5.78e-07 ***

X2 0.4161 0.1856 2.242 0.051687 .

X4 -0.2365 0.1733 -1.365 0.205395

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.309 on 9 degrees of freedom

Multiple R-squared: 0.9823, Adjusted R-squared: 0.9764

F-statistic: 166.8 on 3 and 9 DF, p-value: 3.323e-08

由显示结果看到：回归系数检验的显著性水平很很大提高，但变量X2、X4系数检验的显著性水平仍然不理想。下面该如何处理呢？

R语言还有两个函数可以用来作逐步回归，这两个函数是add1()和drop1()，它们的调用格式为：

其中object是由拟合模型构成的对象，scope是模型考虑增加或去掉项构成的公式，scale是用于计算残差的均方估计值，缺省值为0或NULL。下面用drop1()函数计算。

> drop1(lm.step)

Single term deletions

Model:

Y ~ X1 X2 X4

Df Sum of Sq RSS AIC

<none> 47.97 24.974

X1 1 820.91 868.88 60.629

X2 1 26.79 74.76 28.742

X4 1 9.93 57.90 25.420

从运算结果来看，如果去掉变量X4，AIC值会从24.974增加到25.420，是增加最少的。另外，除了AIC准则外，残差的平方和也是逐步回归的重要指标之一。从直观上来看，拟合越好的方程，残差的平方和应越小，去掉X4，残差的平方和上升9.93，也是最少的，因此从这两项指标来看，应该去掉变量X4。

> lm.opt<-lm(Y~X1 X2,data=cement);

> summary(lm.opt);

Call:

lm(formula = Y ~ X1 X2, data = cement)

Residuals:

Min 1Q Median 3Q Max

-2.893 -1.574 -1.302 1.363 4.048

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 52.57735 2.28617 23.00 5.46e-10 ***

X1 1.46831 0.12130 12.11 2.69e-07 ***

X2 0.66225 0.04585 14.44 5.03e-08 ***

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.406 on 10 degrees of freedom

Multiple R-squared: 0.9787, Adjusted R-squared: 0.9744

F-statistic: 229.5 on 2 and 10 DF, p-value: 4.407e-09

这个结果应该是满意的，因为所有的检验都是显著的，最后得到“最优”回归方程为：

Y= 52.58 1.468*X1 0.6623*X2

,