R でステップワイズ回帰

ステップワイズ回帰は、複数の説明変数の候補がある場合にその中から有効な説明変数を絞り込んでいく方法です。Rのstep()関数はAICを基準として、AICのより小さい説明変数の組合せを選びます。

まず被説明変数と説明変数のデータを準備します。ここではsdata.txtとして1列目に被説明変数、２～11列目に説明変数の候補が入っているデータを作成し、Dドライブの直下に保存しました。これをRのxというオブジェクトに読み込みます。xは仮にxとしているだけで、何でもいいです。

「header=T」というのは1行目に列名のヘッダーがあるという意味、「sep="\t"」はデータがタブで区切られているという意味です。

全部の変数を使った重回帰 lm()を行い、結果を「res1」に入れています。

全変数を使った重回帰の結果 res1のサマリー
********************************************************
summary(res1)

Call:
lm(formula = x[, 1] ~ x[, 2] + x[, 3] + x[, 4] + x[, 5] + x[,
    6] + x[, 7] + x[, 8] + x[, 9] + x[, 10] + x[, 11])

Residuals:
     Min       1Q   Median       3Q      Max
-2.14407 -0.45743  0.04198  0.39069  1.45638

Coefficients:
             Estimate Std. Error t value Pr(>|t|)  
(Intercept) -0.129139   0.116580  -1.108  0.27308  
x[, 2]       0.770721   0.327115   2.356  0.02227 *
x[, 3]      -0.167470   0.276670  -0.605  0.54761  
x[, 4]      -0.041000   0.162924  -0.252  0.80230  
x[, 5]      -0.077808   0.125503  -0.620  0.53799  
x[, 6]       0.053005   0.046636   1.137  0.26093  
x[, 7]       0.055532   0.046713   1.189  0.23993  
x[, 8]      -0.022178   0.054304  -0.408  0.68466  
x[, 9]       0.441864   0.093823   4.710 1.89e-05 ***
x[, 10]     -0.005901   0.090956  -0.065  0.94852  
x[, 11]      0.614947   0.219764   2.798  0.00719 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.7848 on 52 degrees of freedom
Multiple R-squared:  0.7796,    Adjusted R-squared:  0.7372
F-statistic: 18.39 on 10 and 52 DF,  p-value: 8.854e-14
********************************************************
この時点でp値の小さい順に9、11、2列目の変数がsignificantとなっていますね。

res1をstep()関数の引数にしてステップワイズ回帰を行い、結果をres2に返しています。

summary(res2)

res2のサマリー結果を示しています。

結果としては、2，9，11列目の変数、つまり1，8，10番目の変数を用いた重回帰モデルが選ばれています。

途中ではp値の大きい順に外してAICを計算し、AICが最も小さくなる変数の組合せを選んでいるようです。

以下がRの出力結果の全体の記録です。

*************************************************************

Start:  AIC=-20.62

x[, 1] ~ x[, 2] + x[, 3] + x[, 4] + x[, 5] + x[, 6] + x[, 7] + 

    x[, 8] + x[, 9] + x[, 10] + x[, 11]

          Df Sum of Sq    RSS      AIC

- x[, 10]  1    0.0026 32.033 -22.6111

- x[, 4]   1    0.0390 32.069 -22.5395

- x[, 8]   1    0.1027 32.133 -22.4144

- x[, 3]   1    0.2257 32.256 -22.1738

- x[, 5]   1    0.2368 32.267 -22.1522

- x[, 6]   1    0.7957 32.826 -21.0702

- x[, 7]   1    0.8705 32.901 -20.9269

none   　　             32.030 -20.6162

- x[, 2]   1    3.4194 35.450 -16.2259

- x[, 11]  1    4.8230 36.853 -13.7796

- x[, 9]   1   13.6621 45.693  -0.2356

Step:  AIC=-22.61

x[, 1] ~ x[, 2] + x[, 3] + x[, 4] + x[, 5] + x[, 6] + x[, 7] + 

    x[, 8] + x[, 9] + x[, 11]

          Df Sum of Sq    RSS      AIC

- x[, 4]   1    0.0418 32.075 -24.5290

- x[, 8]   1    0.1002 32.133 -24.4143

- x[, 3]   1    0.2278 32.261 -24.1645

- x[, 5]   1    0.2605 32.294 -24.1008

- x[, 6]   1    0.7977 32.831 -23.0615

- x[, 7]   1    0.9214 32.955 -22.8244

  none   　　            32.033 -22.6111

