ぎょーむ日誌 2008-03-(01-10)

2008 年 03 月 01 日 (土)

• 0730 起床． 本日は大雪山のふもとで 進化植物学研究会． 朝飯． コーヒー． 0835 自宅発． 曇． 0845 北 11 西 4 サンンクス前の集合場所． 0900 すぎ， 斎藤君の車にのせてもらって出発．
• 道央道にのって旭川方面へ． 旭川にちかづくにつれ天候が悪化してくる． 途中で 1020 ごろ砂川 SA で休憩． 1130 ごろ旭川鷹栖 IC からおりる． 市街地をうかいして富良野方面への国道にのる．

[1144 国道 276 号]

旭川から旭川空港・美瑛・富良野方面にむかう道． すでに風雪状態． 視界もわるい．

[1307 国立大雪青少年交流の家]

白金温泉のちょっと奥にある巨大な施設． よーやくついた． この施設は 29 才以下が過半数だと宿泊料は無料 (ただし晩飯・朝飯・その他で 1500 円ぐらい) で使わせてもらえる なかなかすばらしいところではあるんだけど ……

[研修室]

ここで今日・明日と五題の生態学発表． ひとりあたり一時間以上時間をつかうゆったりしたスケジュールなんで， いいかんじで勉強ができた． 研究発表会の前に， 施設のヒトによるオリエンテイション． やはり教育施設なので明朝 0715 よりラジオ体操とか， 朝はきちんと掃除しましょう， といったような ……

[1315 より研究発表会]

本日は荒木さんのスズラン個体群動態・繁殖生態， 福井さんのコウモリの森林ギャップ利用・河川上空利用， 最後に奥田君の (例の) 太平洋岸 5 地域 × 5 海岸 × 5 岩礁 の群集データにみられる群集構造モデリング， といったよい内容でした．

[1700 ゆーすぴあタイム]

この施設の毎日の行事のひとつ． 20 分間ほど． 本日の施設利用する全グループ， すなわちわれわれと， 旭川のとある小学校のミニバスケチームの自己紹介みたいな． この体育館は本日改装後初使用とのことで， 温泉による床暖房になってる． 暖かすぎるのか天井で越冬してるテントウムシがやたらと落ちてきたんだけど．

[1750 晩飯]

この施設はセルフサーヴィスというかビュッフェ形式． 注目していただきたいのは， 内容ではなく大学院生たちの皿数である． 小学生たちでさえお行儀よく一皿に料理もりつけてたべているというのに， このヒトたちは一人あたま二皿の料理を …… さらに猛者は「おかわり」をこころみていた． じつはこの「食い放題」が場所選択の最重要要因だったのである ……

• 晩飯後 & 温泉内でも奥田君たちにベントス話をうかがう． ああいう 岩礁性潮間帯 の固着生物なんかも「集中分布する」 機構と利点がある， といったハナシなんかもあって …… だったらランダム群集からずれるのも当りまえだよねえ， といったような． ランダムでなければ Gibbs 群集， ということになるのだが． はてさて， どう考えればよいものやら ……
• その後は談話室で懇親会． 総参加者 24 名の自己紹介のあと， また研究に関係のあるハナシないハナシがつづく． 福井さんからコウモリ echolocation の詳細などおもしろいハナシを いろいろと教えていただく． 周波数を変えて定位するとかドップラーシフトの利用とか ……
• 2200 で利用時限になったので， おひらき． 当施設は早寝早起きの原則が全員に適用されるのである． ただし宿泊棟内でその後も歓談つづく． 私はぎょーむ日誌など書いて， 早めに寝てしまおう …… と思ったんだけど， その前によびだしをくらってしまった． けっきょく寝たのは 2340 ごろ． なにしろ明朝 0630 起床だからね． むろん， この点をあまり重視してないヒトたちは ……
• 生態学会福岡大会初日 (3/14 金) まであと 13 日．
• [今日の運動]
  • 寝る前にちょっとストレッチ．
• [今日の食卓]
  • 朝 (0750): ヨーグルト． リンゴ．
  • 昼 (1300): くる途中， 昼飯の時間がなかったので， セミナー前とかにこっそりセイコーマートにぎり飯とか食ってしまった ……
  • 晩 (1750): 青少年交流の家晩飯．

2008 年 03 月 02 日 (日)

• 0630 起床． 国立大雪青少年交流の家， 進化植物学研究会 の第二日目． 紺野さんと宿泊棟ロビーでベイズ雑談．

[0715 さわやかタイム]

当施設の朝の行事． このあと小学生たちといっしょにラジオ体操． ただし， こちらのメンバーのうち何人かは部屋でまだ寝ていたのだが ……

• 食堂で朝飯． コーヒー． 荒木さんとスズランのモデリング雑談 …… ともうしますか， ベイズ布教． Genet 差と場所差を分離推定するのはむずかしそうだねえ ……
• 宿泊棟のそうじ． その後， またセミナーの研修室に移動．

[0900 本日のセッション開始]

今日は飯島君のトドマツ・エゾマツの実生定着要因の解析 (Meadow による lmer() 実演つき)， そして渡邊さんによる北大構内から発掘された縄文・続縄文・擦文・アイヌ分化 の遺跡からみつかった (おもに) 炭化木材の走査電顕写真から特定された 昔のヒトたちの樹木利用について． トネリコ属 (ヤチダモとか) なんかが竪穴住居の柱でアジサイ属 (ノリウツギ) とかが木製道具になるそうで．

