RARtrials

https://github.com/yayayaoyaoyao/rartrials
Admin
RARtrials is an R package for simulating popular response-adaptive randomization methods in clinical...
Tokens:97,701
Snippets:504
Trust Score:3
License:GPL-3.0
Update:1 week ago
Show doc for...
Tokens:
Raw
### Install RARtrials from GitHub

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/README.md

Install the RARtrials package from GitHub using devtools. The commented line is a prerequisite for installing devtools if it is not already installed.

```r
#install.packages('devtools')
devtools::install_github("yayayaoyaoyao/RARtrials")
```

--------------------------------

### pmax_beta example with 5 arms and upper side

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/pmax_beta.Rd

Example call for a trial with 5 arms, specifying Beta parameters for each arm and side='upper' to compute the probability that arm 1 is the best (maximal).

```R
pmax_beta(armn=5,a1=8,b1=10,a2=5,b2=19,a3=8,b3=21,
a4=6, b4=35, a5=15, b5=4, side='upper')
```

--------------------------------

### pmax_beta example with 4 arms and lower side

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/pmax_beta.Rd

Example call for a trial with 4 arms, specifying Beta parameters for each arm and side='lower' to compute the probability that arm 1 is the worst (minimal).

```R
pmax_beta(armn=4,a1=56,b1=98,a2=25,b2=70,a3=87,b3=107,
a4=106, b4=202, side='lower')
```

--------------------------------

### pmax_beta example with 3 arms and upper side

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/pmax_beta.Rd

Example call for a trial with 3 arms, specifying Beta parameters for each arm and side='upper' to compute the probability that arm 1 is the best (maximal).

```R
pmax_beta(armn=3,a1=60,b1=46,a2=55,b2=46,a3=35,b3=36,side='upper')
```

--------------------------------

### Install RARtrials from CRAN

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/README.md

Install the RARtrials package from CRAN.

```r
install.packages('RARtrials')
```

--------------------------------

### pmax_normal examples for 3, 4, and 5 arms

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/pmax_normal.Rd

Three example calls to pmax_normal with different arm counts (5, 4, 3) and side values ('upper', 'lower', 'upper'). Each call specifies the required mean and sd arguments for the given number of arms.

```r
pmax_normal(armn=5,mean1=0.8,sd1=0.2,mean2=0.5,sd2=0.1,mean3=0.8,
sd3=0.5,mean4=0.6,sd4=0.2,mean5=0.6,sd5=0.2,side='upper')
```

```r
pmax_normal(armn=4,mean1=8,sd1=2,mean2=8.5,sd2=2,mean3=8.3,
sd3=1.8,mean4=8.7,sd4=2,side='lower')
```

```r
pmax_normal(armn=3,mean1=80,sd1=20,mean2=50,sd2=10,mean3=80,
sd3=15,side='upper')
```

--------------------------------

### Integration-based posterior probability calculation with pgreater_beta and pmax_beta

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/vignettes/RARtrials.Rmd

Uses direct integration functions pgreater_beta and pmax_beta to compute the same posterior probabilities as the simulation example, but faster. Requires the same alpha and beta vectors; delta=0.1 with side='upper' for larger effects, delta=-0.1 with side='lower' for smaller effects. Results are printed as comments.

