ぎょーむ日誌 2005-02-(11-20)

2005 年 02 月 11 日 (金)

• 1030 起床． 朝飯． コーヒー．
• Perl なメイルフィルター SpamAssassin の update でぢたばたする． Vine Linux 3.1 の SpamAssassin 一式
  • perl-Mail-SpamAssassin-2.64-0vl2
  • spamassassin-2.64-0vl2
  • spamassassin-tools-2.64-0vl2
  の rpm をインストールしてる状態から， 最新の 3.0.2 に更新した (これは CPAN から直接取得)． インストールしたら，
  • sudo /etc/init.d/spamassassin restart が必要
  • bayes db version 2 is not able to be used, aborting! at... といった警告が出るので， このあたり みながら修正
  他にも .procmailrc のオプション指定の変更とか． このあたりは警告みながら変更すればよい．
• あの R Tips の舟尾さんの R 本 ``The R Tips - データ解析環境 R の基本技・グラフィックス活用集 -'' がもうすぐ出版される， と 霜月さんところ で知る． 自由集会 ペイジの参考文献に追加する．
• たしかに， 生態学会大会の会場で R 本とか販売するのは良いことだと思うし， どれぐらい売れるのかたいへん興味ある …… あ， しかし出店の申し込みしめきりは 昨日 ではないか．
• 1600 自宅発の北大構内走． 天気予報はずれて晴れているけどすごく寒い． 1650 帰宅． 体重 73.8kg． 昼飯 …… ってもう日没時刻なんだが．
• あ， いかん． 怠業．
• 1915 自宅発． 夜． 1925 研究室着．
• 先日からうだうだと一般化線形混合モデルの性質なんぞを調べてるもとになった biometry 投稿質問， よくよく読み直してみたら， 生存時間分布の問題だった …… われながらここしばらく阿呆阿呆度がかなり強まってきているようだ． ということで， これは他ならぬ 生存時間解析 である． えー， つまりまた私がカンちがいしており， 粕谷さんのご指摘のほうが正しかった， ということで．
• このあたり， 私はやったことないんだが …… Cox の比例ハザードモデルとかで， R で計算すんなら survival package の coxph() あたりか． ここで nested な構造を扱うには cluster() なる項を用いるらしい …… いや， ちがう． frailty() を使うほうが 正しいよーで?
• CRAN に kinship なる package あり， これに coxme() という ``Fit a mixed-effects Cox model'' な関数あるけど， これはバグあって動作しない …… 動かすと

  > coxme(Surv(time, status) ~ rx, data=rats, random= ~1|litter)
  Error in coxme(Surv(time, status) ~ rx, data = rats, random = ~1 | litter) : 
          couldn't find function "coxph.wtest"
  

  と止まるわけだが， 問題の coxph.wtest とはユーザーからは呼び出せない関数なのである． やれやれ．
• と思ったら， 作者 Jing Hua Zhao 氏の FAQ ペイジに何やら書いてある …… が， どうもこういう小細工は， 非ゐんどーづ版 survival package (ver.2.16) には通用せんよーで． そもそも coxph.wtest <- function(...) ... が見つからん ……
• ハンパなまま撤退． 2350 研究室発． 2400 帰宅． 晩飯． ちょっと生活周期がずれぎみ．
• [今日の運動]
  • 北大構内走 50 分間．
• [今日の食卓]
  • 朝 (1040): クラッカー．
  • 昼 (1710): うどん． ネギ・卵の炒めもの．
  • 晩 (2420): スパゲッティー． レトルトパウチドカレー．

2005 年 02 月 12 日 (土)

• 0920 起床． 朝飯． コーヒー．
• ふろりだの奥山さんから， 過去数日の GLMM 混乱ハナシに関するコメントいただく (生存解析その他でやるべきでは， という内容でした …… 昨日のぎょーむ日誌をアップロードする前だったもんで)．
• これに関連して， ``Negative binomial loglinear mixed models'' (負の二項分布の対数線形混合モデル) という「ひゃー」というような題名のふろりだ産論文おしえていただいた．

  Booth, J. G., Casella, G., Friedl, H. and Hobert, J. P. (2003). Negative binomial loglinear mixed models, Statistical Modelling, 3: 179-191.

  Poisson 分布をガンマ分布混合でばらつかせたのが 負の二項分布なんだけど， これにさらに正規乱数による random effects (個体差とかブロック差とか) を組み込んだもの． しかも Monte Carlo EM (MCEM) algorithm で計算します …… という何だかスゴいもので． 勉強になりそうだ．