[1150 ごろ終了]

いろいろと勉強になる研究発表のあと， 来年とかどうしようといった議論が 30 分ほどあり， 今回の進化植物学研究会は午前中に終了． 外はまた悪天．

• 撤収準備終了． また斎藤君の車にのせてもらって 1220 雪の中を出発． 美瑛-旭川境界ふきんの丘で昼飯． 旭川から道央道にのって札幌をめざす．

[1530 ごろ岩見沢 SA のトイレ]

北海道における女性の社会的地位の高さを示す注意書きを発見した．

• 札幌付近は晴． 自宅ちかくまで送ってもらって (どうもありがとうございました)， 1640 帰宅．
• ぐだぐだと休む． 夜は研究会とは関係ない道外からのお客さんのご接待．
• 生態学会福岡大会初日 (3/14 金) まであと 12 日．
• [今日の運動]
  • 小学生たちと朝のラジオ体操． ストレッチ．
  • 腹筋運動 30 × 3 回． 腕立ふせ 10 × 3 回． スクワット 100 回．
• [今日の食卓]
  • 朝 (0730): 青少年交流の家の朝飯．
  • 昼 (1310):
  

  [コロッケカレー 750 円]

  ぜるぶの丘なる美瑛町はずれの国道ぞいの丘で． ヤーコン サラダつき．
  
  
  • 晩 (1920): ご接待晩飯． 大学周辺も料理のおいしい居酒屋がふえてきた ……

2008 年 03 月 03 日 (月)

• 0820 起床． よく眠れなかった …… すごくねむい． 朝飯． コーヒー． 0930 自宅発． 曇． 0945 研究室着．
• さすがに減ってきたデータ解析こんさる宿題， さっさとかたづけなくては …… と思っているんだけど， なかなか着手できない．
• 井田君とちょっとタチギボウシデータのモデリング雑談．
• かなり前に 「nested 事後分布図の作図プログラムください」 と依頼されてたので， その準備． 作図プログラムからよけーなところをそぎおとして， 抽象化というか一般化のどあいを少しだけ強めてみる …… という作業に時間とられる． とりあえず， nestedplot20080303.tgz (ゐんどーづなヒトなどはこの中の win ディレクトリ内ファイル を使ったほうがよいかも) というのを作って簡単な説明かき．

[簡単な説明]
とりあえず

・展開したフォルダ中の post.mcmc.RData をダブルクリックして R を起動し
  てください (これは R の起動と作業ディレクトリの変更を同時にやってい
  るだけで，post.mcmc.RData じたいはあとから作図プログラム中からあらた
  めて load() しています)
・とりあえず，R のコマンドラインで source("plot.R")

とすれば図を描くはずです．

ファイルについて解説しますと

・post.mcmc.RData: 例題データ

・plot.R: plot.func.R 内の関数をよびだして post.mcmc に格納されている
  データで作図する R スクリプト

・plot.func.R: 二段に nest している事後分布図を描く関数の定義，このと
  きに nestplot.R 内で定義されている関数を利用している

・nestplot.R: 汎用的な nested 事後分布図を描く関数群の定義，何段 nest
  していても作図できるように定義されている

• bugs オブジェクト (library(R2WinBUGS)) を mcmc オブジェクト (library(coda)) に変換するのは このあたり 参照．
• というふうにまた宿題が一件落着したら， とっくに昼すぎ． 北大生協中央購買まで昼飯調達の旅にでる． 歩いてると， さきほど井田君にかなりでたらめなコトを教えたような不安感がつのってきた． 購買に寄ってから 7F 井田君のところにいって詳細を確認． 結果としては， 当方の思いすごしで問題なさそうだった． いかん． やはり今日はばてぎみだ． 昼飯．
• また別のデータ解析こんさるメイルかき．
• 7F の Annika さんに呼ばれて， 彼女の assistant student の計算機不調の修理． Knoppix を使ってハードウェアや HDD の file system がこわれていないことを確認． 診断: ゐんどーづが自殺的な自己書きかえをしたので OS が起動しなくなりました． 処方: りかばりーディスクみたいなものをさがしてください …… ゐんどーづ一般ユーザーではこれは無理だろうなあ， と危惧しつつ．
• 私が共著者になってる苫小牧でんどろ論文， CJFR に却下くらったのだが …… うーむ， reviewer 二人とも GLM のたぐいの基礎をまったく理解できていない， という状況はなんともはや ……
• ばててきたので早めに撤退． 1750 研究室発． こういう日は散髪でも …… とクラーク会館地下にいくと， 閉店時間がちかいのに先行客が二人も． 気のヨワい私はすごすごと退散． 買いものして 1900 帰宅． 体重 70.6 kg． 晩飯の準備． 晩飯．
• 生態学会福岡大会初日 (3/14 金) まであと 11 日． 最低限のやることリスト:
  • 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  何も着手できていない． 明日はアリポスターをススめたい． 今日いろいろやったからデータ解析こんさる宿題は現時点では一掃されたんで． 嗚呼， データ解析こんさるではない宿題は一件のこっていたなぁ ……
• [今日の運動]
  • うんどう休養日
• [今日の食卓]
  • 朝 (0820): ヨーグルト． リンゴ．
  • 昼 (1400): 研究室． 北大生協サンドイッチ．
  • 晩 (2030): 米麦 0.8 合． 麻婆豆腐． ダイコンサラダ．

