ちなみに、今は暗黙的に離散時間でハザード率を考えていることに注意してください。次のような確率の掛け算で得られます。ハザード率を連続時間に拡張する知っているハザード率の定義と違う!と思った方少々お待ちください。今の話をそのまま連続時間に拡張すると一般的なハザード率の定義になります。用いて次のように表せるはずです。これもさっきの理屈と同じです。解釈もさっきの理屈と同じです。一点違うのは時間が離散から連続になっているという点ですね。これは見方を変えると次のように見ることも出来ますね。もとでそのため、よくハザード率の説明では次のような定義が成されます。これがよく見かける連続時間におけるハザード率の定義です。2種類の観点からハザード率について説明しました。ハザード率だ!と言われても、だからどうした、つまりそれをどう扱えばよいのだとなりかねません。そこで確率分布を何らかの方法で特定しましょう。離散時間でものを考えるのならば、話は中学生の時の確率レベルで構いません。分かれば良いわけです。よって、ハザード率の推定量を次のように置くことにしましょう。このように単純な相対度数を取った推定方法を経験分布による推定と呼びます。確率論の定義に則れば、サンプルサイズ無限大で相対度数は真の分布と定義上等しいので、真っ当な推定方法です。これをもう少し大きな枠組みでノンパラメトリックな推定と呼んだりします。一般的にノンパラメトリックな推定はパラメトリックな推定と比べるとefficientではないという欠点があります。パラメトリック推定ノンパラメトリック推定は分布の仮定をかなり広く取った(正確にはパラメータ空間を広く取った)推定方法でしたが、もう少し分布を限定してしまうことも出来ます。これを片側指数分布と呼びます。時間に依存せず、ハザード率は等しいというような分布の仮定になります。これは生存解析の文脈でハザード率の説明をする際には具体例として必ず取り上げられる特徴的な仮定です。打ち切りデータにおけるハザード率と生存確率ただ、現実はそうそう甘くない。途中で離脱してしまうということも考えられます。この場合、少なくとも〇〇歳まで生きていたことはわかるが、それ以上はよくわかりませんというデータになるわけです。具体例それでは、打ち切りデータの具体例を見てみます。RのパッケージsurvBootOutliersにはget.whas100.dataset()という関数でwhas100という打ち切りデータを取ってくることが出来るので、それを使います。1病気の患者のデータですね。statusが0の場合は生きているか、観察打ち切り(退院、転院等)のどちらかの打ち切りデータです。statusが1の場合はもう亡くなってしまっていると考えます。timesは入院から(死亡or離脱or生存)の期間を表します。survivalパッケージのsurvを利用すると、timesとstatusの関係を整理することが出来ます。数字はtimesです。status=0の場合は数字の横に+がついています。このような打ち切りデータの場合、先ほどのようなハザード率・生存確率の推定では、打ち切りを考慮できていません。2719+ちなみにこのデータはデータ数100です。timesに重複があるのがちょっとややこしいですが、一旦説明の都合上完全に順序があると考えてください。6日で2人亡くなっているわけですが、時間軸上先に一人目が亡くなって、次に二人目が亡くなったと考えてください。考えているものとします。時点での生存確率はノンパラメトリックに(経験分布的に相対度数で)考えると次のようになりますね。この生存確率の推定量をカプランマイヤー推定量と呼びます。生存曲線このカプランマイヤー推定量は次のようにplotして確認してみるのが一般的です。ちなみに下のプロットを生存曲線と呼びます。lty=c(1,2))ちなみに、普通生存曲線は比較する二つの群をプロットして、どちらの生存確率が高いかといったことを確認するものなので、survfit関数では比較する群を識別できるようformulaで記述します。今回は男女で分けてみました。理由は特にありません。信頼区間カプランマイヤー推定量の信頼区間もsurvfitでは求めることができます。plot(survfit(Surv(Data$times,Data$status)~Data$gender),lty=c(1,2))lines(survfit(Surv(Data$times,Data$status)~Data$gender,conf.type="log"),mark.T=F,conf.int=TRUE,col=c("black","red"))legend(locator(1),c("male","female"),lty=1,col=c("black","red"))ただし、この信頼区間の推定、恐らくグリーンウッドの公式を用いた方法が取られていると思われますが、この公式、生存確率が0か1に近い場合過少推定が起こると指摘されています。以下参照https://www.jstage.jst.go.jp/article/dds/30/5/30_474/_pdfパラメトリックアプローチカプランマイヤーはノンパラメトリックな推定方法でしたが、パラメトリックに考えることも当然可能です。Rではsurvreg関数が用意されており、distに分布を指定すればそれに合わせてパラメトリックに生存解析を行うことが可能です。2群の検定(ログランク検定)2群のハザード率・生存確率を比較するにあたり、2群のハザード率が等しいか否かを検定する方法にログランク検定があります。ログランク検定の帰無仮説は次のようになります。(生存確率が等しい)Rではsurvdiff関数によって行えます。結局ハザード率って何なのここで一度整理しておきます。ハザード率の確率論的解釈についてはこの記事の初めでお話しましたが、結局のところノンパラメトリックな観点から見たハザード率は、確率密度の代わりであると考えることも可能です。ちなみにパラメトリックな解析においては、パラメータの推定値を利用して逆算的にハザード率を求めることも可能です。比例ハザードモデル最後に、2群じゃなくてもっと多くの条件からハザードを確認し、多くの共変量の中から特に影響の大きいものを探したいといったモチベーションがあった場合、今までの手法ではうまくいきません。そこでcox比例ハザードモデルというセミパラメトリックな方法を取ります。比例ハザードまず、比例ハザードについて説明します。与えられているとします。時点でのハザード率は次のようになると仮定します。ベースラインハザードと呼ばれます。共変量に全く依存していないことに注意してください。要は、ハザードは比較実験のための情報が同じであれば一定であるという仮定をおいているわけです。医療実験の枠組みで解釈するならベースラインハザードはプラセボ群ということになりますね。無いことにも注意です。このモデルでは、時間によって共変量の影響が変化しないという仮定も入っています。置くことで、これをcox比例ハザードモデルと呼びます。このモデルの尤度を解くことは難しいので、部分尤度から近似推定量を求めます。落とし込めれば基本的には解析的に求めることが出来るわけですが、近似解を求める手法も提案されています。coxphsurvivalパッケージではcoxph関数でcox比例ハザードモデルを推定することが出来ます。efronやbreslowを与えると近似解を求めます。近似解はあくまで近似なので、計算量的な問題を気にしないのであればexactの方が良いと思います。デフォルトはefronです。p=1e-04bmiや年齢が影響を及ぼしていることが確認出来ました。検定によるチェックで、cox比例ハザードモデル、かなり仮定が強い印象があります。この仮定が妥当かどうかを検定でチェックしておく必要があります。survivalではcox.zphによって行えます。この検定については私は何をやっているのか知らないのですが(リサーチ不足)とりあえず先人の知恵を盲信して試しておきます。0.645のため、有意差はありません。cox比例ハザードを棄却出来るほどの情報はないと判断出来ます。今回も少し難しめの話。専門にする大学院生がかいてます。経済学部出身の文系です。時々趣味や日々の生活についても書きます。専門以外のことに関しては、単なる自分のアウトプットのため、不正確な記述があります。勉強友達・勉強ブロガー友達等々も募集中なので、気軽にコメントとかしてください~。読者になるプライバシーポリシープライバシーポリシーお問い合わせお問い合わせフォーム注目記事スポンサーリンクオススメ記事www.bananarian.netはてなブログをはじめよう!bananarianさんは、はてなブログを使っています。ブログを報告するスターをつけました引用をストックしましたストック一覧を見る閉じる引用するにはまずログインしてくださいログイン閉じる引用をストックできませんでした。再度お試しください閉じる限定公開記事のため引用できません。