<< トップページへ

<< 目次へ

最終更新日: 2008.11.22

lm()におけるパラメータ制約条件



 以下に示す内容はglm()でも同じこと。

 Rのlm()では独立変数にカテゴリカル型の変数を指定すると、それはダミー変数として解析されることになるが、その際に標準では最初の変数(列)が削除されて解析される。要するにダミー変数の最初の変数をベースラインとして計算されるわけだ。

 これに対してSASでは最後の変数をベースラインとして計算し、Minitabでは合計を0として(ベースラインとして)計算している。RではSASやMinitabと同じ形式で結果を得ることができる。その方法を以下に示す。

> Y <- rnorm(18)
> Group <- rep(c(1, 2, 3), c(6, 6, 6))
> Group <- as.factor(Group)
> Group
 [1] 1 1 1 1 1 1 2 2 2 2 2 2 3 3 3 3 3 3
Levels: 1 2 3

# Rの標準時
> result <- lm(Y ~ Group)

# SASと同じように
> result2 <- lm(Y ~ Group, contrasts=list(Group="contr.SAS"))

# Minitabと同じように
> result3 <- lm(Y ~ Group, contrasts=list(Group="contr.sum"))

# 係数表の出力
# Rでは最初の変数がベースライン
> dummy.coef(result)
Full coefficients are 
                                                      
(Intercept):      0.06810832                          
Group:                     1            2            3
                 0.000000000 -0.413742210 -0.002882025

# SASでは最後の列がベースライン
> dummy.coef(result2)
Full coefficients are 
                                                      
(Intercept):       0.0652263                          
Group:                     1            2            3
                 0.002882025 -0.410860185  0.000000000

# Minitabでは合計が0になるようにして計算している
> dummy.coef(result3)
Full coefficients are 
                                                 
(Intercept):    -0.07076642                      
Group:                    1          2          3
                  0.1388747 -0.2748675  0.1359927


 ちなみに、Minitabの"合計が0になるように"というのは"ダミー変数に置き換えた計画行列において各列の合計が0になるように"ということである。これは、contr.SAS()とcontr.sum()という関数を使って実際に見ているとよいだろう。

# SASでは最後の列が削除されている
# 数えてみると"3"の列数は5つであることが確認できる
> contr.SAS(Group)
  1 1 1 1 1 1 2 2 2 2 2 2 3 3 3 3 3
1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

# このように各列の合計が0になるような計画行列が作成されている
> contr.sum(Group)
  [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
1    1    0    0    0    0    0    0    0    0     0     0     0     0
1    0    1    0    0    0    0    0    0    0     0     0     0     0
1    0    0    1    0    0    0    0    0    0     0     0     0     0
1    0    0    0    1    0    0    0    0    0     0     0     0     0
1    0    0    0    0    1    0    0    0    0     0     0     0     0
1    0    0    0    0    0    1    0    0    0     0     0     0     0
2    0    0    0    0    0    0    1    0    0     0     0     0     0
2    0    0    0    0    0    0    0    1    0     0     0     0     0
2    0    0    0    0    0    0    0    0    1     0     0     0     0
2    0    0    0    0    0    0    0    0    0     1     0     0     0
2    0    0    0    0    0    0    0    0    0     0     1     0     0
2    0    0    0    0    0    0    0    0    0     0     0     1     0
3    0    0    0    0    0    0    0    0    0     0     0     0     1
3    0    0    0    0    0    0    0    0    0     0     0     0     0
3    0    0    0    0    0    0    0    0    0     0     0     0     0
3    0    0    0    0    0    0    0    0    0     0     0     0     0
3    0    0    0    0    0    0    0    0    0     0     0     0     0
3   -1   -1   -1   -1   -1   -1   -1   -1   -1    -1    -1    -1    -1
  [,14] [,15] [,16] [,17]
1     0     0     0     0
1     0     0     0     0
1     0     0     0     0
1     0     0     0     0
1     0     0     0     0
1     0     0     0     0
2     0     0     0     0
2     0     0     0     0
2     0     0     0     0
2     0     0     0     0
2     0     0     0     0
2     0     0     0     0
3     0     0     0     0
3     1     0     0     0
3     0     1     0     0
3     0     0     1     0
3     0     0     0     1
3    -1    -1    -1    -1


 さらにちなみに、Rのlm()の標準ではcontr.treatment()が採用されている。だから、lm()の引数にcontrasts=list(Group="contr.treatment")と指定したものと同じになる。

# SASとは逆に、"1"の列が5つである
# 先頭の列が削除されている(=ベースラインにされている)
> contr.treatment(Group)
  1 1 1 1 1 2 2 2 2 2 2 3 3 3 3 3 3
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
2 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1

# ためしにわざわざcontr.treatmentとして実行してみる
> result4 <- lm(Y ~ Group, contrasts=list(Group="contr.treatment"))

# 係数を見てみると、resultと全く同じであることがわかる
> dummy.coef(result4)
Full coefficients are 
                                                      
(Intercept):      0.06810832                          
Group:                     1            2            3
                 0.000000000 -0.413742210 -0.002882025

> dummy.coef(result)
Full coefficients are 
                                                      
(Intercept):      0.06810832                          
Group:                     1            2            3
                 0.000000000 -0.413742210 -0.002882025


 まぁ、こんなことを知っていも実務的にはあまり役立たないかもしれない(^_^;)

 さらに、さらにちなみに、contrast()という関数を利用すると簡単にダミー変数行列を作成することができる(上のような計画行列ではなくて)。

> contrasts(Group)
  2 3
1 0 0
2 1 0
3 0 1