- x[, 2]   1    3.7092 35.742 -17.7084

- x[, 11]  1    4.8964 36.929 -15.6500

- x[, 9]   1   13.9350 45.968  -1.8568

Step:  AIC=-24.53

x[, 1] ~ x[, 2] + x[, 3] + x[, 5] + x[, 6] + x[, 7] + x[, 8] + 

    x[, 9] + x[, 11]

          Df Sum of Sq    RSS      AIC

- x[, 8]   1    0.0948 32.170 -26.3430

- x[, 5]   1    0.2443 32.319 -26.0510

- x[, 3]   1    0.3242 32.399 -25.8955

- x[, 6]   1    0.8360 32.911 -24.9081

- x[, 7]   1    1.0115 33.086 -24.5729

none   　　              32.075 -24.5290

- x[, 2]   1    3.9902 36.065 -19.1421

- x[, 11]  1    6.2805 38.355 -15.2632

- x[, 9]   1   14.7259 46.801  -2.7259

Step:  AIC=-26.34

x[, 1] ~ x[, 2] + x[, 3] + x[, 5] + x[, 6] + x[, 7] + x[, 9] + 

    x[, 11]

          Df Sum of Sq    RSS     AIC

- x[, 5]   1    0.2227 32.392 -27.908

- x[, 3]   1    0.2777 32.447 -27.802

- x[, 6]   1    0.7732 32.943 -26.847

- x[, 7]   1    0.9772 33.147 -26.458

none   　　              32.170 -26.343

- x[, 2]   1    4.0497 36.219 -20.873

- x[, 11]  1    6.4312 38.601 -16.861

- x[, 9]   1   28.2913 60.461  11.408

Step:  AIC=-27.91

x[, 1] ~ x[, 2] + x[, 3] + x[, 6] + x[, 7] + x[, 9] + x[, 11]

          Df Sum of Sq    RSS     AIC

- x[, 3]   1     0.468 32.861 -29.004

- x[, 6]   1     0.626 33.018 -28.703

none   　　             32.392 -27.908

- x[, 7]   1     1.236 33.629 -27.549

- x[, 2]   1     3.983 36.375 -22.603

- x[, 11]  1     6.322 38.715 -18.676

- x[, 9]   1    34.545 66.937  15.819

Step:  AIC=-29

x[, 1] ~ x[, 2] + x[, 6] + x[, 7] + x[, 9] + x[, 11]

          Df Sum of Sq    RSS     AIC

- x[, 6]   1     0.326 33.187 -30.382

none   　　             32.861 -29.004

- x[, 7]   1     1.233 34.094 -28.683

- x[, 2]   1     3.533 36.394 -24.570

- x[, 11]  1     5.854 38.715 -20.675

- x[, 9]   1    41.871 74.732  20.759

Step:  AIC=-30.38

x[, 1] ~ x[, 2] + x[, 7] + x[, 9] + x[, 11]

          Df Sum of Sq    RSS     AIC

- x[, 7]   1     0.934 34.121 -30.633

none   　　             33.187 -30.382

- x[, 2]   1     4.661 37.847 -24.103

- x[, 11]  1     7.499 40.686 -19.547

- x[, 9]   1    41.740 74.927  18.923

Step:  AIC=-30.63

x[, 1] ~ x[, 2] + x[, 9] + x[, 11]

          Df Sum of Sq    RSS     AIC

none   　　             34.121 -30.633

- x[, 2]   1     4.136 38.257 -25.425

- x[, 11]  1     6.695 40.816 -21.346

- x[, 9]   1    47.382 81.502  22.222

summary(res2)

Call:

lm(formula = x[, 1] ~ x[, 2] + x[, 9] + x[, 11])

Residuals:

     Min       1Q   Median       3Q      Max 

-2.11422 -0.51526  0.02026  0.42524  1.65108 

Coefficients:

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

(Intercept) -0.13892    0.10188  -1.364  0.17790    

x[, 2]       0.65220    0.24388   2.674  0.00967 ** 

x[, 9]       0.37875    0.04184   9.052 9.38e-13 ***

x[, 11]      0.58000    0.17046   3.403  0.00120 ** 

---

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

Residual standard error: 0.7605 on 59 degrees of freedom

Multiple R-squared:  0.7652,    Adjusted R-squared:  0.7533 

F-statistic:  64.1 on 3 and 59 DF,  p-value: ＜ 2.2e-16