読者です 読者をやめる 読者になる 読者になる

R、統計学、機械学習と競馬予想

学んだことを面白くメモするブログ

緑本3.4.2(あてはめとあてはまりの良さ)のメモ

緑本を読んでいて、3.4.2のポアソン回帰の最大対数尤度のところで、つまづいたのでメモ。

  • ある個体iの平均種子数 λ_i
{
λ_i = exp (β_1+β_2 x_i)
}

であると仮定

  • リンク関数によって
{log λ_i =β_1+β_2 x_i}

  • このモデルの対数尤度は
\displaystyle{logL(β_1 , β_2) = \sum_i log \frac{{λ_i}^y_i exp(-λ_i)}{y_i !}}

  • 本のデータで最大対数尤度を計算すると-235.4ぐらいとわかります。
    ってところが複雑な式だなぁ、実際にどの数字を使って計算すれば計算できるのかなっていうんで、本のデータで計算してみた。
# データの読み込み
>d <- read.csv("http://hosho.ees.hokudai.ac.jp/~kubo/stat/iwanamibook/fig/poisson/data3a.csv")  

# GLMへの当てはめ
>fit<- glm(y~x,data=d,family=poisson)
Coefficients:
(Intercept)            x  
   1.29172      0.07566  

# 最大対数尤度を計算
> logLik(fit)
'log Lik.' -235.3863 (df=2)

# 最大対数尤度を上記の式で計算してみる
>logL <- function(x,y) log({exp(1.29172+0.07566*x)^y}*exp(-(exp(1.29172+0.07566*x)))/factorial(y))
>sum(mapply(logL,d$x,d$y))
[1] -235.3863

ポアソン回帰で切片と傾きの最尤推定値が得られた後に、最大対数尤度を念のため計算してみて、logLik関数の結果と一致したという結果でした。