ぎょーむ日誌 2004-07-03

2004 年 07 月 03 日 (土)

• 0800 起床． 朝飯． コーヒー． 昨日の開花確率問題に関して， 観測誤差を考慮したうまい推定を考えつく． ごく単純にある一定の確率で「よみ違える」として， その値を DBH 依存性と同時に推定してやればよいだけのことだ， と気づいた． 1100 自宅発． 晴． 1110 研究室着．
• 釧路妨害な自爆テロ， いよいよお電話攻撃まで． こういうのって大会阻止には何の役にもたってなくて， たんに私個人とかが迷惑なだけなんだが． 大会を憎んで下っぱを憎まずのキモチで企画やりなおしてもらいたいね． 2 分で撃退． 自爆テロリストは二度死ぬ．
• 苫小牧直径成長問題， やはり複数月の成長同時推定は無理っぽい． 推定モデル・計算方式・観測データの量の限界で． 苫小牧ボスからもとりあえず 6 月分だけ， という指示をいただいたんで， そちらをひととーり片づけますかね．
• 次の一手は …… アサダ推定の過程で光合成気象値の定義を少し変えたので， カエデたちを再計算してみますかね．
• そのついでに推定計算高速化に着手してみる． 尤度関数を optim() から評価するたびに， parameter mapping という操作が必要になる． いままではそれを R の中でやってから，

  # p -> mp (mapped paramters)
  matrix.map <- rbind(
  	c(1:length(map.parameters))[map.parameters == 1],
  	parameters
  )
  mp <- map.parameters
  for (r in data.frame(matrix.map))
  	mp[r[1]] <- as.numeric(r[2])
  

  C で書いた拡張モジュールにわたしていた．
• この処理も R ではなく C でやったら速くなるんでは， と憶測してみる． 実験対象樹種はアサダ． 改造前は気象値セット一個で 23-24 秒を必要としていた (growth6 のみの推定; growth6 & growth7 同時推定だと 230-240 秒ぐらい)． さーて， parameter mapping in C 版だと …… 16-17 秒というところか． なかなかいいぞ． 意外と速くなった． いままで 6 時間ぐらいかかってた計算が 4 時間ぐらいで片づくようになったわけで．
• とりあえずイタヤカエデの再計算を命じておく． 昼飯．
• イタヤカエデ成長の推定計算は順調． ということで， 昨日ほうりだしてしまった カエデ開花確率推定関数つくりのつづき．
• 1600 とりあえず完成． R の optim() ではいろんな計算法が指定可能なんだけど， 頑健な (しかし遅い) Nelder-Mead 法はこういう 1 parameter 推定にはむいていない． そこで BFGS 法を使った． これは勾配関数を指定する必要あり．

[1 parameter model]

もっとも簡単な確率論的開花モデル． 二項分布なんで見た目の「ずれ」は大きい． 青線が ``|'' を吐きだしたモデルで， 赤線がそれを推定してみたもの．

• …… といった説明を書いてる最中に， 一個体に花が (仮想的に) 「4 個」ついてたのを思い出した． われながら阿呆というかまぬけというか．
• 1700 改訂版できた． こうやって図にしてみると， これはホントに rbinom(100, 4, prob.flower(DBH)) で生成したパターンであるのに， 思ってたよりばらつきがでかいように見える． 見た目で overdispersion かどうか判定するのは危険だな．

[1 parameter model, 4 flowers]

あてはまりがよくなってる理由は， 単純に標本数が 4 倍になったから． こちらのモデルはコインを同時に 4 枚なげるような試行に対応している (上のものは 1 枚だけ投げる)． かかる 「樹木ごとの個体差なし」 と仮定したモデルであっても開花数はけっこうばらついてるように見える．

• あ! この推定関数って M2 平林さんから 宿題としてだされてたデータ解析問題にも適用できる (女性大学院生たちは入学と同時に無敵性という属性がそなわり， たとえば給料とりどもに宿題をだす権限も生じる)． 横軸が花序内滞在時間でタテ軸を訪花数とみるならば． そうだそうだ， いまは logistic な関数型よりこっちのほうが妥当だよね．
• 1755-1845 北大構内走． 晴．
• お， イタヤカエデ計算おわったか …… 光合成気象値の計算期間を 6-9 月と限定することで， 気温フィルターのカタチがそれっぽくなり， さらにあてはまりもよくなった． つまり気象データからノイズを除去できた， ということか．
• つぎにオオモミジの再計算をやらせておいて 2000 研究室発． 2020 帰宅． 体重 72.4kg． 晩飯．
• オオモミジの推定結果もイタヤカエデと同じ方向に改善されている．
• [今日の運動]
  • 北大構内走 1755-1845． ストレッチング．
  • 腹筋運動 30 × 3 回． 腕立ふせ 3 × 3 回．
• [今日の食卓]
  • 朝 (0930): 米麦 0.7 合． コマツナ・豆腐のスープ．
  • 昼 (1430): 研究室お茶部屋． 食パン．
  • 晩 (2130): 米麦 1.0 合． ニラ・豆腐のスープ． トマト．