ぎょーむ日誌 2007-09-12
2007 年 09 月 12 日 (水)
-
0850 起床．
うーむ，
あいかわらず，
やられているな
……
朝飯．
コーヒー．
1020 自宅発．
曇．
1035 研究室着．
-
13 時間かかるとおもってた ECOAS 主要プランクトン
20 種データの階層ベイズモデルの計算，
実際のところは 12 時間で終了してた．
収束はぜんぜんよくない
……
やはり 20x20 の分散共分散行列の推定はこのゐんばぐすでは無理なのか，
と考えたんだけど，
じつは当方も一ヶ所モデリングをしくじってることに気づいた．
これって二重に random effects が入ってる部分がある!
-
……
と思ったけど，
モデルをよくよく整理してみると，
log(密度) = (池の効果) + (種の効果) + (池・種の効果)
とモデリングしてよい，
とわかった．
あ，
でも他の環境要因な線形予測子とは区別して書かんといかんわ．
-
ゐんばぐすの意味不明えらー対策．
初期設定ファイル inits*.txt の読みこみで
expected c operator
うんぬんとでてきて計算続行してたら，
これは初期化の部分で一次元の vector をラベルつきで渡している，
ということ．
この場合は初期値 vector に重畳にも as.vector(...)
つけると解決する．
-
昨晩帰宅時に気にしてた「よけーな逆行列計算」
あってもなくても計算時間にはほとんど影響しない
(わずか 300 MCMC step に 880 秒ぐらい)，
とわかった．
20x20 の逆行列計算ぐらいは数百分の 1 秒ぐらいでできる，
ということなのか．
-
ちなみに分散共分散行列はヤメて，
動物プランクトン種間に相関がなく random effects がある，
というモデルに変更すると 300 MCMC step の計算はわずか
24 秒
で終了しちゃうんだよね．
いやはや．
-
ハラがたつので，
この無相関モデルで少し計算してみるか．
-
12000 MCMC step の計算やらせても
640 秒
で終了する．
それはいいんだけど，
収束がむちゃくちゃだ．
うーむ
……
-
それから夜まで階層ベイズモデリングの試行錯誤
ひたすらひたすら続く．
モデルはぢりぢりと改良されてて収束も良くなってはいるんだけど，
数千 MCMC step ぐらいではなかなかキレイにそろってくれない
……
-
ああ，
なんか，
私にとっては
こういう「水中」「組成」モノはけっこう鬼門だなぁ．
推定がしんどい理由はおそらく，
サンプルごとの random effects が大きい，
にもかかわらず (ホントの) 測定時誤差の統計モデリング
をどうすべきなのかよくわからない，
これらの区別が難しくて MCMC 計算がなかなか定常状態に到達しない，
といったあたりか?
-
そして東北地方からはいろいろとメイルが
……
-
へろへろと撤退．
2000 研究室発．
ばてたので外で晩飯．
2130 帰宅．
-
R
こんさる．
apply()
系関数と文字列操作は R 中級者にむけての昇段試験だねえ．
> # たとえばこういう data.frame があったとする
> d
x y
1 1.6 1.1
2 1.1 0.6
3 0.4 -0.7
4 -0.3 0.8
5 0.2 -0.1
> # apply() と sprintf() を組みあわせて新しい label 列が作れる
> apply(d, 1, function(r) sprintf("xy=%g:%g", r["x"], r["y"]))
[1] "xy=1.6:1.1" "xy=1.1:0.6" "xy=0.4:-0.7" "xy=-0.3:0.8" "xy=0.2:-0.1"
> # あるいは paste() を使ってもよい
> d$label <- apply(d, 1, function(r) paste("xy", r["x"], r["y"]))
> d
x y label
1 1.6 1.1 xy 1.6 1.1
2 1.1 0.6 xy 1.1 0.6
3 0.4 -0.7 xy 0.4 -0.7
4 -0.3 0.8 xy -0.3 0.8
5 0.2 -0.1 xy 0.2 -0.1
だめだ，
夜は仕事がススまん
……
[今日の運動]
[今日の食卓]
- 朝 (0940):
米麦 0.6 合．
ワカメ・豆腐の味噌汁．
ブナシメジ・卵の炒めもの．
- 昼 (1530):
研究室．
米麦 0.6 合．
クノールのカップカレースープ．
- 晩 (2030):
グロビュールのスープカレー．
ここの自家製ソーセイジ
(豚肉あるいはエビ)
がたいへんおいしかった
……