• で， 奥山さんからのご指摘読みつつ返信書いてるうちに， 数日前 に「R の glmmPQL() では死亡率に関してブロック差の中に個体差がある， という状況ではうまく推定できん」 と騒ぎまくっていた問題が解決した …… これは 推定できなくて当然 だったのである． 嗚呼， どうしてこんなあたりまえのことに気づくのに えらく時間かかってしまったりするんだろう (答: 統計学をよく理解してないから)．
• まず， そもそもくだんの biometry 質問への回答は生存時間解析となるべきで， {生, 死} 二値パターンを説明する混合ロジスティック回帰 ってのはあまりよろしくない回答だろう． それはおいといて， glmmPQL() とかで「処理の差に加えて， ランダム効果としてブロック差があってその中に個体差ある」 というモデルで生成した乱数データの ブロック差・個体差をうまく推定できない理由は何だったのか， を考えてみよう．
• ブロックというのはこの場合だと生物をいれとく水槽かなんかで， その水槽のなかに個体が何匹かいる， という状況だ． ブロック差というのは， 実験処理を (同じ水準のなかでは) 均質にしてるはずなのに水槽ごとに何か違いが出てしまう， という差のことである． 個体差は， まあサイズとか遺伝子型の影響で， 同じ毒性にさらされても死にやすいとか死にやすくない， とか．
• 同じ水準の毒をあたえた (つもりの) ときに， ある水槽ではぜんぜん死ななくて， べつの水槽ではばたばた死んだ， という実験結果だったとしよう． このような「ブロックごとのばらつき」 は単純な二項分布モデルからは逸脱している， つまり overdispersion (過分散; 過大分散) が生じている状況である． 簡単に言えば， このときにこういう overdispersion の原因が水槽 (ブロック) 差なのか 個体差なのか原理的に見わけがつかないので， あわれな推定関数 glmmPQL() は錯乱した推定結果を出力していたのである． これはブロック数を増やしたり， あるいはブロック内の個体をいくら増やしても変わらない．
• えー， このあたりは数学の問題というより， まあ常識的なこととゆーか何とゆーか …… ある水槽では実験動物が他よりばたばた死んだ， という状況だったときに，
  • ブロック差説の例: これは水槽にまちがってちょっと多めに毒をいれたせいである (あるいは水循環装置の具合がいまいちだったとか， この水槽を置いてる場所がよくなかったとか)
  • 個体差説の例: この水槽には「毒に弱い」 遺伝子型の個体がたまたま集中してしまったせいである
  などともっともらしく並べてみたのは良いものの， 手もとにある {生, 死} データをいくらヒネくってもどちらが相対的により重要だったか， といったことは推定計算できない．
• これは， 処理の水準とブロックを決めてしまえば， そのブロック内個体の {生, 死} パターンは (個体差があろーがなかろーが) 必ず二項分布にしたがうので (オモテでる確率が 0.5 でないコインであっても， コイン投げ独立反復試行においては オモテでる回数は必ず二項分布にしたがう， というのと同じ)， 個体差の有無やその大小をいえないためだ． 個体差の有無大小がいえないので， ブロックごとに見られる死亡率の相違が， ブロックの影響によるものなのか， その中の個体差によるものなのか観察者にはわからない …… というだけの， ぢつにぢつにあたりまえなハナシで．
• 推定計算してる関数 glmmPQL() の視点にたってみると， 「個体差ぜんぜんナシ」 だろーが 「いや個体差とても大きい」 だろーが， ともかく何であってもデータへのあてはまりかたはほとんど変わらない (実際には乱数生成で作ったデータセットのちょっとした 偏りに左右されて「個体差ゼロ」⇔「個体差最大」 の間のどこかにさまよいこむ)， という状況なんで， ああいう無茶苦茶な計算結果を出していたわけだ．
• いやはや …… {生, 死} データはある期間においては同じ個体から一回しか得られない， というところが特徴的なんだよなぁ． これが複数回観測できる事象なら， ブロック差・個体差を同時に推定できる場合もあるだろう．
• たとえば， その実験動物が元気なときには 一日に何度も X 行動をとる (毒でカラダが弱ると回数が減る)， とする． このときはブロック内の個体数と観測日数が十分であるなら， その個体ごとの X 行動のカウントデータを利用することで， nested なランダム効果である ブロック差と個体差 (in ブロック差) を区別して推定できるはず．
• 1640 自宅発． 曇． 1650 研究室着． 途中で奇妙なモノをみてしまった．

[雪壁上のらくがき]

これは無頼な若者の無軌道な表現活動の帰結ではなく， 北海道警察交通課の 「ここに路駐してた車はレッカー移動しました」 表示なのである． うーむ．

