ぎょーむ日誌 2005-02-13

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 合． ハクサイ・ブナシメジ・ワカメのシチュー．