```r

#Example: Sample code using Integration
#Expect the treatment group to have a larger treatment effect compared to the control group
#Calculate results of Pr(p_k>p_{control}+\delta|data_{t-1})
pgreater_beta(a1=alpha[1],b1=beta[1],a2=alpha[2], b2=beta[2],delta=0.1,side='upper')
pgreater_beta(a1=alpha[1],b1=beta[1],a2=alpha[3], b2=beta[3],delta=0.1,side='upper')
#Return results:  0.7951487 0.3477606

#Calculate results of Pr(p_k=max\{p_1,...,p_K\}|data_{t-1})
pmax_beta(armn=3,a2=alpha[3],b2=beta[3],a3=alpha[2],b3=beta[2],a1=alpha[1],
     b1=beta[1],side='upper')
pmax_beta(armn=3,a2=alpha[1],b2=beta[1],a3=alpha[3],b3=beta[3],a1=alpha[2],
     b1=beta[2],side='upper')
pmax_beta(armn=3,a2=alpha[1],b2=beta[1],a3=alpha[2],b3=beta[2],a1=alpha[3],
     b1=beta[3],side='upper')
#Return results: 0.01796526  0.8788907  0.1031441

#Expect the treatment group to have a smaller treatment effect compared to the control group
#Calculate results of Pr(p_{control}>p_k+\delta|data_{t-1})
pgreater_beta(a1=alpha[1],b1=beta[1],a2=alpha[2],b2=beta[2],delta=-0.1,side='lower')
pgreater_beta(a1=alpha[1],b1=beta[1],a2=alpha[3],b2=beta[3],delta=-0.1,side='lower')
#Return results:  0.001093548  0.03348547

#Calculate results of Pr(p_k=min\{p_1,...,p_K\}|data_{t-1})
pmax_beta(armn=3,a2=alpha[3],b2=beta[3],a3=alpha[2],b3=beta[2],a1=alpha[1],
     b1=beta[1],side='lower')
pmax_beta(armn=3,a2=alpha[1],b2=beta[1],a3=alpha[3],b3=beta[3],a1=alpha[2],
     b1=beta[2],side='lower')
pmax_beta(armn=3,a2=alpha[1],b2=beta[1],a3=alpha[2],b3=beta[2],a1=alpha[3],
     b1=beta[3],side='lower')
#Return results: 0.7560864  0.01230027  0.2316133
```

--------------------------------

### Simulate A-Optimal Allocation Trial with sim_A_optimal_known_var

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/sim_A_optimal_known_var.Rd

Calls sim_A_optimal_known_var with patient accrual, outcome delay, enrollment rate, sample sizes, arm count, means, standard deviations, alpha, arm labels, and side. The example shows the required argument order and typical values for a two-arm trial with known variances.

```R
sim_A_optimal_known_var(
  Pats = 10,
  nMax = 100,
  TimeToOutcome = expression(distribution(sample = 1, rate = 0.5)),
  enrollrate = 1,
  N2 = 50,
  armn = 2,
  mean = c(0, 0),
  sd = c(1, 1),
  alphaa = 0.025,
  armlabel = c("Control", "Treatment"),
  side = 1
)
```

--------------------------------

### sim_dabcd_max_power

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/sim_dabcd_max_power.Rd

Documentation for the sim_dabcd_max_power function, including usage, arguments, return value, and examples.

```APIDOC
## sim_dabcd_max_power

### Description
Simulates a trial for doubly adaptive biased coin design with maximal power strategy targeting generalized Neyman allocation (2-5 arms) or generalized RSIHR allocation (2-3 arms) for binary outcomes. The function returns an object of class "dabcd" containing the final decision, Z test statistics, the simulated data set, and participants accrued for each arm.

### Usage
```r
sim_dabcd_max_power(Pats, nMax, TimeToOutcome, enrollrate, N1, N2, armn, armlabel, h, BB, type, side)
```

### Arguments
- **Pats** (numeric) - Number of patients accrued within a certain time frame.
- **nMax** (numeric) - Assumed maximum accrued number of patients with the disease in the population.
- **TimeToOutcome** (expression) - Distribution of delayed response times or a fixed delay time for responses.
- **enrollrate** (numeric) - Patient enrollment rate.
- **N1** (numeric) - First stage sample size.
- **N2** (numeric) - Second stage sample size.
- **armn** (numeric) - Number of arms (2-5 for Neyman, 2-3 for RSIHR).
- **armlabel** (vector) - Labels for the arms.
- **h** (numeric vector) - Allocation probabilities or related parameters.
- **BB** (numeric) - Biased coin parameter.
- **type** (character) - Allocation type: 'Neyman' or 'RSIHR'.
- **side** (character) - Direction of one-sided test: 'upper' or 'lower'.

### Value
Returns an object of class "dabcd", a list containing:
- Final decision based on Z test statistics (1 = selected, 0 = not selected).
- Z test statistics.
- Simulated data set with columns: participant ID, enrollment time, observed time of results, allocated arm, and participant's result.
- Participants accrued for each arm at termination time.

### Details
- For more than two arms, one-sided nominal level is alphaa divided by arm*(arm-1)/2 (Bonferroni correction).
- Methods are not smoothed and based on Hu & Zhang's formula.
- Delay mechanism parameters adapted from IntroBayesianSimulation.

### Examples
```r
sim_dabcd_max_power(Pats=10,nMax=50000,TimeToOutcome=expression(rnorm(length( vStartTime ),30, 3)),enrollrate=0.9,N1=30,N2=300,armn=3,armlabel=c(1,2,3),h=c(0.2,0.3,0.2),BB=0.1,type='Neyman',side='upper')
sim_dabcd_max_power(Pats=10,nMax=50000,TimeToOutcome=expression(rnorm(length( vStartTime ),60, 3)),enrollrate=0.1,N1=50,N2=500,armn=3,armlabel=c(1,2,3),h=c(0.2,0.3,0.3),BB=0.15,type='RSIHR',side='upper')
```
```