• あ， またズボンのひざの部分が破れた． これは先日ぼろぼろになったのと同じ時期・同じ場所で買ったものである …… このような「壊れる時期がそろう」のは overdispersion の逆の underdispersion なのかもしれない (現象論的には……メカニズムのほうは不明)． まあ， たぶんこちらも破れるだろうと予感してたので， 新しいのは二本買ってあるわけだが． 私は Sony 製品もってないんでよく知らないんだけど， 世に言う ソニータイマー なんかも underdispersion 現象の例なのか?
• いやいや， このあたり Cox の比例ハザードモデルなら， 時間依存な破損率の増大， といったモデルで説明されるにちがいない． 勉強しなくては ……
• その前に， うだうだと上述の 「glmmPQL() 騒動」 の顛末を書く …… なンか最近， 過去のぎょーむ日誌に 「後記」をつける機会が多いような気がする． うかつな日々．
• いやはや． ともあれそもそもの原因は 「個体差だろーがブロック差だろーが， とにかく random effects はどんどん nest して片っ端から混合モデルにつっこめばよいのだ」 という 完全にまちがった 久保ドクトリンにあった， と結論するのが妥当なところでありましょう． いやはや．
• 月曜日の修論発表会にむけて， 修論生をさぽーとする利他的な活動がおこなわれてるよーで， 私も修論生用晩飯のおすそわけにあずかる．
• 2320 研究室発． 2340 帰宅．
• [今日の運動]
  • いかん， 今日も運動休養日になってしまった ……
• [今日の食卓]
  • 朝 (1010): 米麦 0.7 合． ブナシメジ・海藻・豆腐・煮干の味噌汁．
  • 昼 (1430): 米麦 0.7 合． モヤシ・ニラ・エノキダケの炒めもの． ブナシメジ・海藻・豆腐・煮干の味噌汁．
  • 晩 (2030): 研究室お茶部屋． 修論生用晩飯として他の院生がつくったもの． 米麦 0.5 合． ジャガイモ・ナガイモ・ニンジン・鶏肉の煮物． ヤマクラゲ・豆腐の味噌汁．

2005 年 02 月 13 日 (日)

• 0840 起床． 朝飯． コーヒー．
• いかんなあ， どうも怠業状態だ．
• のんびりもしてられんので， 生存解析時間モデルの勉強．
• 1225 自宅発の北大構内走． 雪がふったり晴れたり． 走ってるうちにアタマの中で勉強したことが整理される． 1320 帰宅． 体重 74.6kg． うわっ， いきなり重くなった． 昼飯． 洗濯．
• 1450 自宅発． 曇． 1510 研究室着．
• 生存時間解析の勉強のつづき． The R Books には生存時間分析 (あるいは event history analysis) の章がある． しかしながら， ここだけ読んでもどういうモデルが使われてんのか， まあまったくわからない．
• …… ということが昨日の時点で判明したんで， 今日の午前中は丹後俊郎 「統計モデル入門」 の「イベント発生までの時間の長さに関するモデル」を読んでいた． この本って生態学方面では「難解」ということで有名なんだけど (ただし推定読者数 < 10)， この章は丁寧な説明がなされており (まあコトバの定義とか一部すっとんでるところもあるんだけど)， 読んでいてたいへんよく理解できるものになっている．
• で， そのあとに The R Books のくだんの章にもどって読んでみると， 書かれていることがすらすらとわかる． この章は R の survival pacakage の使いかたに関しては 丁寧に説明されていて， そのあたりについてはたいへん参考になる．
• で， この生存解析においてブロック差や個体差がある場合はどうすれば良いか， という問題なんだが …… 一番簡単には， おそらく前に書いたようにこれは survival package 中の frailty() 関数を使えばよいのであろう． frailty とは robust の反対語である． で， その「頑健ではない」モデルとはどういうものであろうか． ばらつきが正規分布のとき Gaussian frailty みたいに書いてるから， frailty は「ふらつき」だろうか．
• どーもこのあたりのモデルの解説は Mayo 財団なるところの Terry Therneau なるヒトのものをよく見かけるんだが (R の survival package は Therneau が S-plus 用に作ったものを 移植したそーで) …… help(frailty) で引用されてる 論文，

  Therneau T.M., Grambsch P.M., Pankratz V.S. 2003. Penalized Survival Models and Frailty. Journal of Computational & Graphical Statistics 12: 156-175.

  をつらつらとながめてみる (これは 2000 年に投稿されたのに 2003 年出版という …… 内容は 2000 年時点での原稿 (これもネット上にある) とほとんど変わらんのに)．
• …… これはなかなか難物． Penalized partial likelihood (PPL) ってのがぢつはよくわからん (penalized の部分)． いろいろそのへんひっくりかえしてるうちに， これを持ち出す理由は計算が簡単になる， 推定量に漸近正規性があるんで検定やモデル選択とかに便利， というようなぽりしーで使われてるらしい， ということはわかってきたんだが．
• Agresti の Introduction ぢゃないほうの Categorical Data Analysis 本にそれっぽい説明のってた． Penalized log likelihood である L*(p) が L*(p) = L(p) - a(p) と定義されているのであれば， これはべいぢあん推定のつもりで， パラメーター p の事前分布が exp(-a(p)) に比例するカタチになってると考えればよろしい， か …… ふーむ． これはすごくおもしろい．
• ついでに glmmPQL() の PQL (penalized quasi-likelihood) に関する説明も同書にはあり， このあたりの部分も (いまさらながら!) 少しわかってきた． 正規分布 random effects な GLMM の周辺尤度を Laplace 近似で分解すると， なぜかしら quasi(-log?-)likelihood とおつりの項がでてきて， そのおつりの部分が penalty を表すようになってる， ということのよーで． PQL 推定値は bias とかあるんで注意せよ， とのこと．
• 生存時間解析あれこれ， さっさと一区切りつけたいとは思うんだけど， なかなか思うようにススまん． 2100 研究室発． 2110 帰宅． 晩飯．
• [今日の運動]
  • 北大構内走 50 分間．
