廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

「まちがったモデル」で最尤推定すること

今日の川柳

この文章は『In ALL Likelihood』(p.370) を参考にしました。

In All Likelihood: Statistical Modelling And Inference Using Likelihood

In All Likelihood: Statistical Modelling And Inference Using Likelihood

確率(密度)関数 g(x) を「真の分布」と呼ぶことにする。

 x_1, \ldots, x_n を独立に g(x) に従うサンプルとする。

f_\theta(x) を「モデル」と呼ぶことにする。

ほとんどすべての場合、真の分布は不明なのでモデルは分析者が設定する。

だからほとんどすべての場合、f_\theta(x) g(x) は一致しない。

最尤推定

 \frac{1}{n}\sum_{i=1}^{n} \log f_\theta (x_i)

を最大化することと等しい。

サンプルサイズが十分に大きければ、このことは

 - E_g[ \log f_\theta (X)]

を最小化することと同じである。

ここで E_g は真の分布のもとでの平均を表す。

さらに、このことは

 E_g[ \log g (X)] - E_g[\log f_\theta(X)]

を最小化することに等しい。

なぜならば第一項は \theta の選び方に依存しない未知の定数だからである。

この式は真の分布とモデルの間のカルバック・ライブラ距離である。

つまり最尤推定は真の分布とモデルの間のカルバック・ライブラ距離を最小化することと等しい。

最尤法は「まちがったモデル」を選んでも、選んだ範囲で(カルバック・ライブラ距離の意味で)最適な推定を行う。

例題


真の分布を形状パラメータ4、尺度パラメータ1のガンマ分布とする。

g(x)=\frac{1}{6}x^3\exp(-x)

モデルを正規分布とする。

真の分布とモデルの間のカルバック・ライブラ距離は、

 E_g (-\log6+3\log X-X)-\frac{1}{2}E_g\{\log(2\pi \sigma^2)) + \frac{(X-\mu)^2}{\sigma^2}\}

である。

(途中計算は省略して)

\hat \mu=4
\hat \sigma^2=4

で、真の分布とモデルの間のカルバック・ライブラ距離は最小となる。

R を使って最尤推定のシミュレーションをしてみる。

MLEnorm_sim <- function(i,n){
  x <- rgamma(n,4,1)
  muhat <- mean(x)
  n <- length(x)
  s2hat <-  n*var(x)/(n-1)
  c(muhat,s2hat)
}

library(parallel)
res10 <- t(simplify2array(mclapply(1:10000,MLEnorm_sim,n=10,mc.cores = detectCores())))
res100 <- t(simplify2array(mclapply(1:10000,MLEnorm_sim,n=100,mc.cores = detectCores())))
res1000 <- t(simplify2array(mclapply(1:10000,MLEnorm_sim,n=1000,mc.cores = detectCores())))

boxplot(cbind("n=10"=res10[,1],"n=100"=res100[,1],"n=1000"=res1000[,1]),main=expression(hat(mu)))
abline(h=4,lty=2)

boxplot(cbind("n=10"=res10[,2],"n=100"=res100[,2],"n=1000"=res1000[,2]),main=expression(hat(sigma^2)))
abline(h=4,lty=2)

f:id:abrahamcow:20190117224250p:plain

f:id:abrahamcow:20190117224235p:plain