GGPLOT2找到生存数据帧绘制Kaplan-Meier曲线

Question

我试图用画一些Kaplan-Meier曲线时，空行GGPLOT2代码中的发现：https://github.com/kmiddleton/rexamples/blob/master/qplot_survival.R

我曾与在这个伟大的代码的好成绩不同的数据库然而，在这种情况下，它给了我下面的错误......我仿佛在我的数据框空行：

Error en if (nrow(layer_data) == 0) return() : argument is of length zero.

关于这种类型的错误前面的问题似乎并没有对我有用的，因为在我的情况下数据和功能的类型不一样。

我对使用R的统计分析颇为陌生，而且我没有编程背景，所以我认为这在我的数据中必须是一个“哑巴错误”，但我无法找到它的位置......它绝对似乎ggplot2找不到要绘制的行。请你能以任何方式帮助我，提供线索，建议等等。

这里是我的数据和使用的代码，按顺序，准备好了控制台 - 我试着用knitr脚本 - 。最后，我已为我的sessionInfo：

library(splines) 
library(survival) 
library(abind) 
library(ggplot2) 
library(grid) 

我创建了一个名为acbi30（真实数据）数据帧：

mort28day <- c(1,0,1,0,0,0,0,1,0,0,0,1,1,0,1,0,0,1,0,1,1,1,1,0,1,1,1,0,0,1) 
daysurv <- c(4,29,24,29,29,29,29,19,29,29,29,3,9,29,15,29,29,11,29,5,13,20,22,29,16,21,9,29,29,15) 
levo <- c(0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0) 
acbi30 <- data.frame(mort28day, daysurv, levo) 
save(acbi30, file="acbi30.rda") 
acbi30 

然后，我粘贴下面的命令来创建一个功能使用GGPLOT2：

t.Surv <- Surv(acbi30$daysurv, acbi30$mort28day) 
t.survfit <- survfit(t.Surv~1, data=acbi30) 


#define custom function to create a survival data.frame# 
createSurvivalFrame <- function(f.survfit){ 

#initialise frame variable# 
f.frame <- NULL 

#check if more then one strata# 
if(length(names(f.survfit$strata)) == 0){ 

#create data.frame with data from survfit# 
f.frame <- data.frame(time=f.survfit$time, n.risk=f.survfit$n.risk, n.event=f.survfit$n.event, n.censor = f.survfit 
$n.censor, surv=f.survfit$surv, upper=f.survfit$upper, lower=f.survfit$lower) 
#create first two rows (start at 1)# 
f.start <- data.frame(time=c(0, f.frame$time[1]), n.risk=c(f.survfit$n, f.survfit$n), n.event=c(0,0), 
n.censor=c(0,0), surv=c(1,1), upper=c(1,1), lower=c(1,1)) 
#add first row to dataset# 
f.frame <- rbind(f.start, f.frame) 
#remove temporary data# 
rm(f.start) 
} 
else { 
#create vector for strata identification# 
f.strata <- NULL 
for(f.i in 1:length(f.survfit$strata)){ 
#add vector for one strata according to number of rows of strata# 
f.strata <- c(f.strata, rep(names(f.survfit$strata)[f.i], f.survfit$strata[f.i])) 
} 
#create data.frame with data from survfit (create column for strata)# 
f.frame <- data.frame(time=f.survfit$time, n.risk=f.survfit$n.risk, n.event=f.survfit$n.event, n.censor = f.survfit 
$n.censor, surv=f.survfit$surv, upper=f.survfit$upper, lower=f.survfit$lower, strata=factor(f.strata)) 
#remove temporary data# 
rm(f.strata) 
#create first two rows (start at 1) for each strata# 
for(f.i in 1:length(f.survfit$strata)){ 
#take only subset for this strata from data# 
f.subset <- subset(f.frame, strata==names(f.survfit$strata)[f.i]) 
#create first two rows (time: 0, time of first event)# 
f.start <- data.frame(time=c(0, f.subset$time[1]), n.risk=rep(f.survfit[f.i]$n, 2), n.event=c(0,0), 
n.censor=c(0,0), surv=c(1,1), upper=c(1,1), lower=c(1,1), strata=rep(names(f.survfit$strata)[f.i], 
2)) 
#add first two rows to dataset# 
f.frame <- rbind(f.start, f.frame) 
#remove temporary data# 
rm(f.start, f.subset) 
} 
#reorder data# 
f.frame <- f.frame[order(f.frame$strata, f.frame$time), ] 
#rename row.names# 
rownames(f.frame) <- NULL 
} 
#return frame# 
return(f.frame) 
} 


