Single_analysis_workflow

library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(rbmi)
library(rbmiUtils)
library(beeca)

This is an example of the {rbmi} only workflow. The analysis is as follows: Rubin’s Rules after standardization: Following multiple imputation, we would fit 100 logistic regression models and for each, obtain the marginal estimates using beeca. We would then pool these 100 beeca outputs using Rubin’s Rules.

This example exists to compare with the workflow which saves an ADMI data set, then uses {rbmiUtils} to complete the intended analysis.

set.seed(1974)
data("ADEFF")

## Global constants
N_IMPUTATIONS <- 100 

# Prepare the data
ADEFF$TRT <- factor(ADEFF$TRT01P, levels = c("Placebo", "Drug A"))
ADEFF$USUBJID <- factor(ADEFF$USUBJID)
ADEFF$AVISIT <- factor(ADEFF$AVISIT)

Define function to compare last visit using a responder analysis

## function to compare last visit using a responder analysis
gcomp_responder <- function(data, ...) {

   data <- data |>
    dplyr::mutate(
        CRIT1FLN = ifelse(CHG > 3, 1, 0), 
        CRIT1FL = ifelse(CRIT1FLN == 1, "Y", "N"),
        CRIT = "CHG > 3"
    )
  
    ## fit the GLM model
    model <- glm(CRIT1FLN ~ TRT + BASE + STRATA + REGION, data = data, family = binomial)

    ## Calculate the marginal treatment effect estimate and associated variance
    ## for a difference contrast using the Ge et al. method
    ## with robust HC0 sandwich variance
    marginal_fit <- get_marginal_effect(model,
                                        trt = "TRT",
                                        method = "Ge",
                                        type = "HC0",
                                        contrast = "diff",
                                        reference = "Placebo")

    ## Extract results
    marginal_results <- marginal_fit$marginal_results
    diff_est <- marginal_results[marginal_results$STAT == "diff", "STATVAL"][[1]]
    diff_se <- marginal_results[marginal_results$STAT == "diff_se", "STATVAL"][[1]]

    # Return a named list with the estimate and standard error
    res <- list(
        trt = list(
            est = diff_est,
            se = diff_se,
            df = NA
        )
    )
    return(res)
}

STEP 1: Impute missing data

# Define key variables for rbmi
vars <- rbmi::set_vars(
    subjid = "USUBJID",
    visit = "AVISIT",
    group = "TRT",
    outcome = "CHG",
    covariates = c("BASE", "STRATA", "REGION")
)

# Define which imputation method to use (here: Bayesian multiple imputation)
method <- rbmi::method_bayes(
    n_samples = N_IMPUTATIONS,
    burn_in = 200,
    burn_between = 5
)

# Create a dataset with the necessary columns
dat <- ADEFF |>
    select(USUBJID, STRATA, REGION, REGIONC, TRT, BASE, CHG, AVISIT)

# Fit the imputation model
draws_obj <- rbmi::draws(
    data = dat,
    vars = vars,
    method = method,
    quiet = TRUE
)

# Impute the datasets
impute_obj <- rbmi::impute(
    draws_obj,
    references = c("Placebo" = "Placebo", "Drug A" = "Placebo")
)

STEP 2: Analyze each imputed dataset

ANCOVA analysis

ana_obj_prop <- analyse(impute_obj,
                         vars = vars,
                         fun = ancova,
                         delta = NULL)

STEP 3: Pool analyzes and apply Rubin’s Rules

pool_obj_prop <- pool(ana_obj_prop)

STEP 4: Report analysis

pool_obj_prop
#> 
#> Pool Object
#> -----------
#> Number of Results Combined: 100
#> Method: rubin
#> Confidence Level: 0.95
#> Alternative: two.sided
#> 
#> Results:
#> 
#>   ========================================================
#>       parameter      est     se     lci     uci     pval  
#>   --------------------------------------------------------
#>      trt_Week 24    -2.177  0.182  -2.535  -1.819  <0.001 
#>    lsm_ref_Week 24  0.077   0.131  -0.18   0.334   0.557  
#>    lsm_alt_Week 24   -2.1   0.125  -2.347  -1.854  <0.001 
#>      trt_Week 48    -3.806  0.256  -4.308  -3.303  <0.001 
#>    lsm_ref_Week 48  0.044   0.185  -0.32   0.407   0.814  
#>    lsm_alt_Week 48  -3.762  0.175  -4.106  -3.417  <0.001 
#>   --------------------------------------------------------

Responder analysis

ana_obj_prop <- analyse(impute_obj,
                         vars = vars,
                         method = method,
                         fun = gcomp_responder,
                         delta = NULL)

STEP 3: Pool analyzes and apply Rubin’s Rules

pool_obj_prop <- pool(ana_obj_prop)

STEP 4: Report analysis

pool_obj_prop
#> 
#> Pool Object
#> -----------
#> Number of Results Combined: 100
#> Method: rubin
#> Confidence Level: 0.95
#> Alternative: two.sided
#> 
#> Results:
#> 
#>   ==================================================
#>    parameter   est     se     lci     uci     pval  
#>   --------------------------------------------------
#>       trt     -0.063  0.012  -0.087  -0.039  <0.001 
#>   --------------------------------------------------