--------------------------------

### sim_RSIHR_optimal_known_var

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/sim_RSIHR_optimal_known_var.Rd

Documentation for the sim_RSIHR_optimal_known_var function, including usage, arguments, return value, details, examples, and references.

```APIDOC
## sim_RSIHR_optimal_known_var

### Description
Simulates a trial for continuous endpoints with known variances and fixed allocation probabilities, aiming to minimize the criteria \(\sum_{i=1}^{K}n_i\Psi_i\) with constraints \(\frac{\sigma_1^2}{n_1}+\frac{\sigma_k^2}{n_k}\leq C\).

### Usage
```r
sim_RSIHR_optimal_known_var(Pats, nMax, TimeToOutcome, enrollrate, N1, N2, armn, mean, sd, alphaa, armlabel, cc, side)
```

### Arguments
- **Pats** (numeric) - Number of patients accrued within a certain time frame.
- **nMax** (numeric) - Assumed maximum accrued number of patients with the disease in the population.
- **TimeToOutcome** (expression) - Distribution of delayed response times or a fixed delay time for responses.
- **enrollrate** (numeric) - Enrollment rate.
- **N1** (numeric) - First interim sample size.
- **N2** (numeric) - Second interim sample size.
- **armn** (numeric) - Number of arms.
- **mean** (numeric vector) - Mean for each arm.
- **sd** (numeric vector) - Standard deviation for each arm.
- **alphaa** (numeric) - One-sided nominal level for each test (with Bonferroni correction for more than two arms).
- **armlabel** (numeric vector) - Labels for each arm.
- **cc** (numeric) - Fixed constant C in the constraint.
- **side** (character) - Direction of a one-sided test, with values 'upper' or 'lower'.

### Value
Returns an object of class "RSIHRoptimal", a list containing:
- Final decision based on the Z test statistics (1 for selected, 0 for not selected).
- Z test statistics.
- The simulated data set and participants accrued for each arm at the time of termination of that group in one trial.

The simulated data set includes 5 columns: participant ID number, enrollment time, observed time of results, allocated arm, and participants' result.

### Details
This function aims to minimize the criteria \(\sum_{i=1}^{K}n_i\Psi_i\) with constraints \(\frac{\sigma_1^2}{n_1}+\frac{\sigma_k^2}{n_k}\leq C\), where \(k=2,...,K\) for some fixed C. It is equivalent to generalized RSIHR allocation for continuous endpoints with known variances. With more than two arms, the one-sided nominal level of each test is \(\alpha_a\) divided by \(arm \times (arm-1)/2\); a Bonferroni correction. Considering the delay mechanism, parameters \(Pats\), \(nMax\), and \(TimeToOutcome\) are adapted from \url{https://github.com/kwathen/IntroBayesianSimulation}.

### Examples
```r
# Run the function with delayed responses follow a normal distribution with
# a mean of 30 days and a standard deviation of 3 days under null hypothesis
# in a three-armed trial
sim_RSIHR_optimal_known_var(Pats=10,nMax=50000,TimeToOutcome=expression(
rnorm(length( vStartTime ),30, 3)),enrollrate=0.9,N1=9,N2=132,armn=3,
mean=c(9.1/100,9.1/100,9.1/100),sd=c(0.009,0.009,0.009),alphaa=0.025,
armlabel=c(1,2,3),cc=mean(c(9.1/100,9.1/100,9.1/100)),side='lower')