• [今日の食卓]
  • 朝 (0930): 米麦 0.7 合． モヤシ・ニラ・エノキダケの炒めもの． ブナシメジ・海藻・豆腐・煮干の味噌汁．
  • 昼 (1350): 米麦 0.7 合． ブナシメジ・海藻・豆腐・煮干の味噌汁．
  • 晩 (2220): 米麦 0.7 合． ハクサイ・ブナシメジ・ワカメのシチュー．

2005 年 02 月 14 日 (月)

• 0830 起床． 朝飯． コーヒー． 0920 自宅発． 雪． 0930 研究室着．
• R の survival pacakge を調査するプログラム作りあれこれ．
• 1140 どーにかこーにかできた． 下の図は plot.survfit で作図したもので 毒性によって生存曲線が異なる， という 架空のデータ (これは自作プログラム で生成したいわばシミュレイションの結果)． + は「うちきり」 (censored)， つまり実験期間 30 日間を生きのびたことをあらわす．

• この図のもとデータはこういうかんぢで．

  > Surv(data$time, data$status)
    [1] 30+ 30+ 30+ 30+ 30+ 30+ 30+ 30+ 18  30+ 30+ 30+ 30+ 30+ 30+ 30+ 30+ 30+
   [19] 30+ 10  30+ 30+ 30+ 30+ 30+ 27  21  30+ 30+ 30+ 30+  6   9  28  30+ 10 
   [37] 30+ 30+ 27  27  12  30+  4  18   7  27  30+ 20  30+ 30+ 30+ 30+ 14   6 
   [55] 30+ 10  30+ 30+ 27  14  13  30+ 30+ 30   8  30+ 13  20   2  30+ 30+  9 
   [73] 18  30+ 14  30+ 29   6  30+ 30+ 30+ 10  30+ 12   6  30+ 22  12   2  11 
   [91] 12  10  30+  7   9   4  16   1  20  30+  8   1   4  15   2  30+  2   3 
  [109]  2  27  16   8   5  30+  6  30+  8  21  30+ 30+ 20  17  30+  7  11   8 
  

  6 水準 × 7 ブロック × 3 個体． ここでも + は「うちきり」． ブロック間の差がそれなりにあり， 個体差もびみょーにある (どちらもガンマ乱数)．
• これを Cox の比例ハザードモデルへのあてはめ (関数 coxph() による) はこういうかんぢで． method = "exact" を指定してるけど， これは無視されて Breslow 近似になっとるような．

Call:
coxph(formula = Surv(time, status) ~ toxicity + frailty(plot), 
    data = data)
  n= 126 
              coef   se(coef) se2   Chisq DF p      
toxicity.L     2.090 0.478    0.478 19.13 1  1.2e-05
toxicity.Q    -0.874 0.452    0.452  3.74 1  5.3e-02
toxicity.C     0.177 0.374    0.374  0.22 1  6.4e-01
toxicity^4    -0.647 0.315    0.315  4.22 1  4.0e-02
toxicity^5    -0.214 0.297    0.297  0.52 1  4.7e-01
frailty(plot)                        0.00 0  9.3e-01
           exp(coef) exp(-coef) lower .95 upper .95
toxicity.L     8.083      0.124     3.168     20.62
toxicity.Q     0.417      2.396     0.172      1.01
toxicity.C     1.194      0.837     0.574      2.49
toxicity^4     0.524      1.910     0.282      0.97
toxicity^5     0.808      1.238     0.451      1.45
Iterations: 6 outer, 18 Newton-Raphson
     Variance of random effect= 5e-07   I-likelihood = -283.3 
Degrees of freedom for terms= 5 0 
Rsquare= 0.259   (max possible= 0.992 )
Likelihood ratio test= 37.7  on 5 df,   p=4.36e-07
Wald test            = 27.9  on 5 df,   p=3.83e-05

• この推定結果は， こちらがデータ生成プログラムで与えてたのと …… まあ， だいたい傾向としては正しいけど， 推定値はぜんぜん正確ではないように見える …… 例の選点直交多項式のおかげでわかりにくいものになっているが．
• しかし， ちょっとこれはズレが大きすぎるよーな …… なんかまた私がしくじってるんだろうか． いろいろと調べてみる．
• 上記 Cox 比例ハザードモデルはセミパラメトリックモデル である． じゃあ， パラメトリックモデル (指数関数的な減衰) を仮定した推定を survreg() 関数でやらせてみると ……

