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('http://biostat.mc.vanderbilt.edu/wiki/pub/Main/TatsukiRcode/TatsukiRcodeKMplot.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.labxaxis.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.labelbtyなども指定できる。

図だけで、下の数字はいらないや、という場合はsimple=TRUEで消せる。

他にも ggplot を使って似たようなことをした人は何人かいるけれど、ggplot なんか使わないで済むならそれにこしたことはない、と思う。

平成171819202122232425262728293031323310111210111213141516171819202122232425262728293031日 土曜日|統計学コメント(2)

コメント

○ ○

ブラヴォー!素晴らしいです

平成171819202122232425262728293031323310111210111213141516171819202122232425262728293031日 【日】 午前1時午前2時午前3時午前4時午前5時午前6時午前7時午前8時午前9時午前10時午前11時正午午後1時午後2時午後3時午後4時午後5時午後6時午後7時午後8時午後9時午後10時午後11時午前0時頃|Anonymous
○ ○

サンキュー!ありがとうございます

平成171819202122232425262728293031323310111210111213141516171819202122232425262728293031日 【木】 午前1時午前2時午前3時午前4時午前5時午前6時午前7時午前8時午前9時午前10時午前11時正午午後1時午後2時午後3時午後4時午後5時午後6時午後7時午後8時午後9時午後10時午後11時午前0時頃|びい

コメントの投稿

ブログ検索

カレンダー

JanFebMar AprMay JunJulAugSepOctNovDec 2005200620072008200920102011
- - - - 1 2 3
4 5 6 7 8 9 10
11 12 13 14 15 16 17
18 19 20 21 22 23 24
25 26 27 28 29 30 -

最近のコメント

履歴

2018年11月 (1)

2018年10月 (1)

2018年09月 (1)

2018年08月 (1)

2018年07月 (1)

2018年06月 (1)

2018年05月 (1)

2018年04月 (2)

2018年03月 (1)

2018年01月 (1)

2017年12月 (1)

2017年11月 (1)

2017年10月 (1)

2017年09月 (1)

2017年08月 (1)

2017年07月 (1)

2017年06月 (1)

2017年05月 (1)

2017年04月 (1)

2017年03月 (1)

2017年02月 (1)

2017年01月 (1)

2016年12月 (1)

2016年11月 (1)

2016年09月 (1)

2016年08月 (1)

2016年07月 (1)

2016年05月 (1)

2016年04月 (1)

2016年03月 (1)

2016年02月 (1)

2016年01月 (1)

2015年12月 (1)

2015年07月 (1)

2015年06月 (2)

2015年05月 (2)

2015年04月 (3)

2015年03月 (2)

2015年02月 (1)

2015年01月 (3)

2014年12月 (1)

2014年11月 (1)

2014年10月 (2)

2014年09月 (2)

2014年08月 (2)

2014年07月 (1)

2014年06月 (1)

2014年05月 (2)

2014年04月 (1)

2014年03月 (2)

2014年02月 (1)

2014年01月 (2)

2013年12月 (1)

2013年11月 (3)

2013年10月 (2)

2013年09月 (1)

2013年08月 (4)

2013年07月 (1)

2013年06月 (2)

2013年05月 (2)

2013年04月 (3)

2013年03月 (1)

2013年02月 (1)

2013年01月 (4)

2012年12月 (1)

2012年11月 (3)

2012年10月 (1)

2012年09月 (1)

2012年08月 (3)

2012年07月 (3)

2012年06月 (2)

2012年05月 (6)

2012年04月 (2)

2012年03月 (8)

2012年02月 (2)

2012年01月 (1)

2011年12月 (6)

2011年11月 (5)

2011年10月 (4)

2011年09月 (6)

2011年08月 (9)

2011年07月 (5)

2011年06月 (5)

2011年05月 (5)

2011年04月 (6)

2011年03月 (17)

2011年02月 (6)

2011年01月 (10)

2010年12月 (10)

2010年11月 (4)

2010年10月 (6)

2010年09月 (5)

2010年08月 (11)

2010年07月 (8)

2010年06月 (8)

2010年05月 (3)

2010年04月 (8)

2010年03月 (11)

2010年02月 (4)

2010年01月 (8)

2009年12月 (6)

2009年11月 (6)

2009年10月 (6)

2009年09月 (7)

2009年08月 (6)

2009年07月 (10)

2009年06月 (10)

2009年05月 (10)

2009年04月 (6)

2009年03月 (7)

2009年02月 (9)

2009年01月 (12)

2008年12月 (6)

2008年11月 (10)

2008年10月 (8)

2008年09月 (9)

2008年08月 (12)

2008年07月 (8)

2008年06月 (12)

2008年05月 (12)

2008年04月 (12)

2008年03月 (11)

2008年02月 (10)

2008年01月 (10)

2007年12月 (12)

2007年11月 (14)

2007年10月 (13)

2007年09月 (11)

2007年08月 (16)

2007年07月 (10)

2007年06月 (10)

2007年05月 (6)

2007年04月 (10)

2007年03月 (13)

2007年02月 (10)

2007年01月 (8)

2006年12月 (13)

2006年11月 (15)

2006年10月 (9)

2006年09月 (8)

2006年08月 (18)

2006年07月 (14)

2006年06月 (16)

2006年05月 (23)

2006年04月 (20)

2006年03月 (12)

2006年02月 (14)

2006年01月 (21)

2005年12月 (20)

2005年11月 (17)

2005年10月 (18)

2005年09月 (16)

2005年08月 (10)

2005年07月 (6)