## ----function-signature, eval=FALSE-------------------------------------------
# MRStdCRT_fit(
#   formula, data, cluster, trt, trtprob = rep(0.5, nrow(data)),
#   method, family = gaussian(link = "identity"), corstr, scale,
#   alpha = 0.05
# )

## ----transformed-covariate-formula, eval=FALSE--------------------------------
# data$x_sq <- data$x^2
# formula <- y ~ x + x_sq

## ----treatment-probabilities--------------------------------------------------
library(MRStdCRT)
library(dplyr)
data(ppact)

prob <- rep(0.5, nrow(ppact))

## ----estimated-treatment-probability------------------------------------------
cluster_assignments <- ppact %>%
  distinct(CLUST, INTERVENTION)

estimated_prob <- rep(
  mean(cluster_assignments$INTERVENTION == 1),
  nrow(ppact)
)

## ----hypothetical-block-probabilities-----------------------------------------
clusters <- sort(unique(ppact$CLUST))
n_clusters <- length(clusters) 
block_info <- data.frame(
  CLUST = clusters
) %>%
  mutate(
    # Create a 'block' identifier for each cluster
    block = case_when(
      row_number() <= 25 ~ 1,
      row_number() <= 66 ~ 2,
      TRUE ~ 3
    )
  ) %>%
  mutate(
    # Assign the desired hypothetical treatment probability to each block
    hypothetical_prob = case_when(
      block == 1 ~ 0.3,
      block == 2 ~ 0.5,
      block == 3 ~ 0.6
    )
  )
prob_lookup <- setNames(block_info$hypothetical_prob, block_info$CLUST)
probs_vector <- prob_lookup[as.character(ppact$CLUST)]

## ----fit-ppact----------------------------------------------------------------
example <- MRStdCRT_fit(
  formula = PEGS ~ AGE + FEMALE + comorbid + Dep_OR_Anx + pain_count + PEGS_bl +
    BL_benzo_flag + BL_avg_daily + satisfied_primary + n,
  data = ppact,
  cluster = "CLUST",
  trt = "INTERVENTION",
  trtprob = prob,
  method = "GEE",
  corstr = "independence",
  scale = "RR"
)
## To view the summary, use the following command
summary(example)
## One may also extract specific statistical feature, such as table for point and interval estimates, and standard errors
example$estimate

