[(株)数理システム]

(S-PLUS) S-PLUS mini course 第31回


S-PLUS mini course の31回目をお届けします。

さて、前回は関数 step を使った変数選択についてご説明しました。
今回はその step の細かい点について、触れたいと思います。

・step は総当りではない
前回もちょっと触れましたが、総当りをして、AIC をすべて比較するということは行いません。
例えば、前回のメールでのコマンド

options(contrasts=c("contr.treatment", "contr.poly"))
# ダミー行列を分かりやすいように、通常の1,0 にかえる
fitall <- lm(Fuel ~(Weight + Disp. + Type)^2, data=fuel.frame)
# Weight, Disp. Type(因子)とその2次の交互作用項までを含む回帰を行う

であれば、可能なモデルは
主効果のみ:
変数1個 3モデル
変数2個 3*2/2=3モデル
変数3個 1モデル、合計7モデルが考えられる

主効果と交互作用:
上の全モデル(7)に、2次の交互作用の組み合わせ(多分、7通り)をかけた、最大49モデルについて、回帰を行い、AIC を求めるということをしなければなりません。変数が3でこの程度ですから、変数が増えると大変面倒なことになってしまいます。

そこで、S-PLUS は、次にどの変数を減らしたらいいか、増やしたらいいか、を統計量 Cp を用いて調べ、その変数を増やした(あるいは減らした)ものと、元のモデルのAICを比べて、その変数の増減によって新しいモデルのAICが大きくなるようなら、その1つ前のモデルでストップします。

では、

fits <- step(fitall)

では、何をしているでしょうか?これは、まず、2次のフルモデルなので、変数増減が許されていても、減らす方しか、考えられません。そのとき、落とすことができるのは、この場合ならば、2次の交互作用項のみです(主効果から落ちて、交互作用だけが残る変数は想定していないためです)。
そのため、まず、Weight:Type(交互作用)を落とす、ということが試みされ元のモデルより AIC が小さくなりますから、さらにモデル増減を進めます。

次に落とすことができるのは、残りの2つの交互作用項で、Weight:Disp.を落とす、ということがやはり試されます。
すると、また AIC が小さくなりますから、さらに進めます。
今度は、可能性として、残った交互作用 Type:Disp. と、先ほどの変数選択で交互作用が全部落ちたWeight の主効果も候補になります。
そして、主効果 Weight を落としたらよいという結論に落ち着き、続いてnone と Disp.:Type を比較すると、none のほうがよいと推測されるので、もう変数を落とすことをしません。つまり、ここで収束します。

デフォルトで変数増減法(direction="both") なのに、どうして、変数を増やすほうがされていないか、この場合ならば、「主効果のない、交互作用項が考えられていない」ためです。お手元のデータで変数を増やして試していただけると、多分、増減が発生していると思います。

次回はこの変数選択の続きと、他のモデル(樹形モデル)を使った変数選択についてご紹介したいと思います。
お楽しみに!

----------------------------------------------------------------

? お問い合わせ/資料請求はお気軽に。
(株)数理システム《S-PLUS》グループ <splus-info@msi.co.jp>