# モデル選択_実践編 # 馬場真哉 # 2015年9月2日 # シミュレーションデータの作成 set.seed(0) N <- 100 # サンプルサイズ Intercept <- 5 # 切片 B1 <- 10 # 係数1 B2 <- 5 # 係数2 x1 <- sort(rnorm(N, sd=2)) # 説明変数1 x2 <- rnorm(N, sd=2) # 説明変数2 e <- rnorm(n=N, sd=3) # 誤差 y <- Intercept + B1*x1 + B2*x2 + e # 普通に実行すると、毎回値が変わる rnorm(10) rnorm(10) rnorm(10) # set.seed(0)を設定すると値が固定される set.seed(0) rnorm(10) set.seed(0) rnorm(10) set.seed(0) rnorm(10) # モデルの作成 # 説明変数が1つだけ model1 <- lm(y ~ x1) model1 # 説明変数2つ model2 <- lm(y ~ x1 + x2) model2 # x1とx2の交互作用も入れてモデリング model3 <- lm(y ~ x1*x2) model3 # summary関数による検定 summary(model3) # anova関数による検定 anova(model2, model3) # モデル2とモデル3の比較 anova(model1, model2) # モデル1とモデル2の比較 # AICによるモデル選択 AIC(model1) AIC(model2) AIC(model3) # step関数によるモデル選択 model.best1 <- step(model3) model.best1 # MuMInパッケージを用いた総当たり法によるモデル選択 library(MuMIn) options(na.action = "na.fail") kekka.AIC <- dredge(model3, rank="AIC") kekka.AIC # 最適なモデルを一つ選ぶ best.model <- get.models(kekka.AIC, subset = 1)[1] best.model # Model Averaging avg.model <- model.avg(get.models(kekka.AIC, subset = delta < 4)) avg.model summary(avg.model)