9月2日（火）午後 - 2014年度マイコースプログラム

生存曲線　3

まず、カプラン・マイヤー法について簡単にまとめる。

カプラン・マイヤー法とは、時刻 ${t}$ までの生存率を

${\widehat{S}(t) = \displaystyle \prod_{t_i \le t}^{}(1-\frac{\displaystyle d_i}{\displaystyle r_i})}$

と推定するものである。

（ちなみに上の式書くのかなり苦労した。）

ここで、式の中の ${d_i}$ は、時点 ${t_i}$ における死亡数、 ${r_i}$ は時点 ${t_i}$ の直前まで死亡しうる可能性のある観察対象者の数である。

これにより、脱落者や追跡不可能になってしまった人がいても、それが非常に多くない限りは、生存率を近似して求めることができる。

本当はsurvfit()なる関数でこの計算ができるのだが、今回はこれを自作してみた。

library(survival);library(MASS)
data(gehan) 
dim(gehan)

> dim(gehan)
[1] 42  4

gehan[1:42,]

> gehan[1:42,]
   pair time cens   treat
1     1    1    1 control
2     1   10    1    6-MP
3     2   22    1 control
4     2    7    1    6-MP
5     3    3    1 control
6     3   32    0    6-MP
7     4   12    1 control
8     4   23    1    6-MP
9     5    8    1 control
10    5   22    1    6-MP
11    6   17    1 control
12    6    6    1    6-MP
（以下略）

c1<-seq(1,41,2)
c2<-seq(2,42,2)

mygehan.Ctrl<-gehan[c1,2:3]
mygehan.6MP<-gehan[c2,2:3]
#データを使いやすく加工

rt<-function(d,t){
length(d[,1][d[,1]>=t])
}#rの値を出す関数

di<-function(d,t){
cnt<-0
for(i in 1:nrow(d)){
if(d[i,1]==t&&d[i,2]==1)cnt<-cnt+1
}
return(cnt)
}#diを出す関数

S<-function(d,t){
x<-1
for(i in 0:t){
x<-x*(1-di(d,i)/rt(d,i))
}
return(x)
}#カプラン・マイヤー法に適用

Sx<-function(d){
x<-0:max(d[,1])
return(x)
}
Sy<-function(d){
y<-c()
for(i in 0:max(d[,1])){
y<-c(y,S(d,i))
}
return(y)
}

 
plot(Sx(mygehan.Ctrl),Sy(mygehan.Ctrl),xlim=c(0,35), ylim=c(0,1), 
ann=F,type="s", lty=2) 
   
par(new=T)                                          
plot(Sx(mygehan.6MP),Sy(mygehan.6MP),xlim=c(0,35), ylim=c(0,1),
xlab="日数" ,type="s", ylab="生存率" , main="生存曲線") 

legend(20,0.95,c("6-MP","Control"),lty=1:2)

出来たグラフが下の通り。

f:id:kouri_don:20140902172753p:plain

ちなみにsurvfit()を使って作ったものは下の通り。

library(survival);library(MASS)
data(gehan) 

ge.sf<-survfit(Surv(time,cens)~treat, data=gehan)

plot(ge.sf,lty=1:2)

f:id:kouri_don:20140902173808p:plain

打ち切りの箇所に「＋」がつけられている以外は全く同じである。

過程はともかく、出力している内容は同じだと思われる。

実は生存率の近似方法には他にも何種類かあるらしいが、とりあえず今これ以上時間をかけるほどではないかと思うので、生存曲線については一旦これで満足したことにする。

グラフを操作するいい練習になったと思う。

以上

生存曲線 3

生存曲線　3