統計計算のフリーソフトRは、生命医学研究のさまざまなところで使われています。
この記事では、Rを使った生存時間解析を紹介します。
生存時間解析とは
生存時間解析は、あるイベント (event) が起きるかどうかと、それが起きるまでの時間の関係を調べる統計手法のことです。
例えばその「イベント」が死亡であれば生存時間ということになりますが、「イベント」は他にも、例えば製品の故障でも、経営学の会社の倒産でもいいわけです。
そういったものをまとめて「生存時間解析」と呼んでいます。
研究の終了時点で、死亡に関するデータが入手できないとか、治療の中止や転院などでデータが欠落する場合も多く、これを打ち切り (censoring) があるといいます。
打ち切りに対応した統計手法が生存時間解析なのです。
Rには生存時間解析用のパッケージがいくつもありますが、そのうちで最もよく使われているものの1つはsurvivalパッケージです。
生存時間データ
今回はRにもとから入っているkidneyデータを使ってみましょう。これは,ポータブル透析装置の利用が各種腎臓病の患者の生存時間にどう影響するかを検討するためのデータです。
いくつかデータを見てみます。
library(survival)
data(kidney)
head(kidney)
id time status age sex disease frail
1 1 8 1 28 1 Other 2.3
2 1 16 1 28 1 Other 2.3
3 2 23 1 48 2 GN 1.9
4 2 13 0 48 2 GN 1.9
5 3 22 1 32 1 Other 1.2
6 3 28 1 32 1 Other 1.2
各項目の意味はこちらに書かれています。
生存分析の手法は、時間以外の共変量 (複数の要因、説明変数) が生存時間に影響を与えるパラメーターとして考えるかどうかや、生存時間の分布にある特定の確率分布を仮定するかどうかで、
セミノンパラメトリックモデル:共変量を導入する、分布を仮定しない
パラメトリックモデル:共変量を導入する、分布を導入する
と分類できます。このうちパラメトリックモデルはかなり強い制約となり、実際に使われることはかなり少ないので、この記事では残り2つを紹介します。
カプランマイヤー法
ノンパラメトリックモデルの代表がカプランマイヤー法です。surivivalパッケージには、ノンパラメトリック法による生存時間を当てはめる関数survfitがあります。
diseaseで群を分けて、カプランマイヤー法を当てはめてみます。
KM <- survfit (Surv (time, status) ~ disease, data = kidney)
カプランマイヤー法による生存曲線を書くにはこれをplot関数に渡すだけでできます。
plot(KM, xlab='Time', ylab='Survival Rate', col=c('red', 'green', 'blue', 'black'))
legend(300, 0.8, legend=c('Other', 'GN', 'AN', 'PKD'), lty=1, col=c('red', 'green', 'blue', 'black'))
群ごとの差を検定する打表的な手法がLog-Rank (ログランク) 検定で、これにはsurvdiff関数が使えます。
例えば、生存時間に男女差があるのかをログランク検定で調べてみます。
survdiff(Surv(time, status)~sex, data=kidney)
p値は0.004なので、生存時間に男女差がある (女性の方が長生き) と言えるでしょう。
コックス比例ハザードモデル
セミノンパラメトリックモデルの代表が、コックス比例ハザードモデルです。
病気の種類による影響をコックス比例ハザードモデルで調べてみます。
ここでは生存時間を性別と病気の種類で説明するモデルを作ってみました。
kidney.cox <- coxph (Surv(time, status)~sex + disease, data = kidney)
kidney.cox
構築したモデルによる生存時間の当てはめには、survfit関数を使うと便利です。
kidney.fit <- survfit (kidney.cox)
summary (kidney.fit)
survfitで推定された生存曲線と信頼区間をplot関数で図示することもできます。
plot (kidney.fit)
コックス比例ハザードモデルは、ハザード比が時間によらず一定であることを前提としています。survivalパッケージにはこの比例性を分析するcox.zph関数があります。
cox.zph (kidney.cox)
仮設が棄却されると、比例ハザードの仮定が満たされていない可能性があることを示唆しています。
関連図書
この記事に関連した内容を紹介している本はこちらです。
今日も【生命医学をハックする】 (@biomedicalhacks) をお読みいただきありがとうございました。