Call:
survreg(formula = Surv(time, status) ~ toxicity + frailty(plot), 
    data = data, dist = "exponential")
             Value Std. Error      z         p
(Intercept)  3.812      0.159 24.007 2.35e-127
toxicity.L  -2.068      0.476 -4.348  1.37e-05
toxicity.Q   0.880      0.451  1.949  5.14e-02
toxicity.C  -0.175      0.374 -0.469  6.39e-01
toxicity^4   0.637      0.314  2.028  4.25e-02
toxicity^5   0.207      0.297  0.698  4.85e-01
Scale fixed at 1 
Exponential distribution
Loglik(model)= -293.7   Loglik(intercept only)= -312.6
        Chisq= 37.81 on 5 degrees of freedom, p= 4.1e-07 
Number of Newton-Raphson Iterations:  6 23 
n= 126 

• うーむ …… どちらも推定法でも toxicity の影響がほぼ対数線形ってのは (これは乱数データ生成のときにそう設定したわけなんだが) ちゃんと推定できてるわけだが ……
• 昼飯． 窓の外は吹雪．
• 1300 から 1600 ごろまで生態科学専攻化の修論発表会． まあ， みなさん無事に終了， かな． 某研究隠匿師匠から知恵をさずけられて 「あれ? 液晶プロジェクターにつながらない」 と時間かせぎを試みたのがいたけど， けっきょくその 5 分間は発表時間に含まれなかったので ムダにおわった次第で．
• 研究室にもどると， 以前の査読の件であれこれと編集部から注文が …… いやはや， 私これ読まされるの三度めなんですが．
• とりあえず最小限必要な論点だけ作文して送信． 撤退． 1850 研究室発． 雪． 1900 帰宅． 運動． 晩飯．
• また R の生存解析 package 調査のつづき． 乱数生成シミュレイションにちょっとまずいところが見つかったので， 修正してみる． 上に示しているのは修正済のデータ・推定結果である． しかし， 推定結果があまり良いものに見えない， というところは変化ナシ． なんかまちがえてんのかなぁ．
• [今日の運動]
  • エアロバイク 50 分間． 右に続いて， 左ペダルのボールベアリングもイカれぎみ． 教訓: エアロバイクはペダルで選べ．
• [今日の食卓]
  • 朝 (0840): 米麦 0.7 合． ハクサイ・ブナシメジ・ワカメのシチュー．
  • 昼 (1230): 研究室お茶部屋． 米麦 0.7 合． ハクサイ・ブナシメジ・ワカメのシチュー．
  • 晩 (2210): 米麦 0.7 合． ホウレンソウ・ハクサイ・ブナシメジ・ワカメのシチュー．

2005 年 02 月 15 日 (火)

• 0700 起床． 朝飯． コーヒー． 0850 自宅発． 晴． 0900 研究室着．
• 午前中は専攻科の修論発表会場． 生態遺伝講座と低温研の院生たちの発表， 皆さんなかなか水準高くてしかも楽しめる．
• その修論発表聞きに来ていた苫小牧の村上さん・平尾君と 「きゃら亭」 で昼飯． そのまま苫小牧 leafminer 相談会となる．
• 1500 講座セミナー， 今日は M1 牛原さんで屋久島で調べた樹木の葉寿命が 何に影響されていそうか， というハナシ． 私は多種系がニガ手とは前に書いたとーりなんだけど， この問題に関しては先週のうちに牛原さん・塩寺さんから 教えていただいて予習してたんで， まあそれほど動揺せずにすんだ． いわば Leigh しんどろーむ とでも言うような， 葉寿命と N 含量その他との相関関係をひねくる一分野が存在し， このわくぐみのなかで大雑把なハナシをしてるぶんには， あまり種間差とかをまぢめに考えなくても …… などと楽観してるとあとで苦闘することになるかもしれんな．
• セミナー後， このあたりに関連する院生たちが議論してるのを 横で傍聴してまた勉強する．
• 大阪大会の次の生態学会大会はたいへんなことになりそう …… とりあえず， 私は server 管理者 (OS など非こんてんつ系のみ?) でよい， ということのよーなんだが． 学会事務局， 世間の常識からはかなり逸脱してるんでは．
• よくわからぬ申請まきこまれ事件， そのようなものに関与してしまった愚行の代償として， 「何をやるべきか全くわからない研究あぴーる」 準備の分担を申しつけられる． 数時間を費す． で， そのあげくに下のような図をねつぞうすることに．