# Run the function with delayed responses follow a normal distribution with
# a mean of 30 days and a standard deviation of 3 days under alternative hypothesis
# in a three-armed trial
sim_RSIHR_optimal_known_var(Pats=10,nMax=50000,TimeToOutcome=expression(
rnorm(length( vStartTime ),30, 3)),enrollrate=0.9,N1=9,N2=132,armn=3,
mean=c(9.1/100,8.47/100,8.47/100),sd=c(0.009,0.009,0.009),alphaa=0.025,
armlabel=c(1,2,3),cc=mean(c(9.1/100,8.47/100,8.47/100)),side='lower')
```

### References
Biswas (2011) as cited in the RARtrials package.
```

--------------------------------

### Select a_U by calling brar_select_au_binary 2000 times

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/vignettes/RARtrials.Rmd

Runs the brar_select_au_binary function 2000 times with a fixed seed to generate the test statistic distribution under the intermediate hypothesis c(0.3,0.3). The example sets the alternative hypothesis to c(0.2,0.4) and uses the intermediate hypothesis to avoid excessively large or small power.

```R
set.seed(12345)
#Example: Select a_U by calling brar_au_binary 2000 times
simnull3<-lapply(1:2000,function(x){
  set.seed(x)
  brar_select_au_binary(Pats=10,nMax=50000,TimeToOutcome=0,enrollrate=0.9,N1=24,armn=2,
        h=c(0.3,0.3),N2=224,tp=1,armlabel=c(1, 2),blocksize=4,alpha1=1,beta1=1,
        alpha2=1,beta2=1,minstart=24,deltaa=-0.07,tpp=0,deltaa1=0.1,side='upper')     
})
```

--------------------------------

### Simulate cut-off and power for FLGI with known variances

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/vignettes/RARtrials.Rmd

R code that first selects a critical cut-off value under the null hypothesis using flgi_cut_off_known_var(), then calculates power using sim_flgi_known_var(). The example uses continuous outcomes with known variances and the 'FLGI PM' rule. It recommends running more simulations for accuracy and notes the heavy computational burden.

```r
#Example: The forward-looking Gittins Index rule applied in continuous outcomes with known 
#variances
#Select cut-off values under null hypothesis for continuous outcomes with known variances.
#The delayed responses follow a normal distribution with a mean of 30 days, sd of 3 days 
#and enrollment rate of 1. 
set.seed(123452)
stopbound5<-lapply(1:1000,function(x){flgi_cut_off_known_var(
Gittinstype='KV',df=0.995,Pats=10,nMax=50000,TimeToOutcome=0,enrollrate=1,armn=3,
noRuns2=100,N2=120,block=8,rule='FLGI PM',prior_n=rep(1,3),prior_mean=rep(0,3),
mean=c(-0.05,-0.05,-0.05),sd=c(0.346,0.346,0.346),side='upper')})
stopbound5a<-do.call(rbind,stopbound5)
sum(stopbound5a[,1]>2.074 | stopbound5a[,2]> 2.074)/1000 #0.05
#The selected cut-off value is 2.074 with an overall upper one-sided type I error of 0.05.
#It is recommended to run more simulations to obtain a more accurate cut-off value.

#Calculate power based on 1000 simulations
set.seed(123452)
sim5<-lapply(1:1000,function(x){sim_flgi_known_var(Gittinstype='KV',
df=0.995,Pats=10,nMax=50000,TimeToOutcome=0,enrollrate=1,armn=3,
noRuns2=100,N2=120,block=8,rule='FLGI PM', prior_n=rep(1,3), 
prior_mean=rep(0,3), mean=c(-0.05,0.07,0.13),sd=c(0.346,0.346,0.346),
stopbound =2.074,side='upper')})
h1decision<-t(sapply(sim5, "[[", 1))
sum(h1decision[,1]==1)/1000 #0.074
sum(h1decision[,2]==1)/1000 #0.244
sum(h1decision[,1]==1 | h1decision[,2]==1)/1000 #0.267
#Marginal power of rejecting H02 is 0.074
#Marginal power of rejecting H03 is 0.244
#Power of rejecting H02 or H03 is 0.267
```

--------------------------------

### Select a_U and estimate power with brar_select_au_unknown_var

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/vignettes/RARtrials.Rmd