2008 年 03 月 04 日 (火)

• 0830 起床． 朝飯． コーヒー． 苫小牧樹木直径成長問題， もし Bayes 的にやりなおすとしたらどうなるか …… という問題を検討してしまう． 定式化は簡単． 問題はゐんばぐすごときでは計算にすごく時間がかかるかもしれない， という点だ． これは実験してみないとわからない． もしダメだった場合は自作 MCMC プログラミングが必要になる． Metropolis-Hastings 法は実装は簡単だけど， 良い計算やらせるまでがたいへん． Adaptive rejection sampler (ARS) は実装がたいへんそう． となると slice sampler あたりがお手軽 …… ではあるんだけど， 非負の値しかとらない確率分布からのサンプリングはたいへんそう， と気づいてしまった． うーむ．
• あらら， もうこんな時間． 1000 自宅発． 晴． 1015 研究室着．
• 北大ネット， 大学-外界のあたりで切れてる …… 30 分ほどで復旧． あてにならぬ 「運転状況」 ペイジではいつものごとく「なかったこと」 あつかいのようで．
• 来年度も地環研のネットワーク委員会担当せよ， という指示メイル． 委員というか地環研さーばー管理者ですな．
• いくつかメイルかき．
• アリポスター作業， よーやく再開．
• 昼飯． アリポスター作業， つづく．
• 院生たちが 「spam メイルがきたんだけど何か問題でも?」 とのことで …… うーむ， この メンバーペイジ はそれなりの対策をしているので， そう簡単には収集ロボットに仕事をさせてないはずなんだが． 状況をうかがうと， 何でも一通きただけのようであり， (同じ条件におかれている) 他の院生はあいかわらず spam free な世界にいるようなので， ともかくもうちょい様子をみることに．
• そろそろ髮がぼさぼさになってきたので散髪するべく 1730 研究室発． 1740 クラーク会館地下というか 1F の理髪店． 今日はぜんぜん待たずにすんだ． 1835 終了． 2625 円． 買いものなどして 1930 帰宅． 晩飯の準備． ちょっとうんどう． 晩飯．
• 夜は作業がススまない …… いくつかメイルかき， ぐらい．
• 生態学会福岡大会初日 (3/14 金) まであと 10 日． 最低限のやることリスト:
  • 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  本日のポスター進捗はこんなかんぢで．

• 発表準備作業おそい久保にしては けっこうススんでるんでは， と思われるかもしれないがこれは 昨年 7 月末 とかにじたばた作ってた投影資料の図をうつしかえてるだけ， なんだよね …… 例によって \documentclass[portrait,a0b,final,colordvi]{a0poster} なポスターなんで電気鼠いらずのレイアウト作業なんだが …… うーむ， やってるときはあまり気にならないんだけど (いや LaTeX 使ったレイアウトは難しいんですが)， やっぱり時間を費してしまうのである．
• なんだか些細なコトなんだけど， LaTeX ファイルとかなおす → make → evince で reload して確認 → もとにもどる， というるーちんの中で evince での確認に時間かかってるような． この viewer は何でもむりやり図化してしまう強力なソフトウェアなんだけど， 動作がトロいんだよねえ …… 上述のごういんなる A0 ポスター生成法では gv だの xdvi だのは何の役にも立たないし．
• 今さらながらふと思いついたんだけど， \usepackage{pstricks,pst-grad} 使わぬ地味なポスターにしておけばよかったのかな? そもそも， この package をぜんぜんうまく使いこなせていないわけだし．
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 10 × 3 回． スクワット 100 回．
• [今日の食卓]
  • 朝 (0820): ヨーグルト． リンゴ． ダイコンのコンソメスープ．
  • 昼 (1320): 研究室お茶部屋． 米麦 0.6 合． ジャガイモ・ニンジン・タマネギ・豚肉の煮物．
  • 晩 (2130): 米麦 0.8 合． ニラ卵． キャベツ・ニンジン・ネギ・ハムのサラダ． ダイコンのコンソメスープ．

2008 年 03 月 05 日 (水)

• 0730 起床． なんだか本日もねむい． もっと早寝しなくては． 朝飯． コーヒー． 0910 自宅発． 晴． 0925 研究室着．
• 夕方までぽすたー作りのつづき． LaTeX a0poster class をつかった作業でコマかい点がいくつかわかった．
  • 用紙指定は a0b ではなく a0 にしたほうがよさそう (a0b はタテ:ヨコ比が A 用紙にあってない)
  • preview 問題: Ghostscript viewer (gv) はダメ (そもそも何も表示しない) で evince は動作が遅い (レンダリングが丁寧すぎる) という難点は ggv を使うと多少マシになるとわかった
  などなど．
• とちゅう， いろいろとわりこみが． なんだか良くわからない群集生態学作文の査読やれとかいう 問答無用なかんじのメイルとか (当方がひきうけるのは当然といった内容のすごく奇妙なしろものだった)， いまも修論作文の書きなおしで困っている院生からの相談とか ……
• ともあれ夕方になってひととーりできたので， やや苦闘して PDF 変換 (postscript 修正が必要だった)， ファイルアップロード， 共著者たちへのメイル連絡， でこの件はとりあえずいったん終了． 時刻はすでに 1700 すぎ．