• えー， 樹冠に個体別の なんだかキモチ悪い色がついていますけど， これはまあ 「カラーひよこ」 みたいなもんで， 中身はくろーなるなぱいぷ樹木 そのものでございます．
• カラーひよこばて． 2200 研究室発． 2220 帰宅． 体重 74.6kg． 重い． 晩飯．
• [今日の運動]
  • うんどう休養日．
• [今日の食卓]
  • 朝 (0730): 米麦 0.7 合． ホウレンソウ・ハクサイ・ブナシメジ・ワカメのシチュー．
  • 昼 (1250): 「きゃら亭」 ランチ． 和風かにたま．
  • 晩 (2250): パン． ハクサイ・ショウガ・ワカメ・豆腐のスープ．

2005 年 02 月 16 日 (水)

• 0730 起床． 朝飯． コーヒー． 0910 自宅発． 晴． 0920 研究室着．
• 一昨日 のつづきで， R の survival pacakge の推定結果を解明する試行錯誤のあれこれ． われながら理解がのろっちいと思うんだが …… どうにかこーにか推定値の意味をとれた． 図で表現するとこういうかんぢだ．

• 破線は比例ハザードモデル (coxph()) ではなく， 指数関数モデル (survreg()) で推定させたものである． まあ， どちらのモデルでも (時間依存項) * exp(- 線形予測子) となることには変わりないんだが．
• 蛇足ながら， 私が乱数データ生成用に与えた fixed effects だけで生残曲線を描くと下の実線のようになる． 実際の生残曲線とのずれは random effects の与えかたが不適切であったため， のよーで …… これは推定値が違ってくるはずだ．

• ともかく， ハザード関数を構成する線形予測子の推定値と 観測データの対応関係はこれで理解できた． 生残データ生成のシミュレイションのほうはまだ工夫せねばならんが．
• オホーツク海に面する 興部 (おこっぺ) から講座出身者がやってきたのでオホーツク雑談．
• ともあれ， 先週から気になってた biometry への質問， 生残時間問題は survival pacakge で解決できるんでは， といった教唆投稿． ひとくぎりついたんで昼飯．
• 生残解析に関する質問メイルなどが来てたんで， とりあえず知ってることについては返信．
• 修論発表後も研究に邁進する M2 たちの一人， 平林さんから (人工) 花序内ハチ訪花パターンを 実験データから自動作図するプログラムの改良依頼． とりあえず， 安直な改良版を作ってアップロード． R はこういう自動作図にも有用です．

• 昨日のとある修論発表に関して， コメントのメイルを送ったところ， 反撃の返信いただいてひるむ． ぢりぢりと後退しつつご説明メイルを返信．
• とりあえず目の前の懸案はざっと片づいたんで， アカマツ問題とかにとりくむべきなんだが …… なんとなく統計学勉強方面に流されてしまう．
• 昨日， 苫小牧の平尾君から， 最近の Ecology Letters に Bayes 統計学の Review がのってると教えてもらったんで， そのあたり調べてみる．

  Ellison, Aaron M. (2004) Bayesian inference in ecology. Ecology Letters 7 (6), 509-520. (abstract)

  Clark, James S. (2005) Why environmental scientists are becoming Bayesians. Ecology Letters 8 (1), 2-14. (abstract)

  解説としては後者のほうがわかりやすいかも． しかし， この Jim Clark のほうのタイトル， Why environmental scientists are becoming Bayesians ってのは …… (少なくとも本朝では) 最尤推定法すらぜんぜん普及してないのに， いっきにべいぢあんな時代にされてしまってもなぁ．
• さらに Jim Clark は Models for Ecological Data なる本をちかぢか出版するつもりらしく …… 目次を見るに， これはぢつにすごい教科書だ!
• 1850 研究室発． いろいろな人々のすけぢゅーりんぐの都合で， 本日「追い出しコンパ」である． 一次会は北 12 の「ひだまり庭」， 二次会は北 18 の「がらん」． 楽しくうるさい場でありました． 2440 帰宅．
• [今日の運動]
  • 今日もまたうんどう休養日か ……
• [今日の食卓]
  • 朝 (0750): パン． ハクサイ・ショウガ・ワカメ・豆腐のスープ．
  • 昼 (1330): 研究室お茶部屋． 米麦 0.7 合． ホウレンソウ・ハクサイ・ブナシメジ・ワカメのシチュー．
  • 晩 (1900): おいだしこんぱ． 一次会は鴨キムチ鍋など． たくさんたべてしまった．

2005 年 02 月 17 日 (木)

• 0750 起床． 朝飯． コーヒー． ちょっとメイル書き． 昨日の生存解析がらみの続きで． 0920 自宅発． 雪． 0930 研究室着．
• 小林さんからアカマツに関する新しい (!!) データおくっていただいたので， とりあえず作図してみる． 作図．

[広島県某所のアカマツたち]

• しかし調査場所 (色) と相関のありそーな形質がけっこうありそうな気もする (下の図は R の pairs() によるもの)． まあ， 野外調査データではいろいろな原因による 交絡 (confounding) が回避しがたいわけだが．

