S-PLUS mini course の23回目をお届けします。
まず、お知らせです。
*******************************************************
新しいバージョン、S-PLUS 6.1 for Windows がリリースされました!
詳しくは
http://www.msi.co.jp/splus/
をご覧ください。
(ユーザーの皆様には順次、バージョンアップお申込書をお送りしています)
*******************************************************
前回は線形回帰結果を用いて、予測を行う方法をご紹介しました。そのときには、クラス lm が持っている、メソッド predict.lm が呼び出されるのでした。
今回は、説明変数に因子データがあったらどうなるか、ご説明したいと思います。
例えば、ずっと取り上げているサンプルデータセット fuel.frame には
> names(fuel.frame)
[1] "Weight" "Disp." "Mileage" "Fuel" "Type"
という変数があり、このうち、Type という変数は
> class(fuel.frame$Type)
[1] "factor"
と、因子(factor)です。例えば、
> table(fuel.frame$Type)
Compact Large Medium Small Sporty Van
15
3 13 13
9 7
とすると、各因子水準の度数分布を調べることができます。(単に因子水準を調べるには、levels(fuel.frame$Type) としてください)
では、この因子を説明変数にしたいときはどうしたらよいでしょうか?ソフトによっては、ダミー変数行列を作る必要がありますね。
(if Type="Compact" then mat1=1 else mat1=0 みたいなプログラムを作っていた記憶があります)
S-PLUS では、この変数が因子であれば、説明変数に追加された段階でダミー行列化されます。
例えば
> fit3 <- lm(Fuel ~ Weight + Type, data=fuel.frame)
# 変数 Type は因子
とすると、結果は
> fit3
Call:
lm(formula = Fuel ~ Weight + Type, data = fuel.frame)
Coefficients:
(Intercept) Weight Type1 Type2
Type3
1.686876 0.0008846407 0.02159053 0.02687706 -0.1111044
Type4
Type5
-0.02557048 0.1029998
Degrees of freedom: 60 total; 53 residual
Residual standard error: 0.3634024
と、6水準の因子 Type に対して、5個の係数がありますから、確かにダミー変数が使われていることが分かります。
ですが、Type1 がどの水準に対応しているのか(Compact か、Largeか、・・・)はここでは分かりません。これを調べるには
> fit3$contrasts
$Type:
[,1]
[,2] [,3] [,4] [,5]
Compact -1 -1 -1 -1 -1
Large 1 -1 -1
-1 -1
Medium 0 2 -1
-1 -1
Small 0 0
3 -1 -1
Sporty 0 0
0 4 -1
Van 0
0 0 0 5
とします。すると、おなじみの1,0 で表現される行列ではありませんね。何も指定しないと、互いに直交するヘルマート対比(Helmert Contrats) が使われます。
これは、Type1 は「第1水準と第2水準の差」、Type2 は「第3水準と第1第2水準の平均との差」、を求めるものです。
(因子水準の順序は、levels(fuel.frame$Type) でわかります。特に
指定しなかった場合は、アルファベット順となります)
なんか、わかりにくいなー、と感じた方、次回をお楽しみに。
次回は因子変数のコーディング方法を変更する方法をご紹介します。