• で， メイルかきとかあれこれやってるうちにばててきたので， 1855 研究室発． 1910 帰宅． 晩飯の準備． ちょっとうんどう． 晩飯．
• 生態学会福岡大会初日 (3/14 金) まであと 9 日． 最低限のやることリスト:
  • (ひとまず終了) 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  今日はとりあえずポスターができた． さてさて ……
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 10 × 3 回． スクワット 100 回．
• [今日の食卓]
  • 朝 (0750): ヨーグルト． リンゴ． キャベツ・ニンジン・ネギ・ハムのサラダ．
  • 昼 (1320): 研究室お茶部屋． 米麦 0.7 合． キャベツ・ニンジン・ネギ・ハムのサラダ．
  • 晩 (2010): 米麦 0.7 合． キャベツ・ニンジン・ホウレンソウのシチュー． 豆腐・卵の炒めもの．

2008 年 03 月 06 日 (木)

• 0750 起床． 朝飯． コーヒー． 0910 自宅発． 雪． 0925 研究室着．
• なんだかいろいろなメイル書きにまた時間とられる．
• Agricultural and Forest Meteorology に提出した母子里 Gibbs 林冠モデル 論文原稿の 2nd revision が 3ヵ月ほど塩づけ状態なので， といあわせのメイルを書いてみる．
• 生態学会福岡大会初日 (3/14 金) まであと 8 日． 最低限のやることリスト:
  • (ひとまず終了) 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  とりあえず， 占部さん計算したうけをやることに． 昨年 9 月上旬 に試行錯誤してたアレをデータを少し変えて再計算せよ， というもの．
• またゑくせるファイルとか送られてきてイヤだなぁ， と思っていたのだが …… よく見ると全 86 湖沼のうち「泉」になってる 6 湖沼をのぞいて 再計算せよ， というだけのことだった． 今回のぞかれるのかこんなサイトで． 泉でもサカナがけっこういるぢゃん．

> d[c(1, 9, 24, 72, 77, 84),]
   site site.name year date latitude altitude     area max.depth w.temp   pH      DIN      DTP
1     1  御在所沼 2005    2   39.949      890 0.010000       3.1  20.10 2.90  9.79029 0.153215
9     9      菅沼 2005   11   40.593      465 0.121170       1.0  13.50 8.90  3.41061 0.276334
24   24    太郎湖 2005   39   43.433      419 0.030000       8.7  20.30 8.20  0.50861 0.426813
72   72    鶴間池 2006   10   39.056      815 0.017700       3.8   8.53 6.62 13.75769 0.608117
77   77      空沼 2006   31   43.625     1460 0.012800       0.5   8.64 6.84  7.90282 0.689601
84   84  姿見の池 2006   36   43.657     1630 0.002902       4.5  18.46 3.89  0.94985 0.098842

   sestonC sestonN  sestonP  totalN  totalP   Chl.a bacteriaC  biomassZP fish latitude.1 longitude
1  46.2281  6.2841 0.080322 16.0744 0.23354 0.93010        NA 3.1426e-01    0     39.949    140.94
9  25.5378  4.4798 0.265815  7.8904 0.54215 8.34960    16.488 1.1081e+00    1     40.593    140.95
24 53.6178  8.3922 0.420943  8.9008 0.84776 6.17700    43.582 1.1014e+00    1     43.433    144.14
72  9.0628  1.6433 0.045590 15.4009 0.65371 0.13579        NA 4.6903e-06    1     39.060    140.06
77 19.5992  2.9069 0.160286 10.8098 0.84989 0.61060        NA 0.0000e+00    0     43.630    142.91
84 33.0476  4.7156 0.119099  5.6654 0.21794 1.84777        NA 4.3665e+00    0     43.661    142.83

   rain.fall   c.area cathchment.area
1     845.39   0.6387          0.6321
9     793.94   2.2995          2.2771
24   1156.29 162.8525        147.0418
72   1830.81   2.8430          2.8194
77    668.75   0.6115          0.5950
84    668.75   0.8613          0.8587

• 計算時間は ThinkPad X31 でも 300 秒ぐらい …… しかしながら， 収束しない． 逆 Wishart 分布に与える対角要素を 2 から 10 に変更してみる．
• やはり収束しない． ひょっとしたら勘ちがいしてるのかなと不安になったので 0.1 にするとこれはぜんぜんダメ．
• しょうがないので 20 にしてみる …… こういういきあたりばったりではダメのようで． 基本的に， 「魚の密度」 の部分が発散ぎみだ． これはわりととうぜんというか， 昨年 9 月の計算結果をみてもばらつきが大きいわけで．

• ということで， v.diag <- rep(2.0, N.level) にもどして， 計算 step 数を 12000 から 22000 にしてみる …… しかし， これもダメ． 魚密度に関してもう少しモデリングの工夫が必要というかんじだ．
• うまくいかないまま昼飯．
• どうも魚の密度だけが問題なのではなく， 水温モデリングも改善すべきだ， という気がしてきた．
• 試行錯誤の計算まち時間中に， 年度末的なめんどー事務書類にてきとーに対処してみる． はてさて， どうなることやら．
• 超事前分布をより無情報化することで， よーやく結果が得られた． 12000 MCMC step 300 秒 (事後分布表, BUGS code)． まあ， 動物プランクトン密度まわりの相関が見えにくくなったのはしょうがない， というかんぢだ …… 水温まわりは改善されたかな?

