そして，
またあまり使わない glm(..., family = "gaussian")
の「AIC をどう評価しているのか?」調査のつづき
……
writeLines(as.character(body((glm.fit))), con = "glm.fit.R")
で glm.fit
の内部を
glm.fit.R
出力してみると，
aic.model <- aic(y, n, mu, weights, dev) + 2 * rank
という行がある．
この aic()
なる関数は「上」のほうで，
aic <- family$aic
と定義されているので，
family
を調べないと行けない．
そこで，
今回はなぜかしら
family = gaussian
の場合について調べないといけないので，
writeLines(as.character(body(gaussian)), con = "gaussian.R")
というふうにファイル
gaussian.R
に出力してみると，
aic = function(y, n, mu, wt, dev) {
nobs <- length(y)
nobs * (log(dev/nobs * 2 * pi) + 1) + 2 - sum(log(wt))
}
という定義があったので，
このように AIC が評価されているのがわかる．
もとの glm.fit()
ないで aic()
を呼ぶときの residual deviance はやはり family
中の dev.resids() 関数を使っていて，
これは gaussian の場合，
dev.resids = function(y, mu, wt) wt * ((y - mu)^2)
このように定義されている．