カプランマイヤー曲線をggplot系の図で作りたい【R】【ggsurvplot】【 survminerパッケージ】
Rでカプランマイヤー曲線の図を作る機会があった。
せっかくなので、デフォルトのsurvivalパッケージだけではなく, survminerパッケージのggsurvplotというggplot系コードを使って鮮やかにしてみたので、備忘録的に
今回使うデータとは異なるデータですが、こんな感じのが作れたりします
使うデータと分析
データはMASSパッケージのgehanを使います。
このページで示されてるのと同じ分析をggsurvplotでやってみます
https://www1.doshisha.ac.jp/~mjin/R/Chap_36/36.html
通常のsurvival パッケージでのカプランマイヤー
#パッケージの読み込み
library(survival) #生存時間の関数を作るために呼びだし
library(survminer) # ggsurvplotのパッケージ
library(MASS) # gehanデータを使うために呼び出し
data(gehan) #データの呼び出し
#変数リスト
#pair: Pair index, time: Remission time (weeks), cens: 0=censored, treat: Treatment: control or 6-MP
ge.sf <- survival::survfit(Surv(time,cens)~treat,data=gehan) # モデル
plot(ge.sf,lty=1:2) #プロット
legend("topright" ,c("6-MP投与群","対照群"), lty=c(1,2))
ggsurvplotでカプランマイヤーやってみる
シンプルなカプランマイヤー曲線
ggsurvplot(ge.sf,dat=gehan, legend.title = "legend",legend=c(.8,.8), size=1,,xlab="xlab", ylab="ylab",censor=T, ggtheme = theme_gray())
# sizeは曲線の太さ, ggthemeでは背景を様々に変えることができます。
ちょっとオプションを足すと・・・
#信頼区間やLog-rank testの結果を示す pval や conf.int をTRUEに
ggsurvplot(ge.sf,dat=gehan, legend.title = "legend",legend=c(.8,.8), size=1,
pval =T ,conf.int = T, xlab="xlab", ylab="ylab",censor=T, ggtheme = theme_gray())
図の設定をいじりたい
#ggplot系ですので、どんどん後ろにスクリプトを足していくことで、いろいろな設定を変えていけます。
#たとえば さきほどのシンプルな図を temp に代入し,
temp <- ggsurvplot(ge.sf,dat=gehan, legend.title = "legend",legend=c(.8,.8), size=1,,xlab="xlab", ylab="ylab",censor=T, ggtheme = theme_gray())
temp$plot + scale_x_continuous(breaks=seq(0,40,4))
#と付け足すと
x軸を4weekごとの目盛りを入れることができます。
ちょっとした注意点として,tempに代入後,temp$plot の続きに足していく要素を書かないとエラーが出てきてしまいます。
Risk Table?
あと社会学系では見たことないのですが、リスクテーブルなるものも表示できるようです。
ggsurvplot(ge.sf,dat=gehan, legend.title = "legend",legend=c(.8,.8), size=1,,xlab="xlab", ylab="ylab",censor=T, ggtheme = theme_gray(),risk.table=T)
おまけ
最初の図はこんなコードでつくりました
fit <- survfit(Surv(month,event)~age311_cat,data=dat)
temp <- ggsurvplot(fit, data = dat, size=2,legend.title = "世帯主年齢", legend.labs=c("20未満","20代","30代","40代","50代","60代","70代","80歳以上"),legend=c(.2,.85),conf.int = F, ylab="すまい再建率", xlab="時間", censor=F,pval = F, ggtheme = theme_gray(), fun="event")
temp <- temp$plot + theme(axis.text=element_text(size=30))
temp <- temp + theme(axis.title.x = element_text(size=40),axis.title.y = element_text(size=40))
temp <- temp + scale_x_continuous(breaks=c(0,12,24,36,48,60,72), labels=c("2011/3", "2012/3","2013/3","2014/3","2015/3","2016/3","2017/3"))
temp <- temp+ylim(0,1)
temp + theme(legend.title = element_text(size=38),legend.text = element_text(size=38))
# fun= "event" で y軸を逆転できます。 あとはこのデータは,2011年3月観察開始で,仮設住宅を退居するまでの時間を分析する図なので,x軸は12ヶ月後ごとに設定,そのラベルとして年を示しています。