• これであとは報告を書いて， と思っていたら何やらまた PC 修理だの何だのにかりだされてしまった …… うう．
• 1900 研究室発． 1920 帰宅． 体重 70.2 kg． 晩飯の準備． ちょっとうんどう． 晩飯．
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 10 × 3 回． スクワット 100 回．
• [今日の食卓]
  • 朝 (0810): キャベツ・ニンジン・ホウレンソウのシチュー． ヨーグルト． リンゴ．
  • 昼 (1310): 研究室お茶部屋． 米麦 0.7 合． キャベツ・ニンジン・ホウレンソウのシチュー．
  • 晩 (2100): 米麦 0.8 合． 豆腐・ワカメのサラダ． ホタテ・ネギの炒めもの．

2008 年 03 月 07 日 (金)

• 0750 起床． 朝飯． コーヒー． 0920 自宅発． 晴． 0935 研究室着．
• 昨日のつづき， すなわち占部さんたちの高山湖沼の魚-動物プランクトン-植物プランクトン モデルの再計算に関する検討， ついでにモデルの説明図かき．
• 図をとりあえず作りはじめてみると， 環境変数のあつかいとか首尾一貫してないので， なんとかならんかなーと (お絵描きススめつつ) 5-6 回やりなおしてみたんだけど， やっぱり収束がうまく行かなくなったので， 昨日のモデルそのまま使うことにした．
• 図のほうは …… 占部さんが企画集会でそのまま投影資料として使うんで …… まあ， しょうがないので Linux 版 OpenOffice.org で手がき． BUGS code どーりにするとワケわからなくなるので， てきとうに手ぬきする．

• で， 最後に poorpoint ファイルとして「別名で保存」． ファイルアップロード． 計算結果に関して ECOAS wiki 上でまとめて， そこから図のファイルなんかもダウンロードできるようにする． メイルかいて一件落着．
• すでにとっくに正午すぎ． 昼飯． 修士論文の紙たばの事務室提出なんかも終了したようで， A 棟 8F はいたって静かで平和だ．
• 大雪 → 北大の学部二次試験の日程が一日ずれる → 採点したうけ要員の私の帯広出張をキャンセルせざるをえなくなった， で取り消しになった出張， こういう場合の旅費返納は JTB (北大の公費旅行ぎょーむの一部を下うけしている) のヒトがやってきて現金にて返納することになる (私の口座にすでに旅費がふりこまれているから)． 39300 円．
• 生態学会福岡大会初日 (3/14 金) まであと 7 日． 最低限のやることリスト:
  • (ひとまず終了) 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • (ひとまず終了) 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  占部さん下うけはひとまずナンとかなった． 次はベイズ企画集会の「コメント」とかいうやつあたりに着手したいんだけど ……
• その前にいくつかやらねばならぬコトが． まず， アリポスターなんだけど， 藤原君がいろいろとコメントをくれたので， 修正．
• で， 一件落着したと思ったら粕谷さんから R 自由集会 後のうちあげに関する問いあわせ …… うーむ， 今回はなんとなくいきあたりばったりでゴマかそうかと考えていたんだけど， (いつも人数不確定予約はイヤがられるので) …… まあ，金曜日の夜なのでかかる ノープランでコトにのぞんでしまうと 「夜の博多区・中央区あてのない R ユーザー集団散策」 を楽しむことになりそうだ．
• ということで， うちあげ会場さがしにとりくんでみる． 福岡国際会議場からの場所のちかさと 収容力の大きさ (小さい店は人数の急変に対応できない可能性が高い) でさがす． 私の場合， こういうのが食べたいとかあまりないんだけど， 今回はそういや福岡って小さい餃子なんかがうまかったような， とかやっぱり福岡の焼き鳥は他とはちがうような， といった趣味が反映されてしまって …… とりあえず今年も うちあげ店候補一覧 を作ることになってしまった． 粕谷さん・竹中さんに連絡．
• うーむ， 次に処理すべき件に関しては進捗せず …… 1855 研究室発． 1910 帰宅． 体重 70.2 kg． 晩飯の準備． 晩飯．
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 10 × 3 回． スクワット 100 回．
• [今日の食卓]
  • 朝 (0820): ヨーグルト． リンゴ．
  • 昼 (1300): 研究室お茶部屋． 米麦 0.6 合． ホタテ・ネギの炒めもの．
  • 晩 (2100): 米麦 0.8 合． ブナシメジ・サバの電磁波酒蒸し．

2008 年 03 月 08 日 (土)

• 0950 起床． なんだか寝不足ぎみ …… 朝飯． コーヒー． 洗濯． 怠業．
• 怠業ばかりしてられん …… ということで午後になってから放置してたトドマツ宿題にとりくむ． subject editor 指摘する (そして先日すでに計算してしまった) 母子間距離平均の分散なるものだけではおそらくダメだろう， ということでもう少し別の randomization の方法を検討してみる． しかし， ことごとくよろしくない．
• 気分かえるべく 1745 自宅発． 晴． ただしほぼ日没． 1800 研究室着．
• で， いろいろやってみたけどやっぱりダメでしたメイルを書きはじめてみる． さらに「なぜダメだったのか」について書いているうちに …… アタマの中が整理されて， もっとマシな改善策を思いついた．
• やってみたらうまくいった． 子供数 vs 最大母子間距離の単純な線形相関係数 は -0.072 と負の値になった． こんどは稚樹 419 個体の位置を random にい れかえて子供数 vs 最大距離の相関係数を検定統計量とする randomization test をしてみると

  > v.cor <- sapply(
  +    1:1000, function(i) cor(v.n, tapply(sample(d$dist), d$mother, max))
  +)
  > quantile(v.cor, prob = c(0.025, 0.975))
  
  2.5%   97.5%
  0.38021 0.74378
  

  となり， つまり子供の位置がランダムならば子供数が多いほど 最大母子間距離は広くなる，といったことがわかった． 考えてみればあたりまえのことだなぁ．
