---
title: "dtpcrm: Dose Transition Pathways with Continual Reassessment Method"
author: "Christina Yap, Daniel Slade"
date: "`r Sys.Date()`"
output: 
    rmarkdown::html_vignette:
    fig_caption: yes
vignette: >
  %\VignetteIndexEntry{dtpcrm: Dose Transition Pathways with Continual Reassessment Method}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---


```{r setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

# Description

The dtpcrm package presents the design and conduct of a Continual Reassessment Method (CRM) trial with the use of the Dose Transition Pathways (DTP). Building on the R dfcrm package [1], the dtpcrm package provides the DTP to project in advance the doses recommended by a model-based dose-finding design for subsequent subjects (stay, escalate, de-escalate or stop early) using all the accumulated toxicity information [2]. DTP can be used as a design calibration tool and an operational/trial conduct tool. It can be presented as a table or as a flow diagram. The dtpcrm package also provides the modified CRM [3] and Time to event CRM (TITE-CRM) [4] with added practical modifications to allow stopping early when:  

-    there is sufficient evidence that the lowest dose is too toxic and/or
-    there is a sufficient number of subjects dosed at the maximum tolerated dose (MTD)

Here, we will illustrate the use of the DTP as  

(1) a design calibration tool
(2) a trial conduct tool 

when designing a CRM, using the Viola trial as a case study [2,5].


##  Getting started
In order to run the code, you have to first install the package and load it with

```{r}
library(dtpcrm)
```

--------

# Design a CRM trial with DTP as a calibration tool


## CRM Design Considerations

We need to input the parameters required for design considerations for CRM (see Figure 1 in Yap et al [2]). Those could be grouped under the three broad categories of Clinical Parameters, Model Specification Parameters and Practical Considerations. 


### 1. Clinical Parameters

```{r}
number.doses= 7

start.dose.level = 3

max.sample.size  = 21

target.DLT = 0.2

cohort.size = 3

```

### 2. Model Specification Parameters

A one-stage one-parameter Bayesian empiric model is used where the slope parameter is assumed to be normally distributed with mean 0 with a specified prior variance. 

The prior DLT rates for the 7 doses and prior variance are specified as:
```{r}
prior.DLT = c(0.03, 0.07, 0.12, 0.20, 0.30, 0.40, 0.52)
prior.var = 0.75
```


### 3. Practical Considerations  

Using a model-based dose-escalation design, dose recommendations to escalate, de-escalate or stay at the current dose for subsequent cohorts of subjects are directed by the specified model, prior and observed data. In certain cases, those purely model-based recommendations may cause safety concerns, and it may be desirable to introduce practical constraints. 

Practical considerations for Viola include:

- No skipping of doses in escalation (default option)
- Allow skipping of doses in de-escalation. 
- No dose escalation if the overall observed rate of DLT at the current dose level is above the target DLT rate (default option)

```{r}
no_skip_esc = TRUE
no_skip_deesc = FALSE
global_coherent_esc = TRUE
```

- Stop the trial early for toxicity if there is strong evidence that the lowest dose is too toxic based on a Bayesian safety early criterion: 

p(true DLT rate at lowest dose > target DLT rate + $x%$ | current
observed data and any relevant prior information) > $y$

```{r}
y = 0.72
x = 0.1
tox_lim = target.DLT + x
```

- Stop the trial early if there are 12 subjects at the recommended MTD. 

```{r}
n.mtd= 12     
```


Stopping early criterion (for safety and sufficient number of subjects at MTD)

```{r}
stop_func <- function(x) {
  y = stop_for_excess_toxicity_empiric(x, tox_lim = tox_lim, prob_cert = y, dose = 1, nsamps = 100000)
  if(y$stop){
    x <- y
  } else {x = stop_for_consensus_reached(x, req_at_mtd = n.mtd)}
}
```


## Simulations: assess operating characteristics of the proposed CRM design

Simulated operating characteristics (OCs) are an important tool to assess the performance of a proposed adaptive design for informed decision making. 

Let's examine three scenarios, S1, S2 and S3 with assumed true DLT rates at each dose as defined below. S1 and S2 presents scenarios where the true MTD are at d(1) and d(3) respectively. S3 presents a setting whereby all doses are too toxic.


```{r}
S1=c(0.03, 0.05, 0.07, 0.20, 0.36, 0.45, 0.55)
S2=c(0.005, 0.01, 0.03, 0.05, 0.07, 0.20, 0.36)
S3=c(0.45, 0.55, 0.65, 0.7, 0.75, 0.80, 0.85)
```

```{r, echo=FALSE}
scenarios<-rbind.data.frame(S1,
                            S2,
                            S3)