#define custom function to draw kaplan-meier curve with ggplot# 
qplot_survival <- function(f.frame, f.CI="default", f.shape=3){ 
#use different plotting commands dependig whether or not strata's are given# 
if("strata" %in% names(f.frame) == FALSE){ 
#confidence intervals are drawn if not specified otherwise# 
if(f.CI=="default" | f.CI==TRUE){ 
#create plot with 4 layers (first 3 layers only events, last layer only censored)# 
#hint: censoring data for multiple censoring events at timepoint are overplotted# 
#(unlike in plot.survfit in survival package)# 
ggplot(data=f.frame) + geom_step(aes(x=time, y=surv), direction="hv") + geom_step(aes(x=time, 
y=upper), directions="hv", linetype=2) + geom_step(aes(x=time,y=lower), direction="hv", linetype=2) + 
geom_point(data=subset(f.frame, n.censor==1), aes(x=time, y=surv), shape=f.shape) 
} 
else { 
#create plot without confidence intervals# 
ggplot(data=f.frame) + geom_step(aes(x=time, y=surv), direction="hv") + 
geom_point(data=subset(f.frame, n.censor==1), aes(x=time, y=surv), shape=f.shape) 
} 
} 
else { 
if(f.CI=="default" | f.CI==FALSE){ 
#without CI# 
ggplot(data=f.frame, aes(group=strata, colour=strata)) + geom_step(aes(x=time, y=surv), 
direction="hv") + geom_point(data=subset(f.frame, n.censor==1), aes(x=time, y=surv), shape=f.shape) 
} 
else { 
#with CI (hint: use alpha for CI)# 
ggplot(data=f.frame, aes(colour=strata, group=strata)) + geom_step(aes(x=time, y=surv), 
direction="hv") + geom_step(aes(x=time, y=upper), directions="hv", linetype=2, alpha=0.5) + 
geom_step(aes(x=time,y=lower), direction="hv", linetype=2, alpha=0.5) + 
geom_point(data=subset(f.frame, n.censor==1), aes(x=time, y=surv), shape=f.shape) 
} 
} 
} 

绘制全球生存曲线（95％CI）：

它不给任何错误：

# Kaplan-Meier plot, global survival (with CI) 
t.survfit <- survfit(t.Surv~1, data=acbi30) 
t.survframe <- createSurvivalFrame(t.survfit) 
t.survfit 
qplot_survival(t.survframe, TRUE, 20) 

绘制分层生存曲线：

给出了错误上面提到：

# Kaplan-Meier plot, stratified survival 
t.survfit2 <- survfit(t.Surv~levo, data=acbi30) 
t.survframe2 <- createSurvivalFrame(t.survfit2) 
t.survfit2 
qplot_survival(t.survframe2, TRUE, 20) 

绘图而不GGPLOT2结果：

t.survframe2的结构对我来说似乎没问题，没有任何空行，所以它必定是qplot_survival在t.survframe2中读取我的数据的问题。做一个简单的阴谋不会返回任何错误：

t.survframe2 
plot(t.survfit2) 

问题与我的数据框在哪里？该功能创建的工作以及与其他数据集，但与这一个...

谢谢你在前进，

Mareviv

会议信息：

sessionInfo() 