• 1250 昼飯たべるためいったん帰宅． 昼飯． 1400 ふたたび研究室にもどる．
• 一時間ちょい， 講座のプリンター (CANON LBP-2810) 不調問題にとりくむ． カラー印刷するとペイジの「上側」あたりに 過剰なトナーがべったりつく， という問題なんだが …… まあ， 原因はハードウェアの故障か， カートリッジ (再生もの) のセレンドラムの不調だろう． 後者のほうがありがちに思えるんだけど， 三色のうち二つが同時にダメという現状を説明できてるのか， よくわからない．
• CANON のプリンター FAQ みてもよくわからん． 「カートッジ振り」 とかやってみたけど効果ナシ． 保証書・マニュアルのたぐいを探しまわったけど見つからず， そうしてるうちにそのあたりの紙切れを (現在育児休暇中の) 雪野さんにおしつけたような気がしてきた． とりあえず， お尋ねのメイル．
• 雑用ついでに 1F の研究科事務室に郵便物とりにいった …… あ， 来週末の試験監督についての紙きれが． 一点だけ不幸中のさいわいみたいなことは， どうもこれはセンター試験とは異なり， 一日中拘束されるわけでもなさそう， というような ……
• 新アカマツデータの「補償」モデルあてはめの推定計算． といっても， 一種の変則的な折れ線回帰を最尤推定法でやってるだけ なんだが． 意外と時間かけてしまう …… というのも， いぜん実験個体群用につくった R の推定プログラムを モジュール化する， というあたりでぢたばたしてしまってですね． とりあえず， 図のアップロード．
• で， 上の図をみればあきらかなよーに， 点どもは「折れ曲がった」線分上にのってるんで， 単なる直線だとこじつけたモデルと比較した場合， パラメーター増やした「補償」モデルのほうがよい， という結果になった． しかし， この図の縦軸っていったい何なんだろう， といまさらながらわけわからなくなったわけだが (実験個体群の場合だと， 単純に個体全体の葉量なんだけど)． あとでお尋ねしなくては．
• とりあえず撤退． 2040 研究室発． 2100 帰宅． 運動すこし． 晩飯．
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 5 × 3 回．
• [今日の食卓]
  • 朝 (0820): リンゴ．
  • 昼 (1320): 米麦 0.8 合． ハクサイ・ショウガ・ワカメ・豆腐のスープ． 昨晩「ひだまり庭」の残りものとして回収してきた チンゲンサイ・鶏肉の炒めもの．
  • 晩 (2220): 米麦 0.9 合． 海草スープ． ハクサイ・ニラ・コンニャクの炒めもの．

2005 年 02 月 18 日 (金)

• 0700 起床． 朝飯． コーヒー． 0850 自宅発． 晴． 0900 研究室着． 輪読会の予習．
• 1000 から輪読会のハズ …… だったんだが， 担当の M1 北村君が階下のコピー室からもどってこないんで， 様子を見にいったら階段したに倒れて気絶してた． 事故発生は 1005 ころか? 後頭部が階段にアタってるけど， とりあえず命に別状はなさそう． まわりの院生をよんでどうしようかと相談してるうちに， 意識をとりもどした．
• しばらくそこに放置して様子をみてから， セミナー室のソファに運んで 30 分ほど様子をみる． 工藤さんの車で北大病院に運んでもらうことにする． 1055 北大病院着． 「救急」 と書かれた入口から入ったけど， そっちは昼間は受けつけていなかったんで， ふつーの外来の窓口へ． 私が診察申込書みたいなの書いたんだけど， 字が汚いとおこられてしまった．
• 2F の脳神経外科 (neurosurgery) なるところに移動． 1128 最初の診察． とりあえず CT スキャン撮ってこい， とのことでまた 1F に移動． 病院内のオンライン化は北大一般事務関係よりはよほどマシらしく， CT スキャンの結果なんかもオンラインでやりとりされるらしい． 患者の個人データは何やらカードみたいなのに入っていて， ですね．
• CT スキャン部屋の前で一時間ほど待たされる． 1240 撮影開始， ですぐに終った．
• また 2F 脳神経外科の前で 20 分ほど待つ． 1305 すでに自力で歩けるようになった北村君が 診察室に行く． 私は待合室で待つ． 1315 診断終了， でとくに問題なかったらしい． 1330 北大病院発． 考えてみれば北大病院入ったのはじめてだなぁ． 1335 研究室にもどる． 昼飯．
• 今回の事件は
  1. ふせっせいなる大学院生の生活
  2. 貧血で階段下で倒れる (階段から落ちたわけではないらしい)
  3. アタマぶつけて脳震盪
  という顛末だったらしい． 大学院生のみなさん， 世間一般とはあまりズレてない 規則ただしい生活をしてください．