• といった計算結果にもとづき， 本文修正文案・cover letter 文案の作文． メイル送信． これで一件落着してほしい ……
• またべつのメイルかきなどやって， 2230 研究室発． 2250 帰宅． 晩飯の準備． 晩飯．
• 生態学会福岡大会初日 (3/14 金) まであと 6 日． 最低限のやることリスト:
  • (ひとまず終了) 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • (ひとまず終了) 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  さっさと企画集会のこめんととやらを準備したいところなのだが …… なんだかいろいろまた宿題が．
• [今日の運動]
  • うんどう休養日
• [今日の食卓]
  • 朝 (0930): ヨーグルト． シリアル．
  • 昼 (1330): イソップベイカリーのクロワッサン・バゲット．
  • 晩 (2320): 米麦 0.7 合． ネギ卵炒飯． キャベツ・ニンジン・ネギ・ハムのサラダ．

2008 年 03 月 09 日 (日)

• 0750 起床． 朝飯． コーヒー．
• …… いかん， また怠業もーどに落ちこんでいた． ぎょーむ日誌など書いてぎょーむもーどにぢりぢりと移行してみる．
• R 自由集会 のメイルかきとか． 竹中さんがすでに発表の準備されたのでそれに関するコメント． 私が考える「初心者にとっての R 作図難所」:

  ひとつめはデータの入力に関する部分です．端的に言って，テキストファイル
  とか CSV ファイルといった概念がよくわからない，と言いますか「ふつうのゑ
  くせるファイルと何がちがうのか」といったあたりがわからなくて困っている
  ヒトがときどきいます．一番極端な事例ですと，foo.xls というファイル名だ
  けを foo.csv に変更して「read.csv() でうまく読めません」といったことが
  ありました．データが read できないと何もできないので，このあたり「(デー
  タ入力はつねにゑくせるでなされるという前提のもと) R では『名前をつけて
  で保存』された text, csv file でないと読めない」といった単純化した説明
  (R commander などではまた別の手口が準備されてるみたいですが) が良いので
  はないかと思います．ここは文字どおり入口の難所だと思います．
  
  二点目は R の作図関数の (おそらく初心者にはやや奇怪な) col, pch, cex,
  lty, lwd といった属性指定の説明が補足できないだろうか，といったあたり
  です．これはキリがないので，今回はごく簡単にしか説明できないと思うので
  すが……やはり上述の 5 属性ぐらいはひとつのスライドでごく簡単に説明さ
  れてよいのかもしれません．また今回の「ループ重視」な話とすこしずれるの
  ですが，単純な plot でも
  
  http://hosho.ees.hokudai.ac.jp/~kubo/ce/IntroductionPlot1.html
  
  これの上から三番目のような図みたいに複数色を一度に指定できるといった点
  は R 作図初心者の皆さんにはいつも感銘をあたえるようです．
  

• あ， 自由集会といえば， うちあげの予約というか「人数確定してなくてごめんね」 予約フォーム投稿とか． こういう場合， 携帯電話をもっていない客 (すなわち私) なんかもイヤがられたりするんだよねえ ……
• 昼飯． 1445 自宅発． 晴． 1500 研究室着．
• またトドマツメイルかき． その他メイルかき．
• ここの研究院 (研究科) の web 管理委員会の wiki を PukiWiki で作る． PukiWiki のインストールはお手軽で， 5 分もあればとりあえず動くものが設定できる．
• さてさて …… 生態学会福岡大会初日 (3/14 金) まであと 5 日． 最低限のやることリスト:
  • (ひとまず終了) 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • (ひとまず終了) 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  …… それでは， 次の準備， ベイズ企画集会のこめんとのため， WinBUGS に空間モデルの事前分布を提供している GeoBUGS の user manual を熟読してみる．
• 読んでみてだんだんわかってきたのは， この GeoBUGS 世界は巧妙ながらなかなかヘンてこだ， というあたりか． 少なくとも下のふたつの Gaussian な conditional autoregressive (CAR) model の区別はつくようになった．
  • car.normal(): intrinsic Gaussian CAR model を使っているので， 「隣接」サイトとの共分散構造は定義されているけれど， 「全体」のばらつきなんかは定義されていないモデル． つまり空間モデルではない GLMM なんかとの対応は見えにくい．
  • spatial.exp(): こちらはクリギングなんかで使われるような， 近傍と全体の分散をつなぐ semi-variogram モデルを念頭においたもの． ただし nugget が定義できないようで， そうだとするとやはり空間構造なしの GLMM との対応は見えにくい (あえていえば Φ がかなり大)．
  …… といったところか? 明日になったらまた訂正するかもしれないけど．
• 2015 研究室発． 2030 帰宅． 晩飯の準備． 晩飯．
• [今日の運動]
  • またうんどう休養日 ……
