ぎょーむ日誌 2005-08-24

2005 年 08 月 24 日 (水)

• 0710 起床． 朝飯． コーヒー． 0905 自宅発． 晴． 0920 研究室着．
• 昨晩やらせといた計算だが …… やはり「ササ区」と呼ばれる観測データへの 葉群あてはめは難しい， とわかった． 昨日追加した「途中から計算」機能をつかって， 25 万 → 75 万 step (4 時間) を計算させてみた． 対数尤度はほとんど向上しなくなる． しかし葉数 (葉面積) はぢりぢりと変化しつづけているんだよね．

• 欠側値補間に問題ありかも， ということでそのあたりの関数かきかえて実験してみる． これまでは「上のハコ」の値をそのまま引きついでいたんだけど， 今度は上下のハコをみてその幾何平均とるような方式を試みてみる．
• これがまたハコがタテふたつ欠落してたり， とわけのわからん状態になってたり …… 再帰的に等比補間するしくみをつくる．
• やってるうちにエラー表示にこるようになったり．

# Fix_light_NA at (90, 240, 630) ... # changed (824, NA -> 622, 470)
# Fix_light_NA at (90, 300, 630) ... # changed (717, NA -> 274, 105)
# Fix_light_NA at (300, 120, 630) ... # changed (93, NA -> 84, 77)
# Fix_light_NA at (300, 180, 630) ... # Fix_light_NA at (300, 180, 690) ... 
# WARNNING: light_NA at 009-003-006 ...# changed (78, NA -> 70, 58)
# changed (70, NA -> 63, 58)
# Fix_light_NA at (300, 480, 630) ... # changed (141, NA -> 68, 33)
# Fix_light_NA at (300, 540, 630) ... # changed (45, NA -> 42, 41)
# Fix_light_NA at (300, 180, 690) ... 
# WARNNING: light_NA at 009-003-006 ...# changed (78, NA -> 70, 58)
# Fix_light_NA at (300, 240, 690) ... # Fix_light_NA at (300, 240, 750) ... 
# WARNNING: light_NA at 009-004-007 ...# changed (1000, NA -> 355, 45)
# changed (355, NA -> 126, 45)
... (後略) ...

• しかしこの エラー表示 みてわかるのは， 「光の落差」がでかいところが林冠最上部にもあったりして， というですね．
• そもそもこの「ササ区」ってのは， 観測区画の最上部にアヤしげな葉群のかたまりがあったりするんだよね． 下の図は新しい方式で等比補間ずみの初期状態なんだけど， それほど改善されてないよーな ……

• とりあえず A801 室の Dell 機で， 新しい初期状態で計算やらせてみる．
• まだ計算途中なんだけど …… 対数尤度や枚数のあがりかたが 前回 と同じもしくはそれよりもゆっくりしている． どうやら欠側値等比補間のやりかたがまずかったわけではないようだ．
• 昼飯． 昼飯くってもやはり対数尤度の上昇はゆっくりしたまま． 初期状態の葉数 (葉面積) をかえてみるか． ということで再計算．
• この程度では， MCMC 計算の「頑健さ」はびくともしない， とわかってしまった …… つまり結論は， 「ササ区」と呼ばれる観測データにあてはまる葉群分布を推定することは， すごく難しいものである (「除去区」データに比べて，という意味)， と．
• 生態学会さーばー雑用． いつも =NEW= にあふれてる 学会ペイジ どものひっこし準備あれこれ． そしてまたメイルアカウント作り …… のため， また R で乱数ぱすわーど作り．

  moji <-c(letters, LETTERS, 1:9) # 0 は使わないことに
  paste(sample(moji, 12, replace = TRUE), collapse = "")
  

• このレンタルサーヴァーの backup がだんだん心配になってきた．
• 大学院雑用 (さいわいにも私には無関係) メイルくる． 「今年度は教育経費は専攻で管理」 ということで， 使ったら専攻科事務に報告の必要あり． で， その経費ってのは修士課程院生は一人 8.9 万円， 博士課程は一人 12.7 万円， だそーで ……
• 母子里計算， 「ササ区」の計算ひとまず終了． つぎに「除去区」を (いまの欠側値等比補間のもとで) また計算させてみる． こっちは 20 倍以上の速さで対数尤度が改善されていくのになぁ．
• アタマが母子里 MCMC ばてぎみなので書類雑用． 「平成 17 年度 教育・管理運営・社会貢献業績調査票」 とゆーやつなんだが …… 何もぎょーせきがない． べつにいいけど． なになに 「III. 社会貢献」の「7. 専門分野に関連して作成した作品」 …… って ぎょーむ日誌 のこと?
• 雪野さんとアヤしそうでアヤしくないお金のやりとり． 6 月の阿寒実習 の出張費 50200 円がいまごろ大学からふりこまれたんだけど， これは雪野さんに全額おわたしする． えー， つまり「たてかえて」もらってることになってるので．
• A801 室の Dell 機に MCMC 計算駆動させつつ， また 伊庭さん の 「ベイズ統計と統計物理」 (表紙) をぱらぱらとながめる． これは MCMC 計算やるヒト必読の一冊です． たとえば p.52 の

  シミュレーテッド・アニーリング法はあくまで最適化法．
  分布からのサンプリングや期待値計算に使える保証は無い

  といった記述を読むと， 「ササ区」で尤度を押しあげる策なんぞを検討していた アタマが冷静になる． そうだよな． ここで「あてはまり」良い葉群配置だけを探索しても意味ないわけで．
