User Guide: Dynamic Borrowing with psborrow

Time-to-event

The second part involves generating time-to-events from either a piece-wise exponential distribution or weibull distribution. In the example below, we assume the survival time follows the piece-wise exponential distirbution such that baseline hazard rate is constant in months 1-2, and from month 3 onward. The baseline hazard rate (lambdaC) for patients with all covariate equaling 0 in the internal control arm in these time periods are 0.0135 and 0.2, respectively. The baseline hazard rate for the internal treatment and external control arm are the product of lambdaC and HR, driftHR, respectively.

The calculation below details how we draw event time for each patient following an exponential distribution. For each patient, we calculate a specific \(\lambda\) value based on the trial and treatment indicator, along with the covariate values for this patient.

\[ \begin{aligned} f(t) &= v\lambda t^{v-1} e^{-\lambda t^v} \\ F(t) &= 1 - e^{-\lambda t^v} \\ \lambda(t) &= \lambda v t ^{v-1} \\ \\ v &\text{: shape parameter (user-specified)}\\ \lambda &\text{: rate parameter}\\ \lambda_c &\text{: baseline hazard rate at time 0 for internal control patient with covariates } X = 0 \\ \\ \text{Alternative parameterizations}\\ f(t) &= \frac{v}{k} \left( \frac{t}{k} \right) ^{v-1} e^{-\frac{t}{k}^v} \\ &= v \left( \frac{1}{k} \right) ^{v} t^{v-1} e^{- \left(\frac{t}{k} \right) ^v} \\ \lambda &= \left( \frac{1}{k} \right) ^ v\\ \lambda(t) &= \lambda v t ^{v-1} \\ log(k) &= - \frac{1}{v} log(\lambda) \\ k &\text{: scale parameter (input for rweibull in R, one k for each patient)}\\ \\ \text{Internal control}\\ \lambda &= \lambda_c \times e^{X \beta}\\ \lambda(t) &= \lambda_c \times e^{X \beta} \times v t ^{v-1} \\ log(k) &= - \frac{1}{v} \Bigl[ log (\lambda_c) + X \beta \Bigr] \\ \\ \text{Internal treatment}\\ \lambda &= \lambda_c \times \text{HR} \times e^{X \beta} \\ \lambda(t) & = \lambda_c \times \text{HR} \times e^{X \beta} \times v t ^{v-1}\\ log(k) &= - \frac{1}{v} \Bigl[ log(\lambda_c) + log(\text{HR}) + e^{X \beta} \Bigr] \\ \\ \text{External control} \\ \lambda &= \lambda_c \times \text{driftHR} \times e^{X \beta} \\ \lambda(t) &= \lambda_c \times \text{driftHR} \times e^{X \beta} \times v t ^{v-1} \\ log(k) &= - \frac{1}{v} \Bigl[ log(\lambda_c) + log(\text{driftHR}) + e^{X \beta} \Bigr] \\ \end{aligned} \]

Since the covariate data we generated above in sample_cov contains 8 covariates, we also want to make sure the length of beta is 8.

evt <- set_event(
  event = "pwexp",
  lambdaC = c(0.0135, 0.02),
  t_itv = 2,
  beta = c(
    rep(0.2, 3),
    0.3,
    rep(0.5, 2),
    rep(2, 2)
  )
)

Additionally, users can choose to expand the covariates when calculating the hazard rate. In the example below, the hazard rate for patient i from the internal control will be calculated as below. Same as internal treatment and external control arms.

\[ \begin{aligned} \lambda_i = \lambda_C \times e^{0.2 \times x_3 + 0.3 \times (x_1 \times x_3) + 0.5 \times (x_2 + x_3) + 2 \times x_1^3} \end{aligned} \]

set_event(
  event = "pwexp",
  lambdaC = c(0.0135, 0.02),
  t_itv = 2,
  beta = c(0.2, 0.3, 0.5, 2),
  change = list(
    c("cov1", "*", "cov3"),
    c("cov2", "+", "cov3"),
    c("cov1", "^", "3")
  ),
  keep = "cov3"
)

Time-to-events can also follow a weibull distribution. Under this setting, we are expected to enter the shape of the weibull distribution which is 0.9 in this case. We let beta = 0.5 to be the coefficients for all covariates. Users can use ?set_event to see more options.

set_event(
  event = "weibull",
  shape = 0.9,
  lambdaC = 0.0135,
  beta = 0.5
)

