[R program]時間依存性ROC曲線法

ある結果を予測するモデルの良さをチェックするときに有用なツールにROC曲線があります。

主にロジスティック回帰モデルの性能をチェックするために使われますが、生存時間解析でも使われることがあります。

ベースライン検査値と生存時間情報がある場合、ベースライン検査値での経時予測性能を時間依存性ROC曲線法で推定できます。

論文はこれです。

Heagerty, P.J., Lumley, T., Pepe, M. S. (2000) Time-dependent ROC Curves for Censored Survival Data and a Diagnostic Marker Biometrics, 56, 337 – 344

RではsurvivalROCパッケージで推定できます。

mayoというexampleデータがあり、ベースライン検査値はmayo4とmayo5というバイオマーカーです。

追跡期間は約12年。

それぞれのマーカーの、1年後と12年後のROCを描くとこんな感じ。

f:id:isseing333:20100820181836j:image

f:id:isseing333:20100820181835j:image

f:id:isseing333:20100820181834j:image

f:id:isseing333:20100820181833j:image

1年後の予測力はmayo4、mayo5共に同じくらいだけど、12年後の予測力は圧倒的にmayo5の方が良いことが分かる。

という感じで、モデル選択の指標として利用することができる。

ちなみにAUCを経時プロットするとこんな感じになる。

f:id:isseing333:20100820183058j:image

mayo4とmayo5の差が歴然としていることがわかる。

ちなみにこの折れ線の検定を同僚が研究していました（完成したのかな？）。


Rプログラムはこちら↓

library(survivalROC)

data(mayo)
ROC <- function(cut){
	cutoff <- cut
	Mayo4.2= survivalROC(Stime=mayo$time,status=mayo$censor, 
		marker = mayo$mayoscore4,predict.time = cutoff, method="KM")
	plot(Mayo4.2$FP, Mayo4.2$TP, type="l", xlim=c(0,1), ylim=c(0,1), 
		xlab=paste( "FP", "\n", "AUC = ",round(Mayo4.2$AUC,3)), ylab="TP", 
		main=paste("Mayoscore 4, Method = KM \n Year =", cutoff/365))
	abline(0,1)
	Mayo4.2$AUC
}

AUC4 <- rep(0, 12)
AUC4[1] <- ROC(365)
AUC4[2] <- ROC(365*2)
AUC4[3] <- ROC(365*3)
AUC4[4] <- ROC(365*4)
AUC4[5] <- ROC(365*5)
AUC4[6] <- ROC(365*6)
AUC4[7] <- ROC(365*7)
AUC4[8] <- ROC(365*8)
AUC4[9] <- ROC(365*9)
AUC4[10] <- ROC(365*10)
AUC4[11] <- ROC(365*11)
AUC4[12] <- ROC(365*12)


ROC2 <- function(cut){
	cutoff <- cut
	Mayo4.2= survivalROC(Stime=mayo$time,status=mayo$censor, 
		marker = mayo$mayoscore5,predict.time = cutoff, method="KM")
	plot(Mayo4.2$FP, Mayo4.2$TP, type="l", xlim=c(0,1), ylim=c(0,1), 
		xlab=paste( "FP", "\n", "AUC = ",round(Mayo4.2$AUC,3)), ylab="TP", 
		main=paste("Mayoscore 5, Method = KM \n Year =", cutoff/365))
	abline(0,1)
	Mayo4.2$AUC
}

AUC5 <- rep(0, 12)
AUC5[1] <- ROC2(365)
AUC5[2] <- ROC2(365*2)
AUC5[3] <- ROC2(365*3)
AUC5[4] <- ROC2(365*4)
AUC5[5] <- ROC2(365*5)
AUC5[6] <- ROC2(365*6)
AUC5[7] <- ROC2(365*7)
AUC5[8] <- ROC2(365*8)
AUC5[9] <- ROC2(365*9)
AUC5[10] <- ROC2(365*10)
AUC5[11] <- ROC2(365*11)
AUC5[12] <- ROC2(365*12)


AUC <- data.frame(Year = rep(1:12, 2), AUC= c(AUC4, AUC5) , mayo = rep(c("4", "5"), c(12, 12)))

library(ggplot2)
ggplot(AUC, aes(x=Year, y=AUC, group=mayo, col=mayo)) + geom_line() + geom_point() + ylim(0.5, 1)

ページTOPへ