• それにしても， 輪読会で北村君担当の日になると本人に何かがおこるんだよな． 特異現象というか． 交通事故だの， 親類の不幸だの， インフルエンザだの， 脳震盪だの ……
• 本日の天気， 晴れたり吹雪になったり．
• またアカマツ新データを以前のモデルにあてはめてみる作業． これまた意外にも， というか …… シュート内の重量分配をみると， 場所によって違うんでは， というのが最良モデルになってしまった． シュート重量の範囲に偏りがあるんで， これが妥当な推定結果なのかどうか， びみょーではあるんですが．

• というようなことを報告したり， また別のメイル書きをやったり ……
• 統計学勉強とかにとらっぷされる …… 2050 研究室発． 2100 帰宅． 晩飯． いきなり寝る．
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 5 × 3 回．
• [今日の食卓]
  • 朝 (0720): パン． ハクサイ・ニラ・コンニャクの炒めもの．
  • 昼 (1345): 研究室お茶部屋． 米麦 0.8 合． ハクサイ・ニラ・コンニャクの炒めもの．
  • 晩 (2220): 米麦 0.8 合． ハクサイ・ニンジン・豆腐のカレー．

2005 年 02 月 19 日 (土)

• 0740 起床． 朝飯． コーヒー． 0930 自宅発． 晴． 0950 低温研着．
• 1600 すぎまで生態学会の 北海道地区会 の発表をきく． 1000-1500 は修論生の発表で， たいはんはすでに一度聞いたものなんだけど， なかなか楽しめた．
• 1610 低温研発． 吹雪． 帰宅ついでにそのへんの古本屋をうろうろと． 1730 帰宅． 晩飯．
• [今日の運動]
  • 腹筋運動 30 × 3 回． 腕立ふせ 5 × 3 回．
• [今日の食卓]
  • 朝 (0800): 米麦 0.7 合． ハクサイ・ニンジン・豆腐のカレー．
  • 昼 (1240): 北大生協北部店 (土曜日は 1400 まで営業)． サバ生姜煮． 豆腐ビビンバ． ライス M． 豚汁 M． 902 キロカロリー …… うーむ， とりすぎてしまった． 514 円．
  • 晩 (2100): 米麦 1.0 合． ハクサイ・ニンジン・豆腐のカレー．

2005 年 02 月 20 日 (日)

• 0900 起床． 朝飯． コーヒー． 洗濯．
• 怠業． 昨日みつけた SF 文庫本数冊を昨晩以来読んでる． いまさらながら， 「3001 年」 とか． これはクラーク老の いつもの太陽系冒険譚ではあるんだけど， 3001 年には技術が進歩しすぎていて …… セコい「とんち」航宙技法が登場しないところが悲しい． 私にとってこの「2001 年」しりーづ の最大の楽しみはそれだったわけなんだけど． 運動量交換カタパルト (2001 年)， 木星大気圏制動シールド (2010 年)， 位置エネルギー併用ブースター (2010 年)， 「行った先で補給」法 (2010 年 & 2061 年) ……
• 1315 自宅発の北大構内走． 曇． 1410 帰宅． 体重 74.6kg． あいかわらず， 重い． 昼飯．
• Linux の X (XOrg-6.7.0-0vl6.1) 上で窓わくとか仮想画面を管理する window manager が交換できないか， いろいろと試してみる． Gnome とか KDE (これらは WM というより「ですくとっぷ環境」かもしれんが) みたいなごちゃごちゃしたのは使う気になれんので， 「軽め」がウリの xfce4, flubox, afterstep, などなど …… けっきょく， どれも電気鼠操作を重視しすぎていて， 使う気になれなかった． いちいちポインターでちっこいボタンの上にもっていって， かちかちと悠長にやってられるかっ ……
• ということで， 鍵盤操作 (しかもばきばきとカスタマイズした) で何でもできてしまう WindowMaker (WindowMaker-0.80.2-0vl3) をしつこく使っている． 私が鍵盤たたいてるだけなのに 作業スペイス切り換え・ターミナル開閉・窓伸縮 などなどすばやく処理できてしまうの見ると， 電気鼠依存なヒトたちはこれを不思議に思うかもしれない．

[WindowMaker]

…… を特にぢみにカスタマイズしてしまうのがワタクシ流である． なンかさー， Vine Linux なるディストリビューションは， この「鍵盤操作派」には最良な WindowMaker を見捨てつつあるんだよねぇ …… まぁ， 自分でコンパイルすればいいんだけなんだけど．

• あ， いかん． 本日は完全怠業のうちに終わりつつある．
• [今日の運動]
  • 北大構内走 55 分間．
• [今日の食卓]
  • 朝 (0930): パンケイキ．
  • 昼 (1440): うどん．
  • 晩 (2000): スパゲッティー． トマト・ホウレンソウ・ブナシメジ・豆腐のソース．