[R版2.15 。2（2012年10月26日） 平台：I386-W64-的mingw32/I386（32位）

locale: 

[1] LC_COLLATE=Spanish_Spain.1252 LC_CTYPE=Spanish_Spain.1252 
[3] LC_MONETARY=Spanish_Spain.1252 LC_NUMERIC=C     
[5] LC_TIME=Spanish_Spain.1252  

attached base packages: 
[1] grid  splines stats  graphics grDevices utils  datasets 
[8] methods base  

other attached packages: 
[1] ggplot2_0.9.3 abind_1.4-0  survival_2.36-14 knitr_0.8  

loaded via a namespace (and not attached): 
[1] colorspace_1.1-1 dichromat_1.2-4 digest_0.5.2  
[4] evaluate_0.4.2  formatR_0.7  gtable_0.1.2  
[7] labeling_0.1  MASS_7.3-22  munsell_0.4  
[10] plyr_1.8   proto_0.3-9.2  RColorBrewer_1.0-5 
[13] reshape2_1.2.1  scales_0.2.3  stringr_0.6.1  
[16] tools_2.15.2  

Answer 1

我做你的qplot_survival()的作用有点整容手术。主要问题似乎是data =参数geom_point中的子集条件;在t.survframe和t.survframe2中，n.censor的表得出值0,3和12.通过将子集条件改变为n.censor > 0，我设法在所有情况下得到一个图。我也没有看到f.CI = "default"的点，所以我将默认设置为TRUE，并相应地修改if条件。

qplot_survival <- function(f.frame, f.CI= TRUE, f.shape=3) 
{ 
# use different plotting commands depending whether 
# or not strata are given# 
if(!("strata" %in% names(f.frame))) 
{ 
    #confidence intervals are drawn if not specified otherwise# 
    if(isTRUE(f.CI)) 
    { 
     # create plot with 4 layers (first 3 layers only events, 
     # last layer only censored)# 
     # hint: censoring data for multiple censoring events at 
     # timepoint are overplotted# 
     # (unlike in plot.survfit in survival package)# 
    ggplot(data=f.frame) + 
     geom_step(aes(x=time, y=surv), direction="hv") + 
     geom_step(aes(x=time, y=upper), direction ="hv", linetype=2) + 
     geom_step(aes(x=time,y=lower), direction="hv", linetype=2) + 
     geom_point(data=subset(f.frame, n.censor > 0), 
       aes(x=time, y=surv), shape=f.shape) 
    } else { 
    #create plot without confidence intervals# 
    ggplot(data=f.frame) + 
     geom_step(aes(x=time, y=surv), direction="hv") + 
     geom_point(data=subset(f.frame, n.censor > 0), 
       aes(x=time, y=surv), shape=f.shape) 
      } 
} else { 
    if(!(isTRUE(f.CI))){ 
#without CI# 
    ggplot(data=f.frame, aes(group=strata, colour=strata)) + 
    geom_step(aes(x=time, y=surv), direction="hv") + 
    geom_point(data=subset(f.frame, n.censor > 0), 
       aes(x=time, y=surv), shape=f.shape) 
} else { 

#with CI (hint: use alpha for CI)# 
    ggplot(data=f.frame, aes(x = time, colour=strata, group=strata)) + 
     geom_step(aes(y=surv), direction="hv") + 
     geom_step(aes(y=upper), direction="hv", 
        linetype=2, alpha=0.5) + 
     geom_step(aes(y=lower), direction="hv", 
        linetype=2, alpha=0.5) + 
     geom_point(data=subset(f.frame, n.censor > 0), 
       aes(y=surv), shape=f.shape) 
     } 
    } 
} 

下面的情节做了这些改变后，所有的工作对我来说：

qplot_survival(t.survframe2, TRUE, 20) 
qplot_survival(t.survframe2, FALSE, 20) 
qplot_survival(t.survframe, TRUE, 20) 
qplot_survival(t.survframe, FALSE, 20) 

一对夫妇的意见：

1. 在函数内部子集可能是危险的，因为有时，在本情况下，满足条件返回一个零行数据帧。我会考虑geom_point()层是否真的有必要。
2. 在几个地方，您在geom_step()呼叫中有directions = "hv"。这个论点不是多元化的，并且已经在上面改变了。
3. 这可以有效地完成多一点我想，而是从一个survfit对象提取所关注列的一种方式，说t.survfit，是这样的：

（展开谱曲时阶层都存在）

comps <- c(2:6, 8, 10); 
t.fit <- as.data.frame(do.call(cbind, lapply(comps, function(j) t.survfit[[j]]))) 
names(t.fit) <- names(t.survfit)[comps] 

Answer 2

这里是另一个版本时出现在你的数据不删失分（@丹尼斯的版本，在这种情况下仍然无法）也占的情况。这可以变得更高效，可能是通过创建一个变量来存储整个数据框中的检查点数量，然后重新使用它，而不是像我在每种情况下那样进行测试。

# define custom function to draw kaplan-meier curve with ggplot 
qplot_survival <- function(f.frame, f.CI="default", f.shape=3){ 
    # use different plotting commands dependig whether or not strata's are given 
    if("strata" %in% colnames(f.frame) == FALSE){ 
    # confidence intervals are drawn if not specified otherwise 
    if(f.CI=="default" | f.CI==TRUE){ 
     # create plot with 4 layers (first 3 layers only events, last layer only censored) 
     # hint: censoring data for multiple censoring events at timepoint are overplotted 



     # (unlike in plot.survfit in survival package) 
     p <- ggplot(data=f.frame) + geom_step(aes(x=time, y=surv), direction="hv") + geom_step(aes(x=time, 
                          y=upper), directions="hv", linetype=2) + geom_step(aes(x=time,y=lower), direction="hv", linetype=2) 
     if(nrow(subset(f.frame, n.censor > 0)) > 0){ 
     p+geom_point(data=subset(f.frame, n.censor > 0), aes(x=time, y=surv), shape=f.shape) 
     }else{ 
     p 
     } 
    } 
    else { 
     # create plot without confidence intervalls 
     p <- ggplot(data=f.frame) + geom_step(aes(x=time, y=surv), direction="hv") 
     if(nrow(subset(f.frame, n.censor > 0)) > 0){ 
     p + geom_point(data=subset(f.frame, n.censor > 0), aes(x=time, y=surv), shape=f.shape) 
     }else{ 
     p 
     } 
    } 
    } 
    else { 
    if(f.CI=="default" | f.CI==FALSE){ 
     # without CI 
     p <- ggplot(data=f.frame, aes(group=strata, colour=strata)) + geom_step(aes(x=time, y=surv), 
                     direction="hv") 
     if(nrow(subset(f.frame, n.censor > 0)) > 0){ 
     p +geom_point(data=subset(f.frame, n.censor > 0), aes(x=time, y=surv), shape=f.shape) 
     }else{ 
     p 
     } 
    } 
    else { 
     # with CI (hint: use alpha for CI) 
     p <- ggplot(data=f.frame, aes(colour=strata, group=strata)) + geom_step(aes(x=time, y=surv), 
                     direction="hv") + geom_step(aes(x=time, y=upper), directions="hv", linetype=2, alpha=0.5) + 
     geom_step(aes(x=time,y=lower), direction="hv", linetype=2, alpha=0.5) 
     if(nrow(subset(f.frame, n.censor > 0)) > 0){ 
     p + geom_point(data=subset(f.frame, n.censor > 0), aes(x=time, y=surv), shape=f.shape) 
     }else{ 
     p 
     } 
    } 
    } 
}