ぎょーむ日誌 2006-11-16

2006 年 11 月 16 日 (木)

• 0650 起床． 朝飯． コーヒー． 0815 自宅発． 晴． 0830 研究室着．
• どうも R を 2.4.0 を変えて以来， update.packages() に失敗する場合が多い． この version up 時に， なんか CRAN package 依存性に齟齬が生じているような気がするな． install.packages() で手動で補正すればいいだけなんだけど．
• Trendy-ml 経由で今月末のセミナー案内． なんと今回は埼玉の崎尾均さんがみえるのか． ということで， セミナー案内 を更新する．
• 一昨日の井田君の光合成速度「個体差」問題， 魔宮的な nlme 本 を探索するのはあきらめて， 階層ベイズモデル化 → R2WinBUGS で MCMC 計算， という戦略にきりかえる． 「めんどうになったらベイズ」 というのは今後の下うけ生活でありがちな状況になりそうだ．
• 個体差があって， 個体の処理 (状態) が途中でかわって， 同じ個体から異なる時期に三回サンプリングして …… という状況で測定される光合成曲線はどうなるか， というあたりをまずは BUGS 言語 (R みたいなもんです) で書いていく． この問題， けっこうややこしいよ ……
• R2WinBUGS 経由で呼びだす (Wine 実行環境上の) WinBUGS が出力する意味不明ぎみの error message どもを解読しつつばぐとり．
• ようやく走ったと思ったら， 計算にえらく時間かかる． 1000 MCMC step で 80 秒ぐらい? 地衣類繁殖の階層ベイズの数倍， といったところか …… で， よーやく結果でたと思ったら MCMC 計算がちっとも混交してないし． WinBUGS を使ってなおこの状況とは．
• 北大構内走するつもりがシューズを自宅におきわすれてきた． まぬけだ． しかし， いずれにせよ雨． 傘をさして， 郵便局へ． 生態学会松山大会参加費 + 懇親会費 12500 円ふりこむ． 手数料 150 円． 研究室にもどって昼飯．
• 小山さんの R + Tinn-R が呪われぎみのようだったので， どちらも最新版に更新しておはらい．
• 昨日のつづきで平尾君に R の optim() を使った最尤推定法のやりかたの説明． 要点は
  • 対数尤度を計算する関数を作る (これは平尾君が準備・点検ずみだった)
  • help(optim) をよく読む: とくに， control と ... つまり上記の対数尤度関数にデータなどわたすあたり
  • 「最小化ではなく最大化だ」 と optim() に指示する control = list(fnscale = -1) 指定を忘れない (忘れると最尤推定ならぬ最悪尤度推定になるので)
  といったところか．
• そして川合さん研究室雑用 Wiki (PukiWiki) 使いかたの説明．
• 宮田さんに R の plot() の説明．
• …… という旅路のはてに自分の机によーやく戻れた． 井田君の光合成曲線の階層ベイズモデル， これはあきらかに収束してないよね． してるパラメーターもあるけど． 何か工夫が必要だ． いや， それよりも前にまちがい点検を．

• 窓の外は雪．
• で， 夜中まで (計算まち時間に) A 棟 8F あちこちをうろうろしつつ， 「チェインごとにみょーなところに収束してしまう」 問題にとりくむ ……
• いろいろ試行錯誤してみてついにわかった答えはぢつに簡単なもので， (WinBUGS で実装されてる強力かつアヤしげな) Gibbs sampler とて尤度の「劣悪すぎる地形」 にいったんトラップされたら自力脱出は不可能なので， パラメーターの初期値は適切に選んでやりましょう， というだけのことだった．
• もちろん R2WinBUGS を使ってるときでも R の中でパラメーター初期値を生成して， それを WinBUGS にわたすことができる． ただし， 今回はいくつかのパラメーターの初期値が指定できない場合があって， ですね …… s[1], s[2], s[3] といったパラメーターがあったときに BUGS コードの中で

  s[1] <- 0.0;
  s[2] ~ dnorm(0.0, iv.noninformative);
  s[3] ~ dnorm(0.0, iv.noninformative);
  

  てなことをやってると (s[1] だけが定数であとは 無情報事前分布をもつパラメーター)， どの s[] にも初期値が設定できないんだよね (設定するとエラーがでる)．
• で， WinBUGS では「初期値はおまかせ」モードがあるんだけど， おまかせしてしまうとトンでもない初期値をほうりこまれてしまい， さらに Gibbs sampler ではそのトンでもない値周辺から脱出できなくなる， と． いやはや， 考えてみるまでもなくあたりまえのことだ．
• ということで， 解決策としては統計モデルをとりあえず少しばかり変更して， アブないパラメーターにはことごとく無難な初期値を設定してやる．
• ちなみにこういった MCMC 計算， Gibbs sampler ではなく Metropolis-Hastings 法を使っていれば， かくのごとき窮地には陥らないですんだのであろうか? これも M-H 法の「さんぷりんぐ方法」に依存しているはずで， random walk sampler を使っていればやはり 「自力脱出不可」 におちいっていたであろう． 独立サンプラーならだいじょうぶかな?
• あるいは母子里林冠モデルで使ったパラレルテンパリング法 を使ってみるとか． WinBUGS も複数チェイン計算はやってるわけだけど， どれも同じ「温度」なので窮地からの自力脱出には 何の役にも立たないんだよね …… まあ， 複数チェイン計算のおかげで 「WinBUGS おまかせ初期値が面倒の原因」 と特定できたわけだけど．
• とりあえず長めのサンプリングを命じておいて， 撤退． 2240 研究室発． 雪はやんでいる． 2255 帰宅．
• [今日の運動]
  • なんというか ……
• [今日の食卓]
  • 朝 (0720): 食パン．
  • 昼 (1330): 研究室お茶部屋． 食パン．
  • 晩 (2030): 研究室お茶部屋． 食パン …… けっこう最悪な一日だ．