• [今日の食卓]
  • 朝 (0830): ヨーグルト． シリアル．
  • 昼 (1330): スパゲッティー． タマネギ・エノキダケのクリームシチュー．
  • 晩 (2120): 米麦 0.7 合． ダイコン・シュンギク・ネギ・マイタケ・タラすりみの煮物．

2008 年 03 月 10 日 (月)

• 0740 起床． 朝飯． コーヒー． 0905 自宅発． 晴． なぜかとつぜん温暖 (日平均気温がプラス) な日々がはじまった． 0920 研究室着．
• 昨日のつづき， つまり ベイズ企画集会 の予習的な勉強のつづき． より具体的には WinBUGS における空間統計な関数 (もともとは GeoBUGS の関数だったもの) の正体解明について． 企画集会の皆さんにメイルかいてみる．

  CAR に関する説明ですが，最初に読んだときは深澤君の説明がなんだかヘンな のでは……と思ったのですが，川渡でご教示いただいた GeoBUGS user manual 読んでみて，ようやくこのあたり理解できました．「ヘン」なのはむしろ car.normal() で使ってる intrinsic Gaussian CAR model のほうだったようで すね (便利なモデルだとは思いますが)．

  intrinsic Gaussian CAR model は隣接する場所間の local な共分散構造だけ きめて，調査地全体の global なばらつき (random effects 全体のばらつき) は決めない，というちょっとふうがわりなモデルだと今さらながら理解しまし た．この単純化のおかげで，ばらつきパラメーターは一個指定するだけで MCMC 計算できているわけです．

  しかしながら空間構造を考えない GLMM なんかとの対応をきちんと説明するの は難しそうですね．私はこの点を誤解していて，普通の GLMM と car.normal() 使ったモデルが簡単に対応づけられる，と考えていました．

  しかしながら実際は簡単ではないようで，

  • 普通の場所差 GLMM: 調査地全体の random effects のばらつきは考慮している (あてはまりの良さと複雑さのバランスがとられる) が，どの二地点間も無相関であると仮定している (モデルで生成される確率場の variogram は水平線になる)
  • car.normal() を使ったモデル: 近傍地点間の相関は推定されるが，調査地全体の場所差のばらつきをあらわすパラメーターはない (モデルで生成される確率場の variogram の最大値 (sill) が不定 (tau が定数であっても) …… のはず)

  となっているように思います．「対応がつかない」と言ってるのは， car.normal() の weight や tau をどうひねくりまわしても「場所差あり，た だし地点ごと独立同分布の GLMM」にはならない，ということです．

  では，もう少し GLMM との対応がわかりやすいモデルはないかと考えると，お そらくこれが WinBUGS (GeoBUGS) では spatial.exp() を使ったモデルで phi が大きい状況 (距離とともに急速に無相関になる) というのが空間構造なしの GLMM に相当するかと思います．モデルで生成される確率場の variogram は原 点 (0, 0) をとおって上に凸の関数で最大値 (sill) に漸近する (これは spatial.exp() においては tau から決まるはず)，つまり空間統計学の教科書 とかでみかけるもの (ただし nugget はゼロ) になります．

  ベイズ企画集会においてこのあたりをごちゃごちゃ説明する必要はまったくな いと思いますが，今回つかっている空間統計モデルが一般的な CAR というより， (やや簡便版ということになるのでしょうか) intrinsic Gaussian CAR model ですよ，というのはスライドのはしっこにでも記載しておけばよいのかもしれ ません．

  上のように考えてみると，このあたりは意外とめんどうなような気も してきました (他にも，たとえば空間相関のある random effects とは「移動 分散制約の効果なのか，環境不均質性の効果なのか……時間変化がわからない データではその区別がつかない」問題とか)．

• で， メイルかきおわってみればすでにとっくに正午すぎ． 確率場 (random field) を使った統計モデリングはおもしろいんだけど， なかなか難しい．
• 生態学会さーばー雑用， 野外安全管理委員会関連． その他メイル書きもつづく．
• 昼飯． またメイルかき． なんでこんなにメイルばかり ……
• どうでもいいが母子里 Gibbs 林冠モデルは Agricultural and Forest Meteorology の web page によれば依然として With Editor 状態のまま止まっている． 最初の revision のときは督促メイルだしたら対応したんだけど， 今回はまったく動く気配がないなぁ．
• たのまれ WinBUGS 調査． dpois() がポアソン乱数 (整数) を生成してないんでないか (実数になる) 疑惑があって， ですね．
• まずは素朴に

  model
  {
  	Tau.noninformative <- 1.0E-2
  	for (i in 1:N) {
  		Y[i] ~ dpois(lambda)
  		y[i] ~ dpois(lambda)
  	}
  	log(lambda) <- beta
  	beta ~ dnorm(0, Tau.noninformative)
  }
  

  この例では Y[i] は観測データなんだけど， y[i] は (データにしばられない) 単なるポアソン乱数になるはず …… そうなっていると確認できた．