dose.labels= seq(-2,4,by=1)
dimnames(scenarios)[[2]]<-paste("d(", dose.labels,")", sep="")
scenarios<-cbind.data.frame("Scenario" = c("S1", "S2", "S3"), scenarios)
knitr::kable(scenarios)
```

It is useful to bear in mind that simulaton results are unavoidably subject to random variation. Hence in order to obtain adequately stable results, 5000 simulations per setting would be recommended for the final design. For illustrative purposes, we have used 50 here. 

```{r}
set.seed(1000)
No.of.simulations = 50
```

Scenario 1: True MTD at d(1),
```{r message=FALSE, warning=FALSE, results= "hide"}
sim1_1 <- applied_crm_sim(true_tox = S1, prior = prior.DLT, target = target.DLT, max_sample_size=max.sample.size, first_dose=start.dose.level, num_sims = No.of.simulations, cohort_size = cohort.size, dose_func = applied_crm, scale = sqrt(prior.var), stop_func = stop_func)
```

Extracting outputs on probability of selecting each dose as the MTD and the number of subjects allocated to each dose, as well as 
stopping probabilty for toxicity or sufficient number at MTD,

```{r}
stop.fun<-function(x, matchtxt) {
   prob<-round(ifelse( !is.null(grep(matchtxt, names(x))>0) , x[grep(matchtxt, names(x))], 0),2)
   prob[is.na(prob)]<-0
   return(prob)
}

out<- cbind.data.frame(
           rbind.data.frame("S1 True Tox" = sim1_1$summary$true_tox,
                            "S1 Prob Select"= sim1_1$summary$mtd,
                            "S1 No.Subjects"= sim1_1$summary$doses_given),
      Stop.Tox  =   c("", stop.fun(sim1_1$summary$prob_stop, "Tox"), ""),
      Stop.nMTD =   c("", stop.fun(sim1_1$summary$prob_stop, "Con"), "")
)
dimnames(out)[[2]][1:length(dose.labels)]<-paste("d(", dose.labels,")", sep="")
dimnames(out)[[1]] <- paste(c(paste("S1", " True DLT rate", sep = ""), paste("S1",  " Selection Probability", sep = ""), paste("S1",  " Mean Number of Subjects", sep = "")))

```

Repeat the procedure for Scenarios 2 and 3.
```{r, echo= FALSE, results="hide"}
sim1_2 <- applied_crm_sim(true_tox = S2, prior = prior.DLT, target = target.DLT, max_sample_size=max.sample.size, first_dose=start.dose.level, num_sims = No.of.simulations, cohort_size = cohort.size, dose_func = applied_crm, scale = sqrt(prior.var), stop_func = stop_func)
```

```{r, echo= FALSE, results="hide"}
sim1_3 <- applied_crm_sim(true_tox = S3, prior = prior.DLT, target = target.DLT, max_sample_size=max.sample.size, first_dose=start.dose.level, num_sims = No.of.simulations, cohort_size = cohort.size, dose_func = applied_crm, scale = sqrt(prior.var), stop_func = stop_func)
```


```{r, echo=FALSE}
out2<- cbind.data.frame(
           rbind.data.frame("S2 True Tox" = sim1_2$summary$true_tox,
                            "S2 Prob Select"= sim1_2$summary$mtd,
                            "S2 No.Subjects"= sim1_2$summary$doses_given),
      Stop.Tox  =   c("", stop.fun(sim1_2$summary$prob_stop, "Tox"), ""),
      Stop.nMTD =   c("", stop.fun(sim1_2$summary$prob_stop, "Con"), "")
)
dimnames(out2)[[2]][1:length(dose.labels)]<-paste("d(", dose.labels,")", sep="")
dimnames(out2)[[1]] <- paste(c(paste("S2", " True DLT rate", sep = ""), paste("S2", " Selection Probability", sep = ""), paste("S2", " Mean Number of Subjects", sep = "")))

out3<- cbind.data.frame(
           rbind.data.frame("S3 True Tox" = sim1_3$summary$true_tox,
                            "S3 Prob Select"= sim1_3$summary$mtd,
                            "S3 No.Subjects"= sim1_3$summary$doses_given),
      Stop.Tox  =   c("", stop.fun(sim1_3$summary$prob_stop, "Tox"), ""),
      Stop.nMTD =   c("", stop.fun(sim1_3$summary$prob_stop, "Con"), "")
)
dimnames(out3)[[2]][1:length(dose.labels)]<-paste0("d(", dose.labels,")")
dimnames(out3)[[1]] <- paste(c(paste("S3",  " True DLT rate", sep = ""), paste("S3",  " Selection Probability", sep = ""), paste("S3",  " Mean Number of Subjects", sep = "")))


