Sunday, Monday, Tuesday, Wednesday, Thursday, Friday, Saturday, 12-17-2011
Kaplan-Meier estimates をグラフにするプログラムを書いた。Rだ。
仕事上、何度となくグラフにしていて、その都度せくせくコードを書いていたのだけど、今回、昔書いたプログラムを大幅に改善した。生存曲線を書いて、下のマージンに numbers at risk を書く。
例。セットアップから。
library(survival)
kma <- survfit( Surv(time, status) ~ rx + adhere, data=colon )
この kmaの中にグラフにするべき数字が隠されている。
僕が書いた function はkmplot()という名前で
ここ にあるので、
source('https://biostat.app.vumc.org/wiki/pub/Main/TatsukiRcode/RFunctions1.R')
で読み込める。
例の続き:
シンプルにkmplot(kma)で、それなりのグラフは得られるんだけど、細かにいろいろと設定することが可能。
kmplot(kma, mark='',
xaxis.at=c(0,.5,1:9)*365, xaxis.lab=c(0,.5,1:9),
lty.surv=c(1,2), lwd.surv=1, col.surv=c(1,1,2,2,4,4),
col.ci=0,
group.names=c('Obs ','Obs tumor adh','Lev','Lev tumor adh',
'Lev+5FU ','Lev+5FU tumor adh'),
group.order=c(5,3,1,6,4,2),
extra.left.margin=6, label.n.at.risk=FALSE, draw.lines=TRUE,
cex.axis=.8, xlab='Years', ylab='Survival Probability',
grid=TRUE, lty.grid=1, lwd.grid=1, col.grid=grey(.9),
legend=TRUE, loc.legend='bottomleft',
cex.lab=.8, xaxs='r', bty='L', las=1, tcl=-.2
)
これでこの図が描ける。
唯一絶対必要なインプットはsurvfit()で作ったオブジェクト(この例ではkma)だけ。それ以降のインプットは・・・
まずmarkは、打ち切りデータをどう表示するか。特に気にしなくて良い。いらない場合はmark=''。
次にxaxis.at。これはn at riskを計算して表示する x軸の値(つまり時)を指定する。何もしていないと適当な所になる。xaxis.labはxaxis.atで表示するラベル。例えば、時間は日にちで指定されているけれど、年を表示したいという時は、xaxis.at=(0:5)*365, xaxis.lab=0:5とすれば、0日、365日、730日のところに、0,1,2と表示することが出来る。
lty.surv, lwd.surv, col.survは、それぞれ線の種類、太さ、色だ。グループが4つあれば、4つ指定することもできるけれど、4つ未満の場合は繰り返して使われる。Rではいつものことだ。lty.surv=1:2, lty.lwd=1, lty.col=1:4とか。
lty.ci, lwd.ci, col.ciは confidence interval を表示する時に使う。表示したくなければ、lty.ci=0あるいはcol.ci=0とする。
group.namesは別に指定しなくていいんだけど、その場合、kmaオブジェクトで使われているものをそのまま引用する。それだと、たいていの場合うざい。
group.names=c('trt1','trt2','trt3','trt4')とかに書き換えてしまえば良い。 group.orderを指定しないと、ま、そのままグラフになるのだけど、n at riskの表示と、上のグラフがだいたい同じ順序だと見やすい。そこで、表示の順番を変えられる。例えば、group.order=c(2,4,1,3)とすると2番目のグループを1番上に書いて、その次が4番目のグループ、その下に、1番目のグループ、3番目のグループとなる。ま、例を見てください。
group.namesを短く指定してもそれでも見切れる場合はextra.left.marginを指定する。label.n.at.risk=TRUEとすると"Numbers at risk"と表示する。必要ないかな。
grid=TRUEで座標を表示する。細かな設定はcol.gridなどで。
legend=TRUEにすると、グループの名前と線の色、種類なんかを表示できる。場所はloc.legendで変えられる。
その他のインプットはplot()に行くので、cex.label、btyなども指定できる。
図だけで、下の数字はいらないや、という場合はsimple=TRUEで消せる。
他にも ggplot を使って似たようなことをした人は何人かいるけれど、ggplot なんか使わないで済むならそれにこしたことはない、と思う。
平成17
平成18
平成19
平成20
平成21
平成22
平成23
平成24
平成25
平成26
平成27
平成28
平成29
平成30
令和元
令和2
令和3
令和4
令和5
令和6
令和7
令和8年
123456789101112月
12345678910111213141516171819202122232425262728293031日 土曜日|
統計学|
コメント(4)
○ ○
ブラヴォー!素晴らしいです
平成17
平成18
平成19
平成20
平成21
平成22
平成23
平成24
平成25
平成26
平成27
平成28
平成29
平成30
令和元
令和2
令和3
令和4
令和5
令和6
令和7
令和8
令和9年
123456789101112月12345678910111213141516171819202122232425262728293031日 【日】 午前1時午前2時午前3時午前4時午前5時午前6時午前7時午前8時午前9時午前10時午前11時正午午後1時午後2時午後3時午後4時午後5時午後6時午後7時午後8時午後9時午後10時午後11時午前0時頃|Anonymous
○ ○
サンキュー!ありがとうございます
平成17
平成18
平成19
平成20
平成21
平成22
平成23
平成24
平成25
平成26
平成27
平成28
平成29
平成30
令和元
令和2
令和3
令和4
令和5
令和6
令和7
令和8
令和9年
123456789101112月12345678910111213141516171819202122232425262728293031日 【木】 午前1時午前2時午前3時午前4時午前5時午前6時午前7時午前8時午前9時午前10時午前11時正午午後1時午後2時午後3時午後4時午後5時午後6時午後7時午後8時午後9時午後10時午後11時午前0時頃|びい
○ Rコード ○
No. at riskをggplot2を使わずに作成する方法を探しており、このページにたどり着きました。
コードを拝見したいと思ったのですが、リンク切れのようでした。
可能でしたらコードをご教示いただければ幸いです。
平成17
平成18
平成19
平成20
平成21
平成22
平成23
平成24
平成25
平成26
平成27
平成28
平成29
平成30
令和元
令和2
令和3
令和4
令和5
令和6
令和7
令和8
令和9年
123456789101112月12345678910111213141516171819202122232425262728293031日 【日】 午前1時午前2時午前3時午前4時午前5時午前6時午前7時午前8時午前9時午前10時午前11時正午午後1時午後2時午後3時午後4時午後5時午後6時午後7時午後8時午後9時午後10時午後11時午前0時頃|はしもと
○ リンク切れ訂正 ○
リンクを訂正しておきました。リンク先の RFunctions1.R の中ほどに隠れています。
9年前に書いたコードですが、最終アップデートは2018年です。ggplot のもいいんですが、自分でゴテゴテ書くと細かくカスタマイズできる利点はあります。
平成17
平成18
平成19
平成20
平成21
平成22
平成23
平成24
平成25
平成26
平成27
平成28
平成29
平成30
令和元
令和2
令和3
令和4
令和5
令和6
令和7
令和8
令和9年
123456789101112月12345678910111213141516171819202122232425262728293031日 【日】 午前1時午前2時午前3時午前4時午前5時午前6時午前7時午前8時午前9時午前10時午前11時正午午後1時午後2時午後3時午後4時午後5時午後6時午後7時午後8時午後9時午後10時午後11時午前0時頃|びい