• それではもう少し質問に近い状況で試験してみる．

  model
  {
  	Tau.noninformative <- 1.0E-2
  	for (i in 1:N) {
  		Y[i] ~ dpois(lambda[Group[i]])
  	}
  	for (j in 1:N.group) {
  		lambda[j] ~ dpois(mean)
  	}
  	log(mean) <- beta
  	beta ~ dnorm(0.0, Tau.noninformative)
  }
  

  これは奇妙なモデルでポアソン分布の平均値である lambda[j] がポアソン分布にしたがう， といった状況だ …… で，この場合はやはりというべきか

  > post.mcmc[, "lambda[1]"]
    [1] 14.65 16.96 15.10 17.87 15.35 15.08 15.56 16.87 13.78 12.98 12.89 13.72
   [13] 14.21 14.55 13.61 15.93 16.36 15.92 14.20 14.46 16.02 15.52 13.64 15.63
   [25] 14.14 13.24 16.72 11.70 14.97 17.59 15.24 15.64 15.82 13.35 15.93 15.97
   [37] 16.46 17.37 16.46 13.53 14.96 14.85 15.29 15.38 14.56 13.79 15.47 15.34
   [49] 15.07 14.41 13.14 15.10 17.02 14.03 14.25 13.95 13.06 14.41 14.50 14.29
   ....
  

  と Gibbs sampling しているねえ． すなわち， WinBUGS 使い大学院生の (良くも悪くも) 独創的というほかない BUGS coding が WinBUGS 作者たちの想定した状況を超越してしまった， というコトのよーで (ひょっとしたらこれは WinBUGS の仕様どおりの動作， という可能性もゼロではないけど)．
• ちなみに， こういうゴマカシもまったく通用しない， と判明．

  model
  {
  	Tau.noninformative <- 1.0E-2
  	for (i in 1:N) {
  		Y[i] ~ dpois(r.lambda[Group[i]])
  	}
  	for (j in 1:N.group) {
  		r.lambda[j] <- 1.0 * lambda[j]
  		lambda[j] ~ dpois(mean)
  	}
  	log(mean) <- beta
  	beta ~ dnorm(0.0, Tau.noninformative)
  }
  

• 整数乱数しか生成しえない dcat(p[]) を使ってみると …… WinBUGS は「正しく」というべきか incompatibility を指摘しているらしいナゾの trap 窓 (おそらく dpois() のパラメーターに整数いれるな， といった趣旨なんだろう) をだして計算停止， と． なるほどね．
• もしも WinBUGS の MCMC 計算が Metropolis-Hastings 法だけで 駆動しているのであれば， 上のような問題はおこらない (したがってどうしても上のようなモデルで計算したければ， 自分で M-H 法のプログラムを書けばよい)． しかしながら WinBUGS の adaptive rejction sampler, slice sampler その他による Gibbs sampling は上のような ベイズモデルの MCMC 計算を拒絶している， ということだね．
• といった報告メイルを送信． (もとのモデルにはちょっとむりやり感もあるので) まだまだ統計モデリングには工夫の余地があるはずですよ， ということで．
• 1600 すぎに 北大情報基盤センター北館 にポスター印刷にいく． 事前に 昨年 3/16 の大判プリンター作業行程をチェック …… しかしうかつにもアカウント書いたメモをもってくるのを忘れて研究室にもどる． いやはや． ボケてるなぁ．

[情報基盤センター]

地環研の北の道路むこう． この中に， 例のアヤしげな北大ネット機器なんかも入っていたり．

[プリンター部屋の位置は変更された]

1F， 玄関はいってすぐ左の廊下奥に．

• で， 印刷作業． 昨年の不満足な結果 (小さすぎるポスターになった) を教訓に今年は慎重に． まず光沢紙だけどこれは長さにかかわらず「1 枚」 (というよくわからぬ単位) 945 円． 福岡大会 のポスター板は 90 x 210 cm とのことなので …… 私がいったときに装着されてた 36 inch (91 cm) ロール紙でちょうどよい．
• 大判プリンター用 PC の Adobe Acrobat でサイズあわせ． A0 は幅 84.1 cm なので 90 cm ぐらいに拡大するためには拡大縮小率とやらを 106% ぐらいにすればよい． 15 分ほどかかって印刷． 10 分ほどかけてインクを乾燥させる． まち時間は空間統計学本よみ． 右横を 2 cm ぐらいレイルカッターで切除． 今回は大きく印刷できてよかった …… しかし時刻はすでに 1720． 研究室にもどる．

[研究室もどる]

ポスター 完成． めでたしめでたし．

• 私はイヤだと拒否したんだけど， 大学院生に命じられて A 棟 8F 婦人トイレ天井蛍光灯を交換作業に従事させられた (なお野郎トイレは 7F)． A 棟 8F における雑用係という社会的地位なので雑用処理はやるけど， ここだけはちょっと ……
• 2010 研究室発． 晩飯． 2120 帰宅．
• 生態学会福岡大会初日 (3/14 金) まであと 4 日． 最低限のやることリスト:
  • (ひとまず終了) 一般講演: アリポスターの準備
  • R プログラミング自由集会: 前口上的な部分
  • ベイズ企画集会: コメント (まあ 1-2 分で終わるような)
  • (ひとまず終了) 山地生態学企画集会: 占部さんの計算したうけ，PP-ZP-Fish 系の再計算
  自由集会・企画集会ではどちらも 5 分以下のちょっとしたハナシ， をやればよさそう． さっさと準備しなくては ……
• [今日の運動]
  • またうんどう休養日 ……
• [今日の食卓]
  • 朝 (0750): ヨーグルト． リンゴ．
  • 昼 (1310): 研究室お茶部屋． 米麦 0.7 合． 野菜スープ．
  • 晩 (2030): 北 14 Bari Barl のポーク焼きトマトカレー．