SelectionImpactPLOT_gg <-
function(x, ci, ci.bounds, get.F, fixed.values,conf.bands, rho, xlab, ylab, xlim, ylim, main){
risk.t0 <- x$derived.data$fittedrisk.t0
risk.t1 <- x$derived.data$fittedrisk.t1
trt.effect <- x$derived.data$trt.effect
if(x$model.fit$family$family == "time-to-event"){
event = rep(0, nrow(x$derived.data))
event.name = x$treatment.name
#setting event.name to treatment.name: this doesn't matter since we use model
#based estimates of event rates to get the marginal treatment effect for the plots.
}else{
event <- x$derived.data[[as.character(x$formula[[2]])]]
event.name = as.character(x$formula[[2]])
}
trt <- x$derived.data[[x$treatment.name]]
n = length(trt)
F.D <- get.F(trt.effect, event, trt, rho = rho)*100
theta.curve <- EventRateVec(risk.t0, risk.t1, F.D, rho, event, trt)
lty = 1 #define the lty for the main curve
mydata = data.frame(theta.curve, F.D, lty )
mydata = mydata[with(mydata, order(F.D)),]
## need the mean risk given t=0 and t=1. We need to account for subsampling in this.
allMeasures <- x$functions$get.summary.measures( x$derived.data, event.name = event.name,
treatment.name = x$treatment.name,
rho = rho,
d = x$model.fit$thresh)
avglines <- cbind(allMeasures$ER.trt0.mod, sort(F.D), 4)
avglines <- rbind(avglines,
cbind(allMeasures$ER.trt1.mod, sort(F.D), 3))
avglines = data.frame(avglines); names(avglines) = names(mydata)
mydata <- rbind(mydata, avglines)
cen <- mean(c(min(theta.curve, na.rm=TRUE), max(theta.curve, na.rm=TRUE)))
ran <- max(theta.curve, na.rm=TRUE) - min(theta.curve, na.rm=TRUE)
if(substr(ci, 1,1) =="v"){ cen <- mean(c(min(c(theta.curve, ci.bounds), na.rm=TRUE), max(c(theta.curve, ci.bounds), na.rm=TRUE))); ran <- max(c(theta.curve, ci.bounds), na.rm=TRUE) - min(c(theta.curve, ci.bounds), na.rm=TRUE)}
ran <- ran*1.1
mylim <- c(cen-ran/2, cen+ran/2)
if(is.null(xlab)) xlab <- "% population recommended treatment"
if(is.null(ylab)) ylab <- "Rate of outcome under treatment policy"
if(is.null(xlim)) xlim <- c(0,100)
if(is.null(ylim)) ylim <- mylim
if(is.null(main)) main <- "Selection impact curve"
p <- ggplot()
#add x/y labels and main
p <- p + xlab(xlab) + ylab(ylab) + ylim(ylim[1], ylim[2])+ ggtitle(main)
# p <- p + theme( text = element_text(size=18)) #,
p <- p + scale_x_continuous(limits = xlim)
if(!is.null(ci.bounds)){
ci.bounds <- matrix(ci.bounds, ncol=length(fixed.values), nrow = 2)
if(substr(ci, 1,1)=="v"){
index.fix <- (fixed.values<= max(F.D) & fixed.values >= min(F.D))
width = 5
}else{
width = .05
index.fix <- (fixed.values<= max(theta.curve) & fixed.values >= min(theta.curve))
}
p <- shade_gg(p, ci.bounds[,index.fix], fixed.values[index.fix], type = substr(ci, 1, 1), bands = conf.bands, lty, width = width)
}
p <- p+geom_step(data = mydata[(1:(n)),], aes(x = F.D, y = theta.curve), direction = "vh")
p <- p+geom_line(data = mydata[-c(1:(n)),], aes(x = F.D, y = theta.curve, linetype = factor(lty)), size = 0.5)
p <- p+ scale_linetype_manual(name = "Event Rate", breaks = c( "3", "4"), values = c( 3, 4), labels = c("treat all", "treat none"))
print(p)
return(list(p=p))
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.