ぎょーむ日誌 2010-01-05

2010 年 01 月 05 日 (火)

• 0750 起床． 朝飯． 0830 自宅発． 曇． 0850 研究室着． コーヒー．
• 他のコトはほうりだして岩波本作業のつづき． ベイズなモデルのハナシを 来月はじめ (ひゃーもう一ヶ月後なのか!) にしなければならないので， 4, 5, 6 章はトバして， 新設の MCMC 導入な第 7 章にとりくんでみる．

[とりあえずらくがき & はりつけ]

いきなり文章とかがすらすら出るはずもなく， いろいろとらくがきしたり， そのへんにぶらさげてみたり ……

• いろいろうだうだと検討してみる ……
• 13 時よりデータ解析こんさるなので， 早めに北部生協にいってみたら …… すでに学部生でちょー満員． ここは 1-2 年生の講義棟群に隣接してるので． 売店でアヤしげなぢゃんくおにぎりを購入してすばやく撤退． 研究室で昼飯．
• で， 13 時から矢澤さんのデータ解析こんさるだったのだが …… 事前の危惧とはうらはらに， なンということはない死亡確率推定モデリングのハナシだった． そんなのは BUGS code でこうこう書けばいいんですよ， ですばやく終了．
• 死亡確率 (死亡率，死亡速度) に関する重要な注意点: 階層ベイズモデルを使っても 死亡確率に関する「個体差」は推定不可能． なぜかといえば「死亡」は対象個体に一度きりしか発生しない事象であり， たとえばなかなか死なない個体がいたとしてもそいつの死亡確率が低いのか， それとも単なる「幸運」な個体なのか区別がつかないからだ． つまり死亡にかぎらず一度しか発生しえない事象に関しては個体差は推定できない．
• さーて， 岩波本作業のつづきを …… 新設の第 7 章は「マルコフ連鎖モンテカルロ法と GLM のベイズ化」 といった章題にして， 次のようなハナシを段階をふんで説明していく:
  1. 最尤推定していた統計モデル → MCMC による「パラメーター推定値の分布」
  2. 「パラメーター推定値の分布」を事後分布と考えることにより， ベイズ統計モデルの導入 (ただしこの章では階層ベイズモデルはまだあつかわない)
  3. 事前分布の導入: 事前分布には主観的事前分布・無情報事前分布・階層的事前分布 がある
  4. glm() であつかっていたポアソン回帰 → 無情報事前分布なベイズモデル化 & WinBUGS の導入
  …… といったかんぢで． つまり
  • MCMC なんかやらなくても最尤推定で簡単にあつかえる例題
  • ベイズモデル化なんかやらなくても glm() で十分な問題
  をあえて MCMC 化・ベイズモデル化してみます， といった内容だ．
• で， MCMC の例題づくりだが …… 第 2 章の最尤推定法導入用の例題をつかいまわす． 正規乱数な random walk sampler なんてもちだしてもどうせ読者が混乱するだけだろうから， このように離散化してみて ……

• λ = 3.0 からスタートしててけてけと左右に歩きまーす …… といった離散な random walk sampler で説明することに．

• 離散 random walk な Metropolis-Hastings アルゴリズムを R でてきとーに書く． 本文ではなく脚注で説明しなければならぬ「はしっこの処理」 はこんなかんぢで．

  	p12 <- 0.5 # prob(i -> i.next)
  	p21 <- 0.5 # prob(i.next -> i)
  	if (i == i.min | i == i.max) {
  		p12 <- 1
  	} else  if (i == i.min + 1 & i.next == i.min) {
  		p21 <- 1
  	} else  if (i == i.max - 1 & i.next == i.max) {
  		p21 <- 1
  	}
  	p <- exp(vlogL[i.next] + log(p21) - vlogL[i] - log(p12))
  

  で， サンプリングさせてみると …… やっぱ WinBUGS とかに比べると格段に収束がわるいね． 下の例は 100 step とばしに 1000 回サンプリングした結果なのだが．

• 事後分布の平均・中央値・quantile などはまあ妥当なのだけど， どうにもこうにもヒストグラムがぼこぼこしている印象が (下の図の曲線は尤度から算出された事後分布)． じつはおそろしいことに plot(density(事後分布標本)) すると「なめらかに & もっともらしく」見えるようになってしまうんだよね ……
• ばててきたので撤収． 1930 研究室発． 買いもの． 2000 帰宅． 晩飯．
• [今日の運動]
  • るー
• [今日の食卓]
  • 朝 (0810): 米麦 0.6 合． 卵焼き．
  • 昼 (1240): 研究室． 北部生協で買ったコロッケ・ゴボウきんぴらいりおにぎり …… ぢゃんく感つよし．
  • 晩 (2030): 米麦 0.8 合． コマツナ・ショウガ・豚肉の炒めもの． ヨーグルト． リンゴ．