%\VignetteEngine{R.rsp::asis} %\VignetteEncoding{UTF-8} %\VignetteIndexEntry{seqDesign's Report Template} %\VignettePackage{seqDesign} \documentclass[11pt]{article} \usepackage{amsmath, amsthm, amssymb} \usepackage[cp1250]{inputenc} \usepackage[english]{babel} \usepackage{amsfonts} \usepackage{multirow} \usepackage[round,authoryear]{natbib} \usepackage{lscape} \usepackage{float} \usepackage{array} \usepackage{subfig} \usepackage{epsfig, psfrag, graphicx, lscape} \textwidth=6.5in \textheight=8.4in \setlength{\topmargin}{0.21truein} \hoffset-18mm \usepackage{parskip} \newcommand\T{\rule{0pt}{2.4ex}} \newcommand\TT{\rule{0pt}{2.9ex}} \begin{document} \begin{center} {\LARGE Operating Characteristics of the Specified Trial Design} \end{center} \vspace*{5mm} <>= rm(list=ls(all=TRUE)) library(xtable) # THE BELOW USER-SPECIFIED VALUES ARE AUTOMATICALLY APPLIED TO THE ENTIRE REPORT N.pla <- 1900 # number of subjects in the control arm N.vax <- 1700 # number of subjects in each treatment arm N.vax.arms <- 2 # number of treatment arms infRate <- 0.04 # annual incidence in the control arm estimand <- "cuminc" # VE estimand definition laggedMonitoring <- TRUE # should post-6 months non-efficacy monitoring be applied? lowerVEuncPower <- 0 # a numeric vector if unconditional powers to reject null hypotheses H0: VE(0-stage1) <= lowerVEuncPower[i]*100% should be reported; NULL otherwise srcDir <- "./" # a directory storing output from simTrial, monitorTrial, censTrial, rankTrial, and VEpowerPP # END OF USER-SPECIFIED VALUES @ <>= # 'aveVE.vector' determines what results are reported in each row aveVE.vector <- c(-2, -1.5, -1, -0.5, seq(0, 0.8, by=0.1)) # generate a shell for the results res <- as.data.frame(matrix(0, ncol=6, nrow=length(aveVE.vector))) colnames(res) <- c("aveVE", "aveHR", "harm", "noneffInterim", "noneffFinal", "higheff") res[,1] <- c(rep("--",sum(aveVE.vector<0)), paste0(aveVE.vector[aveVE.vector>=0]*100,"%")) res[,2] <- 1 - aveVE.vector if (!is.null(lowerVEuncPower)){ designPower <- as.data.frame(matrix(0, nrow=length(aveVE.vector), ncol=length(lowerVEuncPower))) } colnames(designPower) <- paste0("design",1:length(lowerVEuncPower)) # extract the estimated probabilities of each trial outcome and, if desired, unconditional power to reject H0 defined by 'lowerVEuncPower' from the output of 'monitorTrial' row <- 0 for (aveVE in aveVE.vector){ row <- row + 1 filename <- paste0("monitorTrial_nPlac=", N.pla, "_nVacc=", N.vax, "_aveVE=", aveVE, "_infRate=",infRate,"_",estimand,".RData") load(file.path(srcDir, filename)) bounds <- sapply(out, function(trial) trial[[1]]$boundHit) bounds.freq <- table(bounds, useNA="ifany")/length(out) res[row,3:6] <- bounds.freq[c("Harm","NonEffInterim","NonEffFinal","HighEff")] res[row,] <- ifelse(is.na(res[row,]), 0, res[row,]) if (!is.null(lowerVEuncPower)){ altDetectedMatrix <- do.call("rbind",lapply(out, function(oneTrial){ oneTrial[[1]]$altDetected })) if (is.null(altDetectedMatrix)){ designPower[row,] <- numeric(length(lowerVEuncPower)) } else { designPower[row,] <- colSums(altDetectedMatrix, na.rm=TRUE)/length(out) } } } res$noneff <- res$noneffInterim + res$noneffFinal res <- subset(res, select=c("aveVE","aveHR","harm","noneff","higheff")) if (!is.null(lowerVEuncPower)){ res <- data.frame(res, designPower) } res[,3:NCOL(res)] <- res[,3:NCOL(res)]*100 # prints a table with probabilities of reaching each trial monitoring outcome and unconditional power to reject H0 print(xtable(res, align=rep(c("c","|c","c"), c(6,ifelse(!is.null(lowerVEuncPower),1,0),ifelse(!is.null(lowerVEuncPower),length(lowerVEuncPower)-1,0))), digits=1, caption=paste0("Probabilities ($\\times 100$) of reaching each possible trial monitoring outcome",ifelse(!is.null(lowerVEuncPower)," and unconditional power ($\\times 100$) to reject the specified null hypothesis","")," for a 2-arm study design with ",N.pla," placebo and ",N.vax," vaccine recipients") ), scalebox=ifelse(!is.null(lowerVEuncPower),0.85,1), caption.placement="top", include.rownames=FALSE, include.colnames=FALSE, hline.after=4, add.to.row=list(pos=list(0,13), command=c(paste0("\\hline \\hline & & \\multicolumn{2}{c}{Weed Out at Interim Analysis} & ",ifelse(!is.null(lowerVEuncPower),rep("& ",length(lowerVEuncPower)),""),"\\\\\n \\cline{3-4} Average & Average & Potential Harm & Non-Efficacy & High Efficacy ",ifelse(!is.null(lowerVEuncPower),paste0("& ",paste0("\\multicolumn{",length(lowerVEuncPower),"}{c}{Unconditional Power} ")),""),"\\\\\n VE(0-18)$^{\\ast}$ & HR(0-18) & VE(0-18)$<$0\\% & VE(0-18)$<$40\\% & VE(0-18)$>$60\\% ",ifelse(!is.null(lowerVEuncPower),paste0("& ",paste(paste0(paste0("VE(0-18)$>$",lowerVEuncPower*100),"\\%"),collapse=" & ")),"")," \\\\\n \\hline"), paste0("\\hline \\hline \\multicolumn{",5+length(lowerVEuncPower),"}{l}{{\\footnotesize $^{\\ast}$VE halved in the first 6 months}} \\\\\n \\multicolumn{",5+length(lowerVEuncPower),"}{l}{\\footnotesize N=",N.pla,":",N.vax," placebo:vaccine group} \\\\\n \\multicolumn{",5+length(lowerVEuncPower),"}{l}{\\footnotesize ",infRate*100,"\\% annual incidence in the placebo group} \\\\\n \\multicolumn{",5+length(lowerVEuncPower),"}{l}{\\footnotesize 5\\% annual dropout} \\\\\n \\multicolumn{",5+length(lowerVEuncPower),"}{l}{\\footnotesize ",ifelse(estimand=="combined","Cox \\& cumulative incidence-based non-efficacy monitoring}",ifelse(estimand=="cuminc", paste0("Cumulative incidence-based non-efficacy monitoring",ifelse(laggedMonitoring," incl. post-6 months monitoring",""),"}"), paste0("Cox model-based non-efficacy monitoring",ifelse(laggedMonitoring," incl. post-6 months monitoring",""),"}")))," \\\\\n \\multicolumn{",5+length(lowerVEuncPower),"}{l}{\\footnotesize ",ifelse(estimand=="cox", "Cox model-based high efficacy monitoring", "Cumulative incidence-based high efficacy monitoring"),"} \\\\\n \\multicolumn{",5+length(lowerVEuncPower),"}{l}{\\footnotesize ",ifelse(estimand=="cox","Cox model-based unconditional power","Cumulative incidence-based unconditional power"),"}")))) @ \begin{figure}[H] \begin{center} <>= res <- subset(res, aveHR<=1) par(mar=c(4.5,5,3,1.5), las=1, cex.axis=1.2, cex.lab=1.3, cex.main=1.2) plot(-1, 0, xlim=c(0,0.8), ylim=0:1, xaxt="n", yaxt="n", xlab="Average VE(0-18)", ylab="Probability") axis(side=1, at=seq(0,0.8,by=0.1), labels=paste0(seq(0,80,by=10),"%")) axis(side=2, at=c(0,0.1,seq(0.2,1,by=0.2))) axis(side=2, at=0.85) axis(side=4, at=c(0,0.1,0.85,seq(0.2,1,by=0.2)), labels=FALSE) abline(h=c(0.1,0.85), lty="dotted", lwd=2) with(res, lines(1-aveHR, harm/100, type="b", lwd=2, lty="solid", col="darkorange")) with(res, lines(1-aveHR, higheff/100, type="b", lwd=2, lty="longdash", col="red3")) with(res, lines(1-aveHR, noneff/100, type="b", lwd=2, lty="dotdash", col="blue")) if (!is.null(lowerVEuncPower)){ cols <- c("purple","black") ltys <- c("dashed","solid") for (i in 1:length(lowerVEuncPower)){ lines(1-res$aveHR, res[,paste0("design",i)]/100, type="b", lwd=2, lty=ltys[i], col=cols[i]) } } text(0.2, res[res$aveHR==0.8, "noneff"]/100, "Non-efficacy", cex=1.1, pos=2, offset=1.1) text(0.15, 0.06, "Potential harm", adj=0, cex=1.1) text(0.65, 0.25, "High\nefficacy", adj=0, cex=1.1) text(0.28, 0.93, "Unconditional power\n[VE(0-18)>0%]", cex=1.1, pos=4) #text(0.45, 0.6, "Unconditional power\n[VE(0-18)>25%]", cex=1.1, pos=4) text(0.33, 0.5, paste0("N=",N.pla,":",N.vax," placebo:vaccine group\n",infRate*100,"% annual placebo group incidence\n5% annual dropout\nVE halved in first 6 months\n",ifelse(estimand=="combined","Cox & cum inc-based ",ifelse(estimand=="cuminc","Cum inc ","Cox model ")),ifelse(laggedMonitoring,"non-eff w/ post-6 mo monitoring\n","non-eff monitoring\n"),ifelse(estimand=="cox","Cox model high eff & unconditional power","Cum inc high eff & unconditional power"),"\n18-month accrual; average 309 per month"), cex=1, pos=4) cat("\\caption{Probabilities of reaching each possible trial monitoring outcome",ifelse(!is.null(lowerVEuncPower),", and unconditional power to reject the specified null hypothesis","")," for a 2-arm study design with ",N.pla," placebo and ",N.vax," vaccine recipients}", sep="") @ \end{center} \end{figure} \begin{figure}[H] \begin{center} <>= fullVE <- function(aveVE, vePeriods){ aveVE * (vePeriods[3]-1)/((vePeriods[3]-1)-(vePeriods[2]-1)/2) } pwrPP.filename <- paste0("VEpwPP_nPlac=",N.pla,"_nVacc=",N.vax,"_infRate=",infRate,".RData") load(file.path(srcDir, pwrPP.filename)) par(mar=c(4.5,5,3,1.5), las=1, cex.axis=1.2, cex.lab=1.3, cex.main=1.2) plot(-1, 0, xlim=c(0,1), ylim=0:1, xaxt="n", yaxt="n", xlab="Cumulative VE(6.5-18) = [1-Cumul. Incidence Ratio (V/P) by Mo 18] x 100%", ylab="") mtext("Unconditional Power [VE(6.5-18) > 0%]", side=2, line=3.5, cex=1.3, las=0) axis(side=1, at=seq(0,1,by=0.1), labels=paste0(seq(0,100,by=10),"%")) axis(side=2, at=seq(0,1,by=0.2)) axis(side=2, at=0.85) axis(side=4, at=c(seq(0,1,by=0.2),0.85), labels=FALSE) abline(h=0.85, lty="dotted", lwd=2) colNames <- c("black", "darkred", "blue", "purple") ltyNames <- c("solid", "dashed", "dotted", "longdash") for (i in 1:4){ # the x-coordinates are the true values of VE(6.5-18) x <- fullVE(seq(0,0.8,by=0.1),c(1,27,79)) # the y-coordinates are powers to reject H0: VE(6.5-18)<=0 in favor of H1: VE(6.5-18)>0 y <- sapply(pwList, function(oneAveVEData){ oneAveVEData$VEpwPP[i] }) lines(x, y, type="b", lty=ltyNames[i], lwd=2, col=colNames[i]) } text(0.55, 0.15, paste0("N=",N.pla,":",N.vax," MITT placebo:vaccine\n",infRate*100,"% annual incidence placebo\n5% annual dropout\nVE halved in first 6 months"), adj=0, cex=1.1) legend(0.55, 0.65, lty=ltyNames, col=colNames, lwd=2, title="Missing Vaccination\nProbability", legend=c("0% (MITT)",paste0(c(5,10,15),"%")), bty="n") @ \caption{Unconditional power to reject the null hypothesis H$_0$: VE(6.5--18)$\leq$0\% in per-protocol cohorts with a varying probability of a missing vaccination} \end{center} \end{figure} \newpage \begin{figure}[H] \begin{center} <>= cumProbTrialDuration1VaxArm <- function(simTrial.filename, monitorTrial.filename, dirname, stage2=156){ load(file.path(dirname, simTrial.filename)) trialData <- trialObj[["trialData"]] load(file.path(dirname, monitorTrial.filename)) bounds <- sapply(out, function(trial){ trial[[1]]$boundHit }) stopTimes <- sapply(out, function(trial){ trial[[1]]$stopTime }) for (i in 1:length(bounds)){ if (bounds[i] %in% c("Eff","HighEff")){ data.i <- trialData[[i]] endStage2 <- max(data.i$entry + stage2) trialStop <- min (max(data.i$exit), endStage2) stopTimes[i] <- trialStop } } return(quantile(stopTimes, prob=seq(0,1,by=0.01))) } aveVE <- c(0, seq(0.2,0.5,by=0.1)) # store the quantiles of the trial duration for every 'aveVE' quantiles <- matrix(0, nrow=length(aveVE), ncol=length(seq(0,1,by=0.01))) for (i in 1:length(aveVE)){ simTrial.filename <- paste0("simTrial_nPlac=",N.pla,"_nVacc=",N.vax,"_aveVE=",aveVE[i],"_infRate=",infRate,".RData") monitorTrial.filename <- paste0("monitorTrial_nPlac=",N.pla,"_nVacc=",N.vax,"_aveVE=",aveVE[i],"_infRate=",infRate,"_",estimand,".RData") quantiles[i,] <- cumProbTrialDuration1VaxArm(simTrial.filename, monitorTrial.filename, srcDir) } # convert into months quantiles <- quantiles/(52/12) par(mar=c(4.5,5,3,1.5), las=1, cex.axis=1.2, cex.lab=1.3, cex.main=1.2) plot(-10, 0, xlim=c(0,54), ylim=0:1, xaxt="n", yaxt="n", xlab="Vaccine Regimen's Evaluation Duration (Months)", ylab="Cumulative Probability") axis(side=1, at=c(0,seq(10,50,by=10),54)) axis(side=2, at=seq(0,1,by=0.2)) axis(side=4, at=seq(0,1,by=0.2), labels=FALSE) segments(x0=18, y0=0, y1=0.75, lwd=2, col="gray30") segments(x0=24, y0=0, y1=0.875, lwd=2, col="gray30") segments(x0=30, y0=0, y1=1, lwd=2, col="gray30") linesCol <- c("black", "blue", "red", "green", "darkorange") for (i in 1:length(aveVE)){ lines(quantiles[i,], seq(0,1,by=0.01), lwd=2, col=linesCol[i]) } legend(x=0.5, y=0.975, legend=paste(aveVE*100,"%",sep=""), col=linesCol, lwd=2, cex=1.2, title="Average\nVE(0-18)", bty="n") text(0, 0.4, paste("N=",N.pla,":",N.vax," pla:vac\n",infRate*100,"% annual placebo incidence\n5% annual dropout\nVE halved in first 6 months", sep=""), adj=0, cex=1.2) text(x=18, y=0.7, labels="Accrual\ncompleted", cex=1.2, pos=2) text(x=24, y=0.825, labels="Vaccinations through\nmonth 6 completed", cex=1.2, pos=2) text(x=30, y=0.95, labels="Vaccinations through\nmonth 12 completed", cex=1.2, pos=2) @ \caption{Duration of a vaccine regimen's evaluation ($n=$\Sexpr{N.pla} in the placebo arm and $n=$\Sexpr{N.vax} in the vaccine arm)} \end{center} \end{figure} <>= cumProbTrialDuration <- function(trialDataCens.filename, dirname, stage2=156){ load(file.path(dirname, trialDataCens.filename)) trialStopTime <- numeric(length(trialListCensor)) for (i in 1:length(trialListCensor)){ dataI <- trialListCensor[[i]] trialStopTime[i] <- max(pmin(dataI$entry + stage2, dataI$exit)) } return(quantile(trialStopTime, prob=seq(0,1,by=0.01))) } if (N.vax.arms>1){ aveVE <- matrix(c(0,0,0.4,0,0.4,0.4), nrow=3, ncol=2, byrow=TRUE) cat("\\newpage") cat("\\begin{figure}[H]") cat("\\begin{center}") # store the quantiles of the trial duration for every 'aveVE' quantiles <- matrix(0, nrow=NROW(aveVE), ncol=length(seq(0,1,by=0.01))) for (j in 1:NROW(aveVE)){ trialDataCens.filename <- paste0("trialDataCens_nPlac=",N.pla,"_nVacc=",paste(rep(N.vax,N.vax.arms), collapse="_"),"_aveVE=",paste(aveVE[j,], collapse="_"),"_infRate=",infRate,"_",estimand,".RData") quantiles[j,] <- cumProbTrialDuration(trialDataCens.filename, srcDir) } # convert into months quantiles <- quantiles/(52/12) par(mar=c(4.5,5,3,1.5), las=1, cex.axis=1.2, cex.lab=1.3, cex.main=1.2) plot(-10, 0, xlim=c(0,54), ylim=0:1, xaxt="n", yaxt="n", xlab=paste(N.vax.arms,"-Vaccine-Arm Trial Duration (Months)",sep=""), ylab="Cumulative Probability") axis(side=1, at=c(0,seq(10,50,by=10),54)) axis(side=2, at=seq(0,1,by=0.2)) axis(side=4, at=seq(0,1,by=0.2), labels=FALSE) segments(x0=18, y0=0, y1=0.5, lwd=2, col="gray30") segments(x0=24, y0=0, y1=0.65, lwd=2, col="gray30") segments(x0=30, y0=0, y1=0.8, lwd=2, col="gray30") linesCol <- c("black", "blue", "red", "green", "darkorange")[1:NROW(aveVE)] for (i in 1:NROW(aveVE)){ lines(quantiles[i,], seq(0,1,by=0.01), lwd=2, col=linesCol[i]) } aveVEcomb <- NULL for (k in 1:NROW(aveVE)){ aveVEcomb <- c(aveVEcomb,paste("(",paste(paste(aveVE[k,]*100, "%", sep=""), collapse=", "),")",sep="")) } legend(x=0.5, y=0.975, legend=aveVEcomb, col=linesCol, lwd=2, cex=1.2, title="Average\nVE(0-18)", bty="n") text(0, 0.2, paste0("N=",N.pla,":",paste(rep(N.vax,N.vax.arms),collapse=":")," pla:",paste(rep("vac",N.vax.arms),collapse=":"),"\n","18-month accrual; average 309 per month\n",infRate*100,"% annual placebo incidence\n5% annual dropout\nVE halved in first 6 months"), adj=0, cex=1.2) text(x=18, y=0.45, labels="Accrual\ncompleted", cex=1.2, pos=2) text(x=24, y=0.6, labels="Vaccinations through\nmonth 6 completed", cex=1.2, pos=2) text(x=30, y=0.75, labels="Vaccinations through\nmonth 12 completed", cex=1.2, pos=2) cat("\\caption{Total trial duration for the evaluation of ",N.vax.arms," vaccine regimens ($N=",N.vax,"$ per arm) versus one placebo arm ($N=",N.pla,"$)}", sep="") cat("\\end{center}") cat("\\end{figure}") } @ <>= infThroughWeekVec <- c(78, 104, 156) pMissVax <- 0.05 for (infThroughWeek in infThroughWeekVec){ aveVE <- 0.4 p <- c(0.01, 0.05, 0.25, 0.5, 0.75, 0.95, 0.99) ppname <- paste0("pp",which(pMissVax==seq(0,0.15,0.05))) res <- as.data.frame(matrix(0, nrow=2, ncol=length(p)+1)) trialDataCens.filename <- paste0("trialDataCens_nPlac=",N.pla,"_nVacc=",N.vax,"_aveVE=",aveVE,"_infRate=",infRate,"_",estimand,".RData") load(file.path(srcDir, trialDataCens.filename)) # for each trial, create a vector of Stage 1 infection counts by treatment infecListMITT <- infecListPP <- as.list(NULL) for (i in 1:length(trialListCensor)){ dataI <- trialListCensor[[i]] dataI$futime <- dataI$exit - dataI$entry infecListMITT[[i]] <- as.vector( table( as.factor(dataI$trt)[dataI$event & dataI$futime>26 & dataI$futime<=infThroughWeek] ) ) infecListPP[[i]] <- as.vector( table( as.factor(dataI$trt)[dataI[[ppname]]==1 & dataI$event & dataI$futime>26 & dataI$futime<=infThroughWeek] ) ) } # for each trial, calculate the number of infections diagnosed between 28-infThroughWeek weeks among vaccine recipients NinfMITT <- sapply(infecListMITT, function(vectInfbyTx){ sum(vectInfbyTx[-1], na.rm=TRUE) }) NinfPP <- sapply(infecListPP, function(vectInfbyTx){ sum(vectInfbyTx[-1], na.rm=TRUE) }) res[1,1] <- mean(NinfMITT) res[1,-1] <- quantile(NinfMITT, prob=p) res[2,1] <- mean(NinfPP) res[2,-1] <- quantile(NinfPP, prob=p) res <- round(res, 0) colnames(res) <- c("Mean", paste(p*100,"%",sep="")) print(xtable(res, digits=0, align=rep("c",length(p)+2), caption=paste0("Distribution of the number of month 6.5--",infThroughWeek*12/52," infections per vaccine group for use in the immune correlates analysis, for vaccine regimens with average VE of ",aveVE*100,"\\%, halved in the initial 6 months ($N=",N.pla,"$ in the placebo group, $N=",N.vax,"$ in each vaccine group, and ",pMissVax*100,"\\% conditional probability of having missed a vaccination given HIV-negative and ongoing at the Month 6 [Week 26] visit).") ), table.placement="H", caption.placement="top", include.rownames=FALSE, include.colnames=FALSE, #scalebox=0.9, hline.after=NULL, add.to.row=list(pos=list(-1,0,1,2), command=c(paste0("\\hline \\hline & \\multicolumn{7}{c}{Percentiles of distribution of number} \\\\\n & \\multicolumn{7}{c}{of month 6.5--",infThroughWeek*12/52," infections per vaccine arm} \\\\\n \\cline{2-8}"),paste0("Mean & 1\\% & 5\\% & 25\\% & 50\\% & 75\\% & 95\\% & 99\\% \\\\\n \\hline \\multicolumn{8}{c}{Month 6.5--",infThroughWeek*12/52," infections in MITT cohort} \\\\\n"),paste0("\\hline \\multicolumn{8}{c}{Month 6.5--",infThroughWeek*12/52," infections in per-protocol cohort} \\\\\n"),paste0("\\hline \\hline \\multicolumn{8}{l}{\\footnotesize N=",N.pla,":",N.vax," MITT placebo:vaccine} \\\\\n \\multicolumn{8}{l}{\\footnotesize ",pMissVax*100,"\\% probability of a missing vaccination} \\\\\n \\multicolumn{8}{l}{\\footnotesize ",infRate*100,"\\% annual placebo group incidence} \\\\\n \\multicolumn{8}{l}{\\footnotesize 5\\% annual dropout} \\\\\n \\multicolumn{8}{l}{\\footnotesize Average VE=",aveVE*100,"\\%, halved VE in the first 6 months}"))) ) } @ \newpage <>= aveVEsets=matrix(c(0.4,0, 0.4,0.2, 0.4,0.3, 0.5,0.3, 0.5,0.4, 0.6,0.3, 0.6,0.4), nrow=7, ncol=2, byrow=TRUE) typeVec <- "head" for (type in typeVec){ rankSelectPwVector <- headPwVector <- NULL for (i in 1:NROW(aveVEsets)){ rankSelectPwr.filename <- paste0("rankSelectPwr_nPlac=",N.pla,"_nVacc=",paste(rep(N.vax,N.vax.arms), collapse="_"),"_aveVE=",paste(aveVEsets[i,], collapse="_"),"_infRate=",infRate,"_",estimand,".RData") load(file.path(srcDir, rankSelectPwr.filename)) if (type=="head"){ rankSelectPwVector <- c(rankSelectPwVector, out$rankSelectPw*100) # power to correctly identify the winning efficacious vaccine headPwVector <- c(headPwVector, out$headHeadPw[1,1]*100) # power to detect positive relative VE(0-stage1) of Vx1 vs Vx2 } if (type=="pool"){ rankSelectPwVector <- c(rankSelectPwVector, out$poolRankSelectPw*100) # power to correctly identify the best pool headPwVector <- c(headPwVector, out$poolHeadPw[1,1]*100) # power to detect positive relative VE(0-stage1) of, e.g., Ve3-Vx4 vs Vx1-Vx2 } } aveVEchars <- apply(aveVEsets, 1, function(aveVE){ paste0("(",paste(aveVE*100,collapse=", "),")") }) nRows <- length(aveVEchars) res <- data.frame(aveVEchars=aveVEchars, headPwVector=headPwVector, rankSelectPwVector=rankSelectPwVector) print(xtable(res, digits=1, align=rep("c",4), caption=ifelse(type=="head", "Power ($\\times 100$) to detect that relative VE(0--18) $> 0\\%$ comparing head-to-head vaccine regimens 1 vs. 2 and VE(0--18) $> 0\\%$ for vaccine regimen 1, and probability ($\\times 100$) of correct ranking and selection of the winning most efficacious vaccine regimen", "Power ($\\times 100$) to detect that relative VE(0--18) $> 0\\%$ comparing head-to-head pooled vaccine regimens 3--4 vs. 1--2 and VE(0--18) $> 0\\%$ for the pooled vaccine regimen 3--4, and probability ($\\times 100$) of correct ranking and selection among the pooled pairs of the winning most efficacious regimen" ) ), table.placement="H", caption.placement="top", include.rownames=FALSE, include.colnames=FALSE, #scalebox=0.9, hline.after=0, add.to.row=list(pos=list(-1,0,nRows), command=c(paste0("\\hline \\hline & Power ($\\times 100$) & \\\\\n True average VE (\\%)$^1$ & ",ifelse(type=="head","Vx1 vs. Vx2","Vx3-4 vs. Vx1-2")," \\& & Probability ($\\times 100$) \\\\\n"),paste0("(Vx1, Vx2) & ",ifelse(type=="head","Vx1 vs. Placebo","Vx3-4 vs. Placebo"),"$^2$ & select best ",ifelse(type=="head","vaccine","pooled Vx"),"$^3$ \\\\\n"),paste0("\\hline\\hline \\multicolumn{3}{l}{\\footnotesize $^1$ VE halved in the first 6 months} \\\\\n \\multicolumn{3}{l}{\\footnotesize $^2$ Cumulative incidence-based Wald tests of both ",ifelse(type=="head","Vx1/Vx2","Vx3-4/Vx1-2")," and} \\\\\n \\multicolumn{3}{l}{\\footnotesize \\; \\,",ifelse(type=="head","Vx1/Placebo","Vx3-4/Placebo")," VE(0--18) with 1-sided $\\alpha=0.025$} \\\\\n \\multicolumn{3}{l}{\\footnotesize $^3$ Correct selection $=$ ",ifelse(type=="head","Vx1","pooled Vx3-4")," has highest estimated VE(0--36) and} \\\\\n \\multicolumn{3}{l}{\\footnotesize \\; \\, VE(0--18) significantly $>0\\%$} \\\\\n \\multicolumn{3}{l}{\\footnotesize N=",N.pla,":",paste(rep(N.vax,N.vax.arms),collapse=":")," pla:",paste(rep("vac",N.vax.arms),collapse=":")," group} \\\\\n \\multicolumn{3}{l}{\\footnotesize 18-month accrual; average 309 per month} \\\\\n \\multicolumn{3}{l}{\\footnotesize ",infRate*100,"\\% annual incidence in the placebo group} \\\\\n \\multicolumn{3}{l}{\\footnotesize 5\\% annual dropout} \\\\\n \\multicolumn{3}{l}{\\footnotesize Cumulative incidence-based monitoring}"))) ) } @ \newpage <>= aveVE <- c(0, 0.4) p <- c(0.01, 0.025, 0.05, seq(0.1,0.9,by=0.1), 0.95, 0.975, 0.99) res <- as.data.frame(matrix(0, nrow=length(aveVE), ncol=length(p))) for (i in 1:length(aveVE)){ simTrial.filename <- paste0("simTrial_nPlac=",N.pla,"_nVacc=",paste(rep(N.vax,N.vax.arms), collapse="_"),"_aveVE=",paste(rep(aveVE[i],N.vax.arms), collapse="_"),"_infRate=",infRate,".RData") load(file.path(srcDir, simTrial.filename)) Ninf <- sapply(trialObj$NinfStage1, function(vectInfbyTx){ vectInfbyTx[1] + max(vectInfbyTx[-1], na.rm=TRUE) }) res[i,] <- quantile(Ninf, prob=p) } res <- round(res, 0) res <- cbind(c("0%","40%"), res) colnames(res) <- c("Ave VE (0-18)", paste(p*100,"%",sep="")) print(xtable(res, digits=0, align=c("c","l",rep("c",15)), caption=paste("Distribution of the number of Stage~1 infections pooled over the placebo group and the vaccine group with the maximum number of infections, ignoring sequential monitoring for potential harm, non-efficacy, and high efficacy (n=",N.pla," in the placebo group and n=",N.vax," in each vaccine group)",sep="") ), table.placement="H", caption.placement="top", include.rownames=FALSE, include.colnames=FALSE, scalebox=0.85, hline.after=0, add.to.row=list(pos=list(-1,0,NROW(res)), command=c("\\hline \\hline Ave VE & \\multicolumn{15}{c}{Percentiles of distribution of number of Stage 1 infections} \\\\\n","(0-18)$^{\\ast}$ & 1\\% & 2.5\\% & 5\\% & 10\\% & 20\\% & 30\\% & 40\\% & 50\\% & 60\\% & 70\\% & 80\\% & 90\\% & 95\\% & 97.5\\% & 99\\% \\\\\n",paste0("\\hline \\hline \\multicolumn{16}{l}{{\\footnotesize $^{\\ast}$VE halved in the first 6 months}} \\\\\n \\multicolumn{16}{l}{\\footnotesize N=",N.pla,":",paste(rep(N.vax,N.vax.arms),collapse=":")," pla:",paste(rep("vac",N.vax.arms),collapse=":")," group} \\\\\n \\multicolumn{16}{l}{\\footnotesize 18-month accrual; average 309 per month} \\\\\n \\multicolumn{16}{l}{\\footnotesize ",infRate*100,"\\% annual incidence in the placebo group} \\\\\n \\multicolumn{16}{l}{\\footnotesize 5\\% annual dropout}"))) ) aveVE <- c(0, 0.4) p <- c(0.01, 0.025, 0.05, seq(0.1,0.9,by=0.1), 0.95, 0.975, 0.99) res <- as.data.frame(matrix(0, nrow=length(aveVE), ncol=length(p))) for (j in 1:length(aveVE)){ #read in the censored data trialDataCens.filename <- paste0("trialDataCens_nPlac=",N.pla,"_nVacc=",paste(rep(N.vax,N.vax.arms), collapse="_"),"_aveVE=",paste(rep(aveVE[j],N.vax.arms), collapse="_"),"_infRate=",infRate,"_",estimand,".RData") load(file.path(srcDir, trialDataCens.filename)) # for each trial, create a vector of Stage 1 infection counts by treatment infecList <- as.list(NULL) for (i in 1:length(trialListCensor)){ dataI <- trialListCensor[[i]] dataI$futime <- dataI$exit - dataI$entry stage1ind <- dataI$event & dataI$futime<=78 infecList[[i]] <- as.vector( table( as.factor(dataI$trt)[stage1ind] ) ) } # for each trial, calculate the number of Stage 1 infections pooled over all treatment arms Ninf.all <- sapply(infecList, function(vectInfbyTx){ sum(vectInfbyTx, na.rm=TRUE) }) # for each trial, calculate the number of Stage 1 infections pooled over the placebo arm and the vaccine arm with # the maximum number of infections Ninf <- sapply(infecList, function(vectInfbyTx){ vectInfbyTx[1] + max(vectInfbyTx[-1], na.rm=TRUE) }) res[j,] <- quantile(Ninf, prob=p) } res <- round(res, 0) res <- cbind(c("0%","40%"), res) colnames(res) <- c("Average VE (0-18)", paste(p*100,"%",sep="")) print(xtable(res, digits=0, align=c("c","l",rep("c",15)), caption=paste("Distribution of the number of Stage 1 infections pooled over the placebo group and the vaccine group with the maximum number of infections, accounting for sequential monitoring for potential harm, non-efficacy, and high efficacy (n=",N.pla," in the placebo group and n=",N.vax," in each vaccine group)",sep="") ), table.placement="H", caption.placement="top", include.rownames=FALSE, include.colnames=FALSE, scalebox=0.85, hline.after=NULL, add.to.row=list(pos=list(-1,0,NROW(res)), command=c("\\hline \\hline Ave VE & \\multicolumn{15}{c}{Percentiles of distribution of number of Stage 1 infections} \\\\\n","(0-18)$^{\\ast}$ & 1\\% & 2.5\\% & 5\\% & 10\\% & 20\\% & 30\\% & 40\\% & 50\\% & 60\\% & 70\\% & 80\\% & 90\\% & 95\\% & 97.5\\% & 99\\% \\\\\n \\hline", paste0("\\hline \\hline \\multicolumn{16}{l}{{\\footnotesize $^{\\ast}$VE halved in the first 6 months}} \\\\\n \\multicolumn{16}{l}{\\footnotesize N=",N.pla,":",paste(rep(N.vax,N.vax.arms),collapse=":")," pla:",paste(rep("vac",N.vax.arms),collapse=":")," group} \\\\\n \\multicolumn{16}{l}{\\footnotesize 18-month accrual; average 309 per month} \\\\\n \\multicolumn{16}{l}{\\footnotesize ",infRate*100,"\\% annual incidence in the placebo group} \\\\\n \\multicolumn{16}{l}{\\footnotesize 5\\% annual dropout} \\\\\n \\multicolumn{16}{l}{\\footnotesize Cumulative incidence-based monitoring}"))) ) @ \newpage \begin{figure}[H] \begin{center} <>= harmData <- read.csv(file.path(srcDir, paste0("harmBounds_N=60_alphaPerTest=0.0140_pVacc=",round(N.vax/(N.pla+N.vax),2),".csv"))) harmData <- subset(harmData, N>=10) par(mar=c(5.5,5,1,1), las=1, cex.axis=1.1, cex.lab=1.2) plot(-10,0, xlim=c(10,60), ylim=c(9,40), xaxt="n", yaxt="n", xlab="", ylab="", type="n") mtext("Total Number of Infections\n(Placebo + 1 Vaccine Arm)", side=1, line=3.5, cex=1.2) mtext("Number of Vaccine Infections\n(1 Vaccine Arm)", side=2, line=2.5, cex=1.2, las=0) axis(side=1, at=seq(10,60,by=5)) axis(side=2, at=seq(10,40,by=5)) axis(side=3, at=seq(10,60,by=5), labels=FALSE) axis(side=4, at=seq(10,40,by=5), labels=FALSE) abline(0,1, col="gray30") text(28,32,"Y=X", col="gray30") points(harmData$N, harmData$V, pch=19, col="darkred", cex=0.5) for (i in c(20,40,60)){ segments(x0=i, y0=8, y1=harmData$V[harmData$N==i], lwd=2, lty="dashed", col="darkorange") segments(x0=10, x1=i, y0=harmData$V[harmData$N==i], lwd=2, lty="dashed", col="darkorange") } text(42, 21, "Potential Harm\nStopping\nBoundaries", adj=0, cex=1.1, font=2, col="blue") @ \end{center} \end{figure} <>= VEbinomModel <- function(p, r){ 1 - r*p/(1-p) } UL95 <- function(pHat, n, r, bound){ VEbinomModel(pHat, r) + qnorm(0.975)*sqrt(pHat*(r^2)/(n*((1-pHat)^3))) - bound } LL95 <- function(pHat, n, r, bound){ VEbinomModel(pHat, r) - qnorm(0.975)*sqrt(pHat*(r^2)/(n*((1-pHat)^3))) - bound } getpHatNonEff <- function(nVec, randRatioPV, upperBound){ return(sapply(nVec, function(n){ uniroot(UL95, interval=c(0.01,0.9), n=n, r=randRatioPV, bound=upperBound)$root })) } getpHatHighEff <- function(nVec, randRatioPV, lowerBound){ return(sapply(nVec, function(n){ uniroot(LL95, interval=c(0.01,0.9), n=n, r=randRatioPV, bound=lowerBound)$root })) } randRatio <- N.pla/N.vax # calculate the MC average number of infections (pooled over placebo + 1 vaccine arm) triggering the first non-efficacy IA for trials with average VE=20% load(file.path(srcDir,paste0("monitorTrial_nPlac=",N.pla,"_nVacc=",N.vax,"_aveVE=0.2_infRate=",infRate,"_",estimand,".RData"))) firstNoneffCnt <- round(mean(do.call("c", sapply(out, function(trial){ trial[[1]]$firstNonEffCnt }))),0) nInf <- seq(firstNoneffCnt,240,by=20) pHat <- getpHatNonEff(nVec=nInf, randRatioPV=randRatio, upperBound=0.4) nVax <- round(pHat*nInf,0) nSplitVaxPla <- paste0(nVax,":",nInf-nVax) VEhat <- VEbinomModel(pHat, randRatio)*100 tableInfSplits <- data.frame(nInf, nSplitVaxPla, round(VEhat,0)) @ \newpage \begin{figure}[H] \begin{center} <>= par(mar=c(5.5,5,1,1), las=1, cex.axis=1.1, cex.lab=1.2) plot(-10,0, xlim=c(0,240), ylim=c(-25,20), xaxt="n", xlab="", ylab="", type="n") mtext("Total Number of Infections by Month 18\n(Placebo + 1 Vaccine Arm)", side=1, line=3.5, cex=1.2) mtext("Vaccine Efficacy (%)", side=2, line=2.5, cex=1.2, las=0) axis(side=1, at=seq(0,240,by=20)) axis(side=2, at=seq(-25,20,by=5)) axis(side=3, at=seq(0,240,by=20), labels=FALSE) axis(side=4, at=seq(-25,20,by=5), labels=FALSE) for (i in 1:length(nInf)){ segments(x0=nInf[i], y0=-30, y1=VEhat[i], lwd=2, col="darkred") } points(nInf, VEhat, pch=19, col="darkred", cex=0.6) lines(nInf, VEhat, col="darkred") text(5, 15, "Non-Efficacy\nStopping Boundaries", adj=0, cex=1.1, font=2, col="blue") @ \end{center} \end{figure} <>= print(xtable(tableInfSplits, align=rep("c",4), digits=rep(0,4) ), include.rownames=FALSE, include.colnames=FALSE, hline.after=c(-1,0), add.to.row=list(pos=list(0,NROW(tableInfSplits)), command=c("\\multicolumn{3}{c}{Non-Efficacy Stopping$^{\\ast}$}\\\\\n\\hline Total\\TT & Infections & $\\widehat{\\textrm{VE}}^{\\dagger}$\\\\\n Infections & Split V:P & (\\%)\\\\\n", "\\hline \\multicolumn{3}{l}{\\footnotesize $^{\\ast}$Ave VE$=$20\\%, halved in first 6 mo.}\\\\\n\\multicolumn{3}{l}{\\footnotesize $^{\\dagger}$Based on a binomial model}")) ) @ \end{document}