myscratchpad

dat <- read.table("http://blog.cnobi.jp/v1/blog/user/89d80905c7038b4121822249e9062fba/1255449130",header = T)

## 普通の重回帰
dat$college <- factor(dat$college,levels=c("LA","E","BUS"))
dat$cGPA <- dat$GPA - mean(dat$GPA) # 連続変量の中心化
res.idv <- lm(salary ~ cGPA+college, dat) # 交互作用なしモデル
res.int <- lm(salary ~ cGPA * college, dat) # 交互作用ありモデル
summary(res.idv)
summary(res.int)


## 標準化偏回帰係数を求める。ダミー変数をつくり、交互作用項を計算し、全部標準化して同じ分析をするだけ
# ダミー変数の作成。青木先生のコードを参考にさせていただいた。感謝。
lvs <- levels(dat$college) # 水準名を確認し、オブジェクトに格納する
la.dm <- as.integer(dat$college==lvs[1])
e.dm <- as.integer(dat$college==lvs[2])
bus.dm <- as.integer(dat$college==lvs[3])

dat2 <- data.frame(la.dm, e.dm, bus.dm, dat[3:5])
dat3 <- data.frame(scale(dat2))　# 変数の標準化

## 交互作用なしモデルの標準化偏回帰係数
res.idv2 <- lm(scale(salary)~scale(GPA)+e.dm+bus.dm+la.dm, dat2)
res.idv3 <- lm(salary~GPA+e.dm+bus.dm+la.dm, dat3)
summary(res.idv2) # ダミー変数以外を標準化
summary(res.idv3) # 全て標準化
## 交互作用ありモデルの標準化偏回帰係数
# ダミー変数以外を標準化する。ただし、この場合ダミー変数以外の係数が1を超えることがある
summary(lm(scale(salary)~scale(cGPA)*(e.dm+bus.dm), dat2))
summary(lm(salary~GPA*(e.dm+bus.dm), dat3))# dat3は全て標準化したもの


これは間違いらしい
# 交互作用項を先につくっておく。
gpa.la <- la.dm*dat$cGPA # ここは"cGPA"にしておく。GPAだと間違い(?)
gpa.e <- e.dm*dat$cGPA
gpa.bus <- bus.dm*dat$cGPA
dat3.1 <- data.frame(dat2, gpa.la, gpa.e, gpa.bus)
dat4 <- data.frame(scale(dat3.1))
res.int2 <- lm(salary~GPA+e.dm+bus.dm+gpa.e+gpa.bus, dat4)
summary(res.int2)