out<- rbind.data.frame(
      "Prior DLT" = c(prior.DLT, "", ""),
      out, 
      out2,
      out3
)      
```

### Simulation Results 

The table below presents the simulation results for the selected scenarios, with prior and true DLT rates, probability of selecting a dose as the MTD and the number of subjects allocated to each dose. It also presents the probability of stopping early due to toxicity or having sufficient number of subjects at the MTD.  
```{r, echo=FALSE}
knitr::kable(out)
```



## DTP for initial cohorts

Using the specific parameters above, we can produce the DTP for the initial 4 cohorts of the Viola trial.

```{r}
start.dose.level<-3   #(eg. 1,2,3 etc)

viola_dtp <- calculate_dtps(next_dose = start.dose.level, cohort_sizes = c(cohort.size, cohort.size, cohort.size), dose_func = applied_crm, prior = prior.DLT, target = target.DLT,
stop_func = stop_func, scale = sqrt(prior.var),
no_skip_esc = no_skip_esc, no_skip_deesc = no_skip_deesc,
global_coherent_esc = global_coherent_esc)
```

In Viola, there are 7 dose levels where the starting dose level is 3. For display purposes, we would use the dose-labels as used in Viola, which range from d(-2) to d(4), with starting dose at d(0).

```{r}
# Using dose labels                         
viola_dtp[seq(1, ncol(viola_dtp), by=2)] <- viola_dtp[seq(1, ncol(viola_dtp), by=2)] - start.dose.level

# Indicate when stopping early occurs
indSTOP<-is.na(viola_dtp[seq(1,ncol(viola_dtp), by=2)])
viola_dtp.pretty<-viola_dtp
viola_dtp.pretty[seq(1,ncol(viola_dtp.pretty), by=2)][indSTOP]<-"STOP"
viola_dtp.pretty<-cbind.data.frame(1:nrow(viola_dtp.pretty), viola_dtp.pretty)
dimnames(viola_dtp.pretty)[[2]]<-c('Pathway', 'C1 Dose', 'C1 No.DLT',  'C2 Dose', 'C2 No.DLT','C3 Dose', 'C3 No.DLT', 'C4 Dose')

```




### Displaying DTP in a table format
```{r, echo=FALSE}
knitr::kable(viola_dtp.pretty)
```

The term "NA" denotes not applicable. 

###  Displaying DTP as a flow diagram

Another way to display DTP is via a flow diagram. To obtain DTP flow diagram for Cohorts 1-3,

```{r, fig.width=8, fig.height=6}
dtpflow(viola_dtp)
```

One can utilise the DTP for the initial cohorts to discuss with the investigators on the working of the proposed CRM design to check agreement and highlight unanticpated pathways. One can further calibrate the model by fine-tuning the model parameters to achieve DTP that are agreeable by the investigators. In particular, it would be useful to check that the investigators are in agreement on when a trial should be stopped early due to toxicity if excessive DLTs are observed in the initial cohorts. One should then evaluate average model performance of the fine-tuned model by re-running simulations to ensure that it is still favorable (see Figure 2 in [2]). When it comes to choosing a **_desirable_** model-based design,the final recommended model should be one that not only has promising operating characteristics but achieves DTP which the investigators are happy to adopt in practice [2,5].

In the case for Viola, the investigators were keen to only stop the trial early if they saw 2 or 3 DLTs out of the first 3 subjects at the lowest dose, d(-2), which would indicate that all the doses were too toxic. This guided the choice of $y$ to be around 0.72. 

--------

# Trial Conduct with DTP

No DLTs were observed in the first 2 cohorts who were dosed at dose levels 3 and 4. 

```{r}
DLT.outcomes<-c(0,0,0,0,0,0)
dose.level = c(3,3,3,4,4,4)
```

## Model Updates with accumulated information

```{r}
C2 <- applied_crm(prior=prior.DLT, target = target.DLT, tox = DLT.outcomes, level = dose.level, 
                              stop_func = stop_func, 
                              no_skip_esc = no_skip_esc, 
                              no_skip_deesc = no_skip_deesc,
                              global_coherent_esc = global_coherent_esc,
                              scale = sqrt(prior.var))
