その2

さて、よく解析して、p値出して、その通りでしたー、みたいなことを書いてしまいますが、本当に重要なことは、その解析した際に推定された値、変数だったりするわけではないでしょうか。


ということで、より噛み砕こうというのが今回です。
fitted()は、当てはめたモデル式で予測される値です(その元々のデータの数だけ返される)

fitdata <- xtabs(fitted(fitglm) ~spi + state + state, data3)
fitdata
, , pop = osaka

      state
spi    ok out
  saku 19  17
  ume  39  11

, , pop = tokyo

      state
spi    ok out
  saku 45  13
  ume  27  21

(Rは基本的にアルファベット順で並べる傾向にあるので、大阪が先に来ています)
これは、

 

東京 大阪

健康 病気 健全病気
サクラ 45 13 19 17
ウメ 27 21 39 11

という、元々のデータと完全一致しています。

よって、ピアソンカイ2乗統計量も

sum(residuals(fitglm, type = "pearson")^2) ##または
dif <- data3$count - fitted(fitglm)
res.pear <- dif/sqrt(fitted(fitglm))

で、ほぼ0となる。(まぁ、自由度0の飽和モデルなので当たり前か)。
一方、もうひとつのモデルの適合度を計る指標である逸脱度も

sum(residuals(fitglm)^2) ##residual deviance 
2 * sum(data3$count * log(data3$count / fitted(fitglm))) #手計算

?少し、glmの結果と違う気がするが。。?丸め誤差の問題?

上記2の適合度統計量は、ほぼ同じ結果を示すが、Dobson(2002)によれば、逸脱度の方が小さいセル(つまり分割表に小さい数のセルがあると)の影響を受けるので、ピアソンカイ2乗統計量の方がやや良いみたい。

一般化線形モデル入門 原著第2版

一般化線形モデル入門 原著第2版

ちなみに、交互作用を抜いたモデルだと、

fitdata2 <- xtabs(round(fitted(fitglm2)) ~ pop + spi + state, data3)


比較すると、


となる。