• 嗚呼， しかし， 「除去区」だとわりとすぐに定常状態に到達して， いーぐあいに楽しげにふらふらと尤度がゆれるのになぁ ……

  ...
  # 138386: Increment (-7141.0 vs -7154.8, F -4374.7) -> REJECT, n_leaf = 13823
  ...
  # 145344: Increment (-7081.2 vs -7081.7, F -4374.7) -> ACCEPT, n_leaf = 13870
  ...
  # 155471: Decrement (-7086.7 vs -7086.8, F -4374.7) -> ACCEPT, n_leaf = 13927
  ...
  # 171110: Increment (-7030.6 vs -7030.6, F -4374.7) -> ACCEPT, n_leaf = 13886
  ...
  

  いや， 葉数増減はふらついてるが， 尤度はまだぢりぢりと改善されてんのか?
• なーんか市内からサイレンの音が止まらないわ， 下の階のひとたちはどやどやとあがってくるわ， なんなんだと思ってたら …… ここから真南 1000m ぐらいのところから派手に煙がでてた． かなりでかい火事． 時刻は 1645 ごろ． (ケムリで見えてない) 旧道庁のすぐ北あたり， アスティ 45 の「裏」あたりかな． 下の写真は A 棟 8F から撮影．

• いやー， ここから見ると あたかも JR 札幌駅南面ぜんたいがケムリにまかれてるよーな …… 火事でこわいのは火だけぢゃない， とあらためて思い知った次第．
• しばらくするとケムリはおさまったようだ． 車輌通行規制でクルマがたまってる． ヘリコプターとびまくり． ネット上の報道は 読売新聞 が一番速かったか (1718)． 「建設中の駐車場ビルで爆発」 だそーで． 灯油タンク・ガスボンベが置かれてたとか．
• そういや， 先月はじめ にアスティのクリニックで健康診断うけたときに， 裏手でなんだかビルの工事やってたな． あれが爆発したのか．
• だめ押しの解析として， 「ササ区ではいったい どの部分のあてはまりが悪いのか?」 を調べてみた …… やはり観測区画最上部の孤立葉群付近だったか．

• 上の図の左は「除去区」で右が「ササ区」． 球の大きさが「あてはまりの悪さ」をあらわしている．
  • 赤球: 「観測値あり」地点
  • 紫球: 「欠側 (等比補間値で代用)」地点
  「除去区」はあてはまり悪いところない． いっぽうで 「ササ区」(右) の孤立葉群まわりに 大きな球が集中 してるのがわかる つまり「モデルがすごくダメ」なのはここら付近だけ． そして赤球ばかりなんで欠側値補間と関係なさそうだ．

[ずれの方向]

このあたりの図の葉群配置はどれも MCMC 計算で得られた Gibbs 分布標本 (25 万 step)

• 「ササ区」のずれの方向の上の図のとおり．
  • 黄球: 観測値に比べてモデルが「明るすぎ」 (周囲に葉っぱ少なすぎ)
  • 黒球: 観測値に比べてモデルが「暗すぎ」 (周囲に葉っぱ多すぎ)
  なぜ尤度が改善されないかよくわかる． つまり， 黄球周辺で葉っぱ足りないというんで増やしてやると， その下の黒球周辺はますます暗くなる． 黒球周辺を明るくするため， 上の葉っぱを削除すると今度は黄球周辺が明るくなりすぎる． いやはや．

• とゆーわけで， MCMC 計算やってもかくのごとく明るさの過大推定・過小推定が解消されづ， 尤度は改善されない， と．
• なーんか， 観測データうたがう心境 ……
• 観測データに関してはテのうちようがないんで撤退． 1850 研究室発． 1900 帰宅． NHK ラジオの地方ニュースでくだんの火事報道やってるな． 北 4 西 6 に工事現場で 120 名がその場にいて， 溶接作業やってたとか． 2050 自宅発北大構内走． 2145 帰宅． 体重 72.8kg． 晩飯．
• [今日の運動]
  • 北大構内走 55 分間．
• [今日の食卓]
  • 朝 (0810): 米麦 0.7 合． シイタケ炊きこみ飯． タマネギ・ブナシメジ・豆腐の味噌汁． ミニトマト．
  • 昼 (1250): 研究室お茶部屋． 米麦 0.6 合． シイタケ炊きこみ飯． タマネギ・ブナシメジ・豆腐の味噌汁．
  • 晩 (2220): 米麦 0.7 合． モヤシ・ピーマン・ニラ・ショウガ・豚レバの炒めもの． タマネギ・ブナシメジ・豆腐の味噌汁．