#' @title Simulation of one Phase I study
#' @description Simulation of one Phase I study with classical 3+3 design,
#' BOIN, Keyboard or i3+3 design
#' @param phi target toxicity rate (needed for all designs)
#' @param phi1 for BOIN/i3+3 design: the highest DLT rate that is deemed subtherapeutic (i.e., underdosing),
#' such that dose escalation should be made, see \code{\link{get_BOIN_rules}}
#' @param phi2 for BOIN/i3+3 design: the lowest DLT rate that is deemed overly toxic (i.e., overdosing),
#' such that dose de-escalation is required, see \code{\link{get_BOIN_rules}}
#' @param start_dose which dose to start with (could be for instance the second dose level)
#' @param halfkey for Keyboard design: half length of key in keyboard design, see \code{\link{get_Keyboard_rules}}
#' @param maxtox P(DLT rate>phi|c,n)<=maxtox as value for maxprob in \code{\link{get_Elim_rules}}
#' @param N number of patients (not needed for design="3+3")
#' @param truerate scenario of true DLT rates for each dose level. This also defines the number of dose levels
#' @param cohortsize cohort size, default is 3
#' @param maxN If N(patients) at next dose level >=maxN, then stop algorithm, default is 9
#' @param maxNdec condition for maxN. Two options: "STAY"= stop at maxN if decision is stay;
#' "ANY"=stop at maxN at any decision
#' @param acc_tit logical indicator if algorithm should start with accelerated titration=algorithm starts at dose level with first DLT
#' @param dose_no_titr first dose at which no accelerated titration anymore
#' @param BOIN_add33_rule logical indicator: modify the decision from de-escalation to stay when observing 1 DLT out of 3 patients
#' @param design "BOIN" or "Keyboard" or "i3+3"
#' @param MTD_safer imposes that the MTD should be for 1)BOIN:<lambda_d, 2)Keyboard:<phi+halfkey, 3)i3+3:<phi2
#' @param seed define seed
#' @param sim if "NO", calculate decision rules, define "YES" if done within simulation program for multiple studies
#' @param env parent environment to pass objects from \code{\link{ph1_sim_OC}}
#' @param hist complete history of simulation in results
#' @return a data.frame with elements
#' \itemize{
#' \item dose_level: dose level
#' \item dose: actual dose to be used
#' \item incr_ratio: The ratio of successive doses
#' }
#' @references Liu S, Yuan Y. Bayesian Optimal Interval Designs for Phase I Clinical Trials.
#' J R Stat Soc Ser C Appl Stat. 2015;64:507–23
#' Yan F, Mandrekar SJ, Yuan Y. Keyboard: A Novel Bayesian Toxicity Probability Interval Design for Phase I Clinical Trials. Clin Cancer Res. 2017; 23: 3994–4003
#' @export
#' @examples
#' ph1_1sim(phi=0.3, phi1=0.6*0.3, phi2=1.4*0.3, maxtox=0.95, N=15,
#' truerate=c(0.20 ,0.25 ,0.30 ,0.40 ,0.60 ,0.75), cohortsize=3, maxN=9, acc_tit=0,
#' design="BOIN", MTD_safer=TRUE)
#' ph1_1sim(phi=0.3, maxtox=0.95, N=30,
#' truerate=c(0.20 ,0.25 ,0.30 ,0.40 ,0.60 ,0.75), cohortsize=3, maxN=9, acc_tit=1,
#' design="Keyboard", MTD_safer=TRUE, halfkey=0.05)
#' ph1_1sim(phi=0.3, phi1=0.6*0.3, phi2=1.4*0.3, maxtox=0.95, N=15,
#' truerate=c(0.20 ,0.25 ,0.30 ,0.40 ,0.60 ,0.75), cohortsize=3, maxN=9, acc_tit=0,
#' design="i3+3", MTD_safer=TRUE)
#' ph1_1sim(truerate=c(0.20 ,0.25 ,0.30 ,0.40 ,0.60 ,0.75),acc_tit=1,design="3+3")
# truerate=c(0.1 ,0.15 ,0.2 );phi=0.25; phi1=0.6*0.25 ; phi2=1.4*0.25 ; start_dose=2; maxtox=0.95 ; N=15; cohortsize=3; maxN=15;
# acc_tit=0; dose_no_titr=NULL; design="BOIN"; MTD_safer=TRUE; seed=NULL; hist=1;sim="NO"; BOIN_add33_rule=F
ph1_1sim <- function (phi, phi1=NULL, phi2=NULL, start_dose=1, maxtox=NULL, N=NULL,truerate=NULL,cohortsize=NULL,maxN=9, maxNdec="STAY", acc_tit,
dose_no_titr=NULL, BOIN_add33_rule=F, design, MTD_safer=TRUE, seed = NULL, halfkey=NULL, sim="NO", env=parent.frame(),hist=0){ # MTD_safer: MTD should be <lambda_d
if (hist==1){result_list<-list()}
if (!is.null(seed)) {
set.seed(seed)
}
#-----------------------------------------------------------#
# First dose: Create empty vectors and initialize variables #
#-----------------------------------------------------------#
if(is.null(dose_no_titr)){dose_no_titr<-length(truerate)}
ndose<-length(truerate)
dose_in_the_running <- dose <- 1:ndose
npt <- ndlt <- MTD <- rep(0,ndose)
p_iso<-rep(NA,ndose)
dose_inv <- start_dose # initialize "dose under investigation"
n <- 0 # number of patients in study, for first dose
mtd <- NA # initializing scalar for dose level equal to MTD
#--------------------------------------------------------------------------------------------------#
# if sim=="NO" (standard option, if code not run within simulation function), determine boundaries #
#--------------------------------------------------------------------------------------------------"
if (sim=="NO"){
if (! design=="3+3"){
DLT_STOP<-get_Elim_rules(Nmax=N,phi=phi,maxprob=maxtox)$c_STOP
}
if (design=="Keyboard"){
key_dec <-get_Keyboard_rules(N=N,phi=phi,halfkey=halfkey)
}
if (design=="i3+3"){
i33_dec <-get_i33_rules(N=N,phi=phi,phi1=phi1,phi2=phi2)
}
if (design=="BOIN"){
BOIN_thres <- get_BOIN_rules(phi=phi,phi1=phi1,phi2=phi2)
lambda_e <- BOIN_thres$lambda_e
lambda_d <- BOIN_thres$lambda_d
}
}
if (sim=="YES"){
hist<-0
if (! design=="3+3"){
DLT_STOP<-env$DLT_STOP
}
if (design=="Keyboard"){
key_dec <-env$key_dec
}
if (design=="i3+3"){
i33_dec <-env$i33_dec
}
if (design=="BOIN"){
lambda_e<-env$lambda_e
lambda_d<-env$lambda_d
}
}
#-----------------------------------------#
# First dose: Accelerated Titration Phase #
#-----------------------------------------#
if (acc_tit==1){
repeat{
n <- n+1
dlt <- rbinom(n=1,size=1,truerate[dose_inv])
npt[dose_inv] <- 1
ndlt[dose_inv] <- dlt
if (dlt==0){
if (dose_inv!= ndose & dose_inv!=dose_no_titr){dose_inv <- dose_inv + 1} else # Stay in accelerated titration phase
if (dose_inv== ndose | dose_inv==dose_no_titr){break} # Go to cohort-wise dose-escalation algorithm
}
if (hist==1){result_list<-append(result_list,list(cbind(dose,npt,ndlt)))}
if (dlt==1){break} # Go to cohort-wise dose-escalation algorithm
}
}
#print(cbind(dose,npt,ndlt,MTD)) # for debugging
#------------#
# 3+3 design #
#------------#
if (design=="3+3"){
repeat{
#--------------------#
# Get data by cohort #
#--------------------#
if (npt[dose_inv]==1){ # scenario of accelerated titration with DLT or no DLT at last dose
npt[dose_inv] <- npt[dose_inv] + 2
ndlt[dose_inv] <- ndlt[dose_inv] + rbinom(n=1,size=2,truerate[dose_inv])
} else {
npt[dose_inv] <- npt[dose_inv] + 3
ndlt[dose_inv] <- ndlt[dose_inv] + rbinom(n=1,size=3,truerate[dose_inv])
}
if (hist==1){result_list<-append(result_list,list(cbind(dose,npt,ndlt)))}
p_hat<- ndlt/npt
#---------------------------------------------------#
# Get correct decision for dose under investigation #
#---------------------------------------------------#
phat<-p_hat[dose_inv]
# phat<1/3
#---------
if (phat< 1/3 & (dose_inv+1) %in% dose_in_the_running) {dose_inv <- dose_inv+1} else
if (phat< 1/3 & ! (dose_inv+1) %in% dose_in_the_running) {
if (npt[dose_inv]<=3) {dose_inv <- dose_inv} else
if (npt[dose_inv]==6) {mtd<- dose_inv
MTD[dose_inv]<-1
break}
} else
# phat>1/3
#---------
if (phat> 1/3 & (dose_inv-1) %in% dose_in_the_running){
if (npt[dose_inv-1]<=3) {dose_in_the_running <- dose_in_the_running[dose_in_the_running < dose_inv]
dose_inv <- dose_inv-1} else
if (npt[dose_inv-1]==6) {mtd<- dose_inv-1
MTD[dose_inv-1]<-1
break}
} else
if (phat> 1/3 & ! (dose_inv-1) %in% dose_in_the_running) {break} else
# phat=1/3
#---------
if (phat==1/3 & npt[dose_inv]==3) {dose_inv <- dose_inv} else
if (phat==1/3 & npt[dose_inv]==6 & (dose_inv-1) %in% dose_in_the_running){
if (npt[dose_inv-1]<=3){dose_in_the_running <- dose_in_the_running[dose_in_the_running < dose_inv]
dose_inv <- dose_inv-1} else
if (npt[dose_inv-1]==6) {mtd<- dose_inv-1
MTD[dose_inv-1]<-1
break}
} else
if (phat==1/3 & npt[dose_inv]==6 & ! (dose_inv-1) %in% dose_in_the_running) {break}
#print(cbind(dose,npt,ndlt,MTD,p_hat)) # for debugging
}
#cbind(dose,npt,ndlt,MTD,p_hat)
}
#--------------#
# BOIN design #
#--------------#
if (design=="BOIN"){
while ((N-n)>= cohortsize) {
#--------------------#
# Get data by cohort #
#--------------------#
if (npt[dose_inv]==1){ # scenario of accelerated titration with DLT or no DLT at last dose
npt[dose_inv] <- npt[dose_inv] + cohortsize-1
ndlt[dose_inv] <- ndlt[dose_inv] + rbinom(n=1,size=(cohortsize-1),truerate[dose_inv])
n <- n+(cohortsize-1)
} else {
npt[dose_inv] <- npt[dose_inv] + cohortsize
ndlt[dose_inv] <- ndlt[dose_inv] + rbinom(n=1,size=(cohortsize ),truerate[dose_inv])
n <- n+(cohortsize)
}
if (hist==1){result_list<-append(result_list,list(cbind(dose,npt,ndlt)))}
p_hat<- ndlt/npt
#----------------------------------------------------------------------------#
# Check if current dose and doses above should be eliminated due to toxicity #
#----------------------------------------------------------------------------#
if(ndlt[dose_inv]>=DLT_STOP[npt[dose_inv]]){dose_in_the_running <- dose_in_the_running[dose_in_the_running < dose_inv] }
#---------------------------------------------------#
# Get correct decision for dose under investigation #
#---------------------------------------------------#
phat <-p_hat[dose_inv]
if (phat<=lambda_e) {decision<-"E"} else # Escalate
if (phat>=lambda_d) {decision<-"D"} else # De-escalate
if (lambda_e<phat & phat<lambda_d){decision<-"R"} # Retain
if (!is.null(BOIN_add33_rule)){
if (BOIN_add33_rule==T & npt[dose_inv]==3){if (phat==1/3){decision<-"R"}}
}
if (decision=="E") {
if ( (dose_inv+1) %in% dose_in_the_running){dose_inv <- dose_inv+1} else
if (!(dose_inv+1) %in% dose_in_the_running){decision<-"R" } # If no higher dose available: RETAIN DOSE
}
if (decision=="D") {
if ( (dose_inv-1) %in% dose_in_the_running){dose_inv <- dose_inv-1} else
if (!(dose_inv-1) %in% dose_in_the_running){break } # If no lower dose available: STOP
}
if (decision=="R") {
if ( dose_inv %in% dose_in_the_running){dose_inv <- dose_inv} else {break}
}
if (maxNdec=="STAY"){
if (npt[dose_inv]>=maxN & decision=="R"){break} # if number of patients on "dose under investigation">= max specified (regardless decision E/D/R): STOP
}
if (maxNdec=="ANY"){
if (npt[dose_inv]>=maxN){break} # if number of patients on "dose under investigation">= max specified (regardless decision E/D/R): STOP
}
}
#----------------------------------------------------------------#
# Calculate isotonic estimate for p(DLT) (function BOIN package) #
#----------------------------------------------------------------#
if (MTD_safer==T){
p_iso[npt>0] <- as.numeric(as.character(BOIN::select.mtd(target=phi, npts=npt[npt>0], ntox=ndlt[npt>0])$p_est$phat)) # Get isotonic estimates
# which lower than lambda_d
safe_choice <- which(p_iso<=lambda_d)
# of those, which closest to phi
if (length(safe_choice)!=0){
distance<- abs(p_iso[safe_choice]-phi)
mtd<-max(which(distance==min(distance))) # If more than two (under and above MTD), take highest.
}
if (length(mtd)!=0 & ! is.na(mtd)){
MTD[mtd]<-1 # MTD: isotonic estimate < lambda_e
}
} else {
p_iso[npt>0] <- as.numeric(as.character(BOIN::select.mtd(target=phi, npts=npt[npt>0], ntox=ndlt[npt>0])$p_est$phat)) # Get isotonic estimates
mtd<-BOIN::select.mtd(target=phi, npts=npt[npt>0], ntox=ndlt[npt>0], cutoff.eli=maxtox)$MTD # Note if no MTD, function puts it at '99'
if (!(mtd==99)){MTD[which(npt>0)[mtd]]<-1} else
if ( (mtd==99)){mtd=NA}
}
#cbind(dose,npt,ndlt,p_hat,p_iso,MTD)
}
#--------------------------#
# Keyboard or i3+3 design #
#--------------------------#
if (design=="Keyboard" | design=="i3+3"){
if (design=="Keyboard"){dec<-key_dec}
if (design=="i3+3") {dec<-i33_dec}
while ((N-n)>= cohortsize) {
#--------------------#
# Get data by cohort #
#--------------------#
if (npt[dose_inv]==1){ # scenario of accelerated titration with DLT or no DLT at last dose
npt[dose_inv] <- npt[dose_inv] + cohortsize-1
ndlt[dose_inv] <- ndlt[dose_inv] + rbinom(n=1,size=(cohortsize-1),truerate[dose_inv])
n <- n+(cohortsize-1)
} else {
npt[dose_inv] <- npt[dose_inv] + cohortsize
ndlt[dose_inv] <- ndlt[dose_inv] + rbinom(n=1,size=(cohortsize ),truerate[dose_inv])
n <- n+(cohortsize)
}
if (hist==1){result_list<-append(result_list,list(cbind(dose,npt,ndlt)))}
p_hat<- ndlt/npt # not needed for algorithm
#----------------------------------------------------------------------------#
# Check if current dose and doses above should be eliminated due to toxicity #
#----------------------------------------------------------------------------#
if(ndlt[dose_inv]>=DLT_STOP[dose_inv]){dose_in_the_running <- dose_in_the_running[dose_in_the_running < dose_inv] }
#---------------------------------------------------#
# Get correct decision for dose under investigation #
#---------------------------------------------------#
decision<-dec[dec$n==npt[dose_inv] & dec$x==ndlt[dose_inv],"decision"]
if (decision=="E") {
if ( (dose_inv+1) %in% dose_in_the_running){dose_inv <- dose_inv+1} else
if (!(dose_inv+1) %in% dose_in_the_running){decision<-"R" } # If no higher dose available: RETAIN DOSE
}
if (decision=="D") {
if ( (dose_inv-1) %in% dose_in_the_running){dose_inv <- dose_inv-1} else
if (!(dose_inv-1) %in% dose_in_the_running){break } # If no lower dose available: STOP
}
if (decision=="R") {
if ( dose_inv %in% dose_in_the_running){dose_inv <- dose_inv} else {break}
}
if (maxNdec=="STAY"){
if (npt[dose_inv]>=maxN & decision=="R"){break} # if number of patients on "dose under investigation">= max specified (regardless decision E/D/R): STOP
}
if (maxNdec=="ANY"){
if (npt[dose_inv]>=maxN){break} # if number of patients on "dose under investigation">= max specified (regardless decision E/D/R): STOP
}
}
#----------------------------------------------------------------#
# Calculate isotonic estimate for p(DLT) (function BOIN package) #
#----------------------------------------------------------------#
if (MTD_safer==T){
p_iso[npt>0] <- as.numeric(as.character(BOIN::select.mtd(target=phi, npts=npt[npt>0], ntox=ndlt[npt>0])$p_est$phat)) # Get isotonic estimates
# which lower than lambda_d
if (design=="Keyboard"){safe_choice <- which(p_iso<phi+halfkey)}
if (design=="i3+3") {safe_choice <- which(p_iso<phi2)}
# of those, which closest to phi
if (length(safe_choice)!=0){
distance<- abs(p_iso[safe_choice]-phi)
mtd<-max(which(distance==min(distance))) # If more than two (under and above MTD), take highest.
}
if (length(mtd)!=0 & ! is.na(mtd)){
MTD[mtd]<-1 # MTD: isotonic estimate < lambda_e
}
} else {
p_iso[npt>0] <- as.numeric(as.character(BOIN::select.mtd(target=phi, npts=npt[npt>0], ntox=ndlt[npt>0])$p_est$phat)) # Get isotonic estimates
mtd<-BOIN::select.mtd(target=phi, npts=npt[npt>0], ntox=ndlt[npt>0], cutoff.eli=maxtox)$MTD # Note if no MTD, function puts it at '99'
if (!(mtd==99)){MTD[which(npt>0)[mtd]]<-1} else
if ( (mtd==99)){mtd=NA}
}
#cbind(dose,npt,ndlt,p_hat,p_iso,MTD)
}
#----------------#
# Return results #
#----------------#
result<-as.data.frame(cbind(dose,npt,ndlt,p_hat,p_iso,MTD))
if (hist==1){return(list(result,result_list))} else {
return(result)}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.