results<-cbind.data.frame(Dose=paste0("d(", dose.labels, ")"), summary_crm(C2)[2:ncol(summary_crm(C2))])
dimnames(results)[[2]][2:ncol(results)]<-c("Prior DLT", "No.subjects", "No.DLT", "Posterior DLT (90% PI)")
```

Recommended dose for next cohort is Dose level 
```{r} 
C2$mtd 
```


Or dose label
```{r}
paste0("d(", C2$mtd-start.dose.level, ")"  )
```

The table below presents the CRM updates with prior DLT rates, number of evaluable subjects, number of DLTs and posterior DLT rates with 90% probability interval at each of the dose level.

```{r, echo=FALSE}
knitr::kable(results)
```


## Update of Dose Toxicity Curves


```{r, fig.width=6, fig.height=6}
plot_crm(C2, dose_labels = paste0("d(",dose.labels,")"), cohort_sizes = c(3,3), dose_func =applied_crm, 
        no_skip_esc = no_skip_esc, 
        no_skip_deesc = no_skip_deesc,
        global_coherent_esc = global_coherent_esc, 
        scale = sqrt(prior.var),
        ylim = c(0, 0.7), lwd = 2, cex.axis = 1.5, cex.lab = 1.4, cex = 1, cohort.last = T
        )
```

This figure presents the updates of the dose toxicity curves with evolving estimates of DLT rates at each dose, estimated at the end of each cohort. The black curve represents the prior estimates of DLT rates at each dose.


## Update DTP for Subsequent Cohorts

DTP for subsequent cohorts can be updated regularly. For instance, we can update after Cohorts 1 and 2, giving the DTP for Cohorts 3-6. 

```{r message=FALSE, warning=FALSE, results="hide"}
dtp_c3 <- calculate_dtps(next_dose = 5, cohort_sizes = c(3, 3, 3),
                         prev_tox = DLT.outcomes, prev_dose = dose.level,
                         prior = prior.DLT, target = target.DLT, 
                         no_skip_esc = no_skip_esc, 
                         no_skip_deesc = no_skip_deesc,
                         global_coherent_esc = global_coherent_esc, 
                         stop_func = stop_func, 
                         scale = sqrt(prior.var))

# Using dose labels                         
dtp_c3[seq(1, ncol(dtp_c3), by=2)] <- dtp_c3[seq(1, ncol(dtp_c3), by=2)] - start.dose.level
dtp_c3_pretty<-cbind.data.frame(1:nrow(dtp_c3), dtp_c3)
dimnames(dtp_c3_pretty)[[2]]<-c('Pathway','C3 Dose', 'C3 No.DLT',  'C4 Dose', 'C4 No.DLT','C5 Dose', 'C5 No.DLT', 'C6 Dose')
```


###  DTP table for Cohorts 3 to 6 
```{r echo=FALSE, warning=FALSE}
knitr::kable(dtp_c3_pretty)
```


###  DTP flow diagram for Cohorts 3 to 5

```{r, fig.width=8, fig.height=6}
dtpflow(dtp_c3, cohort.labels = c('C3', 'C4', 'C5') )

```

--------

# References
[1] Ken Cheung (2013). dfcrm: Dose-finding by the continual reassessment method. R package version 0.2-2. https://CRAN.R-project.org/package=dfcrm

[2] Yap, C., Billingham, L. J., Cheung, Y. K., Craddock, C., & O'Quigley, J. (2017). Dose Transition Pathways: The missing link between complex dose-finding designs and simple decision-making. Clinical Cancer Research, 23(24), 7440-7447

[3] O'Quigley, J. O., Pepe, M., and Fisher, L. (1990). Continual reassessment method: A practical design for phase I clinical trials in cancer. Biometrics 46:33-48.

[4] Cheung YK, Chappell R. Sequential designs for phase I clinical trials with late onset toxicities. Biometrics 2000;56:1177–82.

[5] Craddock C, Slade D, Santo C De, Wheat R, Ferguson P, Hodgkinson A, Brock K, Cavenagh J, Ingram W, Dennis M, Malladi R, Siddique S, Mussai F and Yap C (2018) Combination lenalidomide and azacitidine: a novel salvage therapy in patients who relapse after allogeneic stem cell transplantation for acute myeloid leukemia. Journal of Clinical Oncology (2019), JCO-18.

[6] Wheeler, G.M., Mander, A.P., Bedding, A., Brock, K., Cornelius, V., Grieve, A.P., Jaki, T., Love, S.B., Weir, C.J., Yap, C. and Bond, S.J., 2019. How to design a dose-finding study using the continual reassessment method. BMC medical research methodology, 19(1), p.18.