For Weibull distribution, the goal is to find the hazard rate \(\lambda(t)\) for patient with covariates \(X\).

\[ \begin{aligned} f(t) &= \lambda e^{-\lambda t} \\ F(t) &= 1 - e^{-\lambda t} \\ \lambda(t) &= \lambda \\ \\ \lambda_c &\text{: baseline constant hazard rate for internal control patient with covariates } X = 0 \\ \text{HR} &\text{: hazard ratio between internal treatment vs control arm (proportional hazards)} \\ \text{driftHR} &\text{: hazard ratio between external vs internal control arm (proportional hazards)} \\ \text{I}(\text{trt}) &\text{: treatment arm indicator} \\ \text{Internal control}\\ \lambda &= \lambda_c \times e^{log(\text{HR}) \times \text{I}(\text{trt}) + X \beta} \\ & = \lambda_c \times e^0 \times e^{X_{i} \beta} \\ & = \lambda_c \times e^{X \beta} \\ \text{Internal treatment}\\ \lambda &= \lambda_c \times e^{log(\text{HR}) \times \text{I}(\text{trt}) + X \beta} \\ & = \lambda_c \times e^{log(\text{HR}) + X \beta} \\ & = \lambda_c \times \text{HR} \times e^{X \beta} \\ \text{External control} \\ \lambda &= \lambda_c \times \text{driftHR} \times e^{X \beta} \\ \text{Reparameterization for rweibull} \\ f(t) &= \lambda e^{-\lambda t} \\ F(t) &= 1 - e^{-\lambda t} \\ \lambda(t) &= \lambda \end{aligned} \]

Enrollment pattern, drop-out, analysis start time

Next, we specify the parameters for enrollment time, drop-out pattern and analysis start time. Both enrollment time and time-to-exit follow piece-wise exponential distribution.

The parameters are set separately for internal and external trials. For internal trial, we recruit 2 patients per month in the first 5 months, 3 patients from month 6 to month 1, and 16 patients per month after month 10. We assume dropouts at a rate of etaC = 1 per month for the control arm in the first 2 month and the rate increases to 1.2 afterwards.

The analysis is set to start at 45 months after the study begins (enrollment of the first patient) and a warning message will show if any patient enrolled after this cut-off date and they will be excluded from any future analysis.

c_int <- set_clin(
  gamma = c(2, 3, 16),
  e_itv = c(5, 10),
  CCOD = "fixed-first",
  CCOD_t = 45,
  etaC = c(0.02, 0.03),
  etaE = c(0.2, 0.3),
  d_itv = 2
)

For the external trial, both the recruitment rate gamma and the drop-out rate etaC are assumed to be constant. We can skip etaE since the external trial contains data on control arm only. Note that in this case, we do not need to specify e_itv or d_itv. The cut-off date CCOD = "event" means the analysis starts once 150 events have been observed.

c_ext = set_clin(
  gamma = 10,
  CCOD = "event",
  CCOD_t = 150,
  etaC = 0.05
)

With all the information provided above, we simulate enrollment pattern, time-to-events and time-to-dropout using function simu_time(). User can also provide a seed value here too.

sample_time <- simu_time(
  dt = sample_cov,
  eventObj = evt,
  clinInt = c_int,
  clinExt = c_ext,
  seed = 47
)

Example of one simulated dataset with time-to-events information from output sample_time:

head(sample_time[[1]], 5)
#>      driftHR   HR ext trt cov1 cov2        cov3 cov4 cov5       cov6       cov7
#> [1,]       1 0.67   0   1    0    0  1.69867515    0    0  1.4069318  1.4069318
#> [2,]       1 0.67   0   1    1    1 -0.42546841    1    1 -0.3578853 -0.3578853
#> [3,]       1 0.67   0   1    1    1 -0.04927021    0    0  0.1370222  0.1370222
#> [4,]       1 0.67   0   1    1    1 -0.42705058    0    0  0.3709728  0.3709728
#> [5,]       1 0.67   0   1    0    1 -0.10903808    0    0  0.2500303  0.2500303
#>            cov8      time cnsr
#> [1,]  1.4069318 0.1019925    0
#> [2,] -0.3578853 6.4700163    1
#> [3,]  0.1370222 3.5066896    1
#> [4,]  0.3709728 1.5454039    1
#> [5,]  0.2500303 6.8677445    0