Selects the futility threshold a_U by calling brar_select_au_unknown_var 2000 times under the null, then uses the chosen a_U (0.9836) in sim_brar_unknown_var to estimate power. The first block runs the null simulations and extracts test statistics, ensuring about 1% futility stops and a type I error near 0.025. The second block runs 1000 power simulations and reports the proportion of 'Superiorityfinal' decisions (0.882).

```R
set.seed(12345)
simnull4<-lapply(1:2000,function(x){
  brar_select_au_unknown_var(Pats=10,nMax=50000,TimeToOutcome=expression(
  rnorm(length( vStartTime ),30, 3)),enrollrate=0.1, N1=192,armn=3,
  N2=1920,tp=1,armlabel=c(1,2,3),blocksize=6,
  mean=c((9.1/100+8.92/100+8.92/100)/3,(9.1/100+8.92/100+8.92/100)/3,
         (9.1/100+8.92/100+8.92/100)/3),sd=c(0.009,0.009,0.009),minstart=192,
         deltaa=c(0.00077,0.00077),tpp=1,deltaa1=c(0,0),V01=1/2,a01=0.3,m01=9/100,
         b01=0.00001,side='lower')
  })  

#Obtain the data set of test statistics
simf<-list()
for (xx in 1:2000){
  if (any(simnull4[[xx]][192:1919,2]<0.01)){
    simf[[xx]]<-NA
  }  else{
    simf[[xx]]<-simnull4[[xx]][1920,2]
  }
}
simf4a<-do.call(rbind,simf)

simf<-list()
for (xx in 1:2000){
  if (any(simnull4[[xx]][192:1919,3]<0.01)){
    simf[[xx]]<-NA
  }  else{
    simf[[xx]]<-simnull4[[xx]][1920,3]
  }
}
simf4b<-do.call(rbind,simf)

#Ensure that around 1% of the trials stop for futility
sum(is.na(simf4a)) #21
sum(is.na(simf4b)) #22
#Select a_U to make sure that the overall type I error is around 0.025
sum((simf4a[,1]>0.9836| simf4b[,1]>0.9836),na.rm=T )/2000#0.025
#The selected a_U is 0.9836.
```

```R
sim4<-lapply(1:1000,function(x) {
    sim_brar_unknown_var(Pats=10,nMax=50000,TimeToOutcome=expression(rnorm(
      length(vStartTime ),30, 3)),enrollrate=0.1,N1=192,armn=3,
            a_U=c(0.9836,0.9836),N2=1920,tp=1,armlabel=c(1,2,3),blocksize=6,
            mean=c(9.1/100,8.92/100,8.92/100),sd=c(0.009,0.009,0.009),
            minstart=192,deltaa=c(0.00077,0.00077),tpp=1,deltaa1=c(0,0),
            V01=1/2,a01=0.3,m01=9/100,b01=0.00001,side='lower')
    
  })

decision<-t(sapply(sim4, "[[", 1))
sum(decision[,2]=='Superiorityfinal')/1000 #0.882
#The simulated power from 1000 simulations is 0.882.
```

--------------------------------

### Convert Normal-Inverse-Chi-Squared parameters to Normal-Inverse-Gamma

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/convert_chisq_to_gamma.Rd

Converts a Normal-Inverse-Chi-Squared parameter list (mu, kappa, nu, sigsq) to a Normal-Inverse-Gamma parameter list (m, V, a, b). The example uses mu=0.091, kappa=2, nu=1, and sigsq=4e-05.

```r
convert_chisq_to_gamma(list(mu=0.091,kappa=2,nu=1,sigsq=4e-05))
```

--------------------------------

### Convert Normal-Inverse-Gamma parameters to Normal-Inverse-Chi-Squared

Source: https://github.com/yayayaoyaoyao/rartrials/blob/main/man/convert_gamma_to_chisq.Rd

Converts a Normal-Inverse-Gamma parameter list (m, V, a, b) to a Normal-Inverse-Chi-Squared parameter list (mu, kappa, nu, sigsq). The example uses V=1/2, a=0.5, m=9.1/100, and b=0.00002.

```r
convert_gamma_to_chisq(list(V=1/2,a=0.5,m=9.1/100,b=0.00002))
```
Tokens:
Raw