統計講座

ダミー変数を用いた重回帰分析

被験者諸属性と試験結果(架空データ)を用いる。

ある集団から無作為に選んだ326人に学力テストを行い、その点数を記録し、さらに受験者のさまざまな属性と合わせて結果を分析する。得点と身長・年齢・性別の関係について回帰分析を行う。

説明変数に性別(男性・女性の2値)というカテゴリー変数が加わった。

ダミー変数

カテゴリー変数を回帰分析に組み入れるためには、数値化しなければならない。このとき変数内の値(水準)を独立した変数とし、1,0の2値に変換する。新たに作られる変数は水準数 - 1である。すべての水準を変数にしてしまうと一つの変数がその他の変数で決定可能になってしまい、独立性をもてない（分析が出来なくなる）。

信号赤黄
赤10
黄01
青00

「青」は「赤」でも「黄」でもないものとして決定するので、「青」変数は作らない。

今回「性別」は「男」「女」の2水準なので、新たに作る変数は一つ（「男」ないし「女」）である。

ダミー変数の作成は原理は単純だが、Rでは手間がかかるので、関数を利用する。

data <- read.csv("http://kyoto-edu.sakura.ne.jp/weblesson/statistics/data/exam-result.csv", fileEncoding = "utf-8")
summary(data)
source("http://kyoto-edu.sakura.ne.jp/weblesson/statistics/socialStatisticsBasic.R", encoding="UTF-8")
dummy.data <- recode.dummy.variables(data)

     no marks height  age gender.male
1     1    50 122.06  8.6           1
2     2    76 153.21 12.0           0
3     3    43 117.53  7.2           1
4     4    51 127.56  8.4           0

recode.dummy.variables関数(socialStatisticsBasic.Rで読み込まれる自作関数)

recode.dummy.variables(data,ref)

データフレーム中のカテゴリー変数をダミー変数に変換する

  • data=データフレーム
  • ref="元変数名.水準名"（ダミー変数から除外する値=基準値、"none"を指定すると全変数を出力する）

ダミー変数作成例

salesweekevent
1400水
1180木
…
1430木
1890金特売
2010土特売

土曜日を基準にしたい時には現変数名weekにweek内にある値"土"を"."(ドット）で連結して引数とする。

この指定をしない時には初期に設定されている最初のカテゴリを基準カテゴリとする。

sales.data <- read.csv("http://kyoto-edu.sakura.ne.jp/weblesson/statistics/data/weekly-sales.csv", fileEncoding = "utf-8")
summary(sales.data)
dummy.sales.data <- recode.dummy.variables(sales.data,ref=c("week.土"))

   sales week.火 week.金 week.月 week.水 week.日 week.木 event.特売
1   1400       0       0       0       1       0       0          0
2   1180       0       0       0       0       0       1          0
3   1650       0       1       0       0       0       0          0
4   1800       0       0       0       0       0       0          0
5   2140       0       0       0       0       1       0          0
6   1380       0       0       1       0       0       0          0
7   1420       1       0       0       0       0       0          0
8   1320       0       0       0       1       0       0          0
9   1430       0       0       0       0       0       1          0
10  1890       0       1       0       0       0       0          1
11  2010       0       0       0       0       0       0          1

これで回帰分析に掛けると「土曜日と比較して日曜日は、月曜日は、火曜日は、…」と回帰係数を読み取ることになる。

ダミー変数を用いた回帰分析

ダミー変数化されたdummy.dataを使って、得点に関して身長・年齢・性別で重回帰分析を行う。

summary(lm(marks~height+age+gender.male,dummy.data))

Call:
lm(formula = marks ~ height + age + gender.male, data = dummy.data)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.6812 -2.3536  0.0504  2.2116  6.2045 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -3.568386   2.509474  -1.422    0.156    
height      -0.009409   0.031459  -0.299    0.765    
age          6.711781   0.205423  32.673  < 2e-16 ***
gender.male -1.854224   0.318831  -5.816 1.46e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.871 on 322 degrees of freedom
Multiple R-squared:  0.9589,    Adjusted R-squared:  0.9585 
F-statistic:  2501 on 3 and 322 DF,  p-value: < 2.2e-16

※R標準の回帰分析のlm関数はカテゴリー変数を説明変数として投入すると自動的にダミー変数化して分析を行う。

lm(marks~height+age+gender,data)

同じ結果が出力される。ただしカテゴリー変数が3つ以上の水準を持つ時、どれを基準とするかの指定は少し面倒。

# データの例
data <- data.frame(
  outcome = c(10, 20, 15, 25, 30, 40),
  group = factor(c("A", "B", "A", "B", "C", "C"))
)

# 基準変数を "C" に変更
data$group <- relevel(data$group, ref = "C")
身長・年齢・性別を説明変数とする得点に対する重回帰式
y(得点) = -3.568386386 + -0.009409179*x1(身長) + 6.711780834*x2(年齢) + -1.854224321*x3(性別.男)
モデル評価
  • 決定係数(寄与率）Multiple R-squared: 0.9589
  • 調整済み決定係数 Adjusted R-squared: 0.9585
分散分析
  • 分散分析 F(3,322)=2501, p < 0.01
回帰係数
  • 身長の係数が負の値なので身長が高ければ得点は下がる(身長1cmに付き得点は-0.01点下がる)。
  • 年齢の係数が正の値なので年齢が高ければ得点も上がる(年齢1歳に付き得点は6.71点上がる)。
  • 性別.男の係数が負の値なので男性ならば得点は下がる(男性であれば得点は-1.85点下がる)。
回帰係数の有意性
  • 年齢 t(322) = 32.67, p < 0.01で1%水準で有意。
  • 身長 t(322) = -0.30, p = 0.77で5%水準で有意ではない。
  • 性別 t(322) = -5.81, p < 0.01で1%水準で有意。

linest関数(socialStatisticsBasic.Rで読み込まれる自作関数)

標準化回帰係数やVIFは標準のlm関数では出力されない。そこで自作関数linestを用いる。

linest(fomula,data)

重回帰分析を行う（R標準のlm関数の拡張)。

  • formula=回帰式(目的変数~説明変数1 + 説明変数2)
  • data=データフレーム
data <- read.csv("http://kyoto-edu.sakura.ne.jp/weblesson/statistics/data/exam-result.csv", fileEncoding = "utf-8")
summary(data)
source("http://kyoto-edu.sakura.ne.jp/weblesson/statistics/socialStatisticsBasic.R", encoding="UTF-8")
linest(marks~height+age+gender,data)

$summary
            Model Summary
R               0.9792101
R2              0.9588523
Adjusted R2     0.9584690
SE              2.8706897

$anova
                  SS  df           MS        F  P
Regression 61835.072   3 20611.690710 2501.158  0
Residual    2653.557 322     8.240859       NA NA
Total      64488.629 325   198.426550       NA NA

$coefficient
                       b         SE         beta          t             P            F      VIF
(Intercept) -3.568386386 2.50947440           NA -1.4219656  1.560041e-01 2.021986e+00       NA
height      -0.009409179 0.03145873 -0.009016326 -0.2990959  7.650598e-01 8.945836e-02 7.111298
age          6.711780834 0.20542316  0.984178271 32.6729513 2.937583e-104 1.067522e+03 7.100378
gendermale  -1.854224321 0.31883060 -0.065917315 -5.8157038  1.455399e-08 3.382241e+01 1.005322
標準化回帰係数beta
  • 年齢の絶対値が0.98で最大であり、身長の絶対値が0.01で最小である。
VIF
  • VIFはどの説明変数も10未満。