Add finite sample correction function

R/robincar-covhr.R

@@ -17,6 +17,7 @@
 #' @param p_trt Treatment allocation ratio for the reference arm.
 #' @param adj_method Adjustment method (one of "CL", "CSL")
 #' @param interval Interval for uniroot function
+#' @param finite_sample Perform a finite-sample size correction
 #'
 #' @export
 #'
@@ -35,7 +36,8 @@ robincar_covhr <- function(df,
                            p_trt=0.5, ref_arm=NULL,
                            car_scheme="simple",
                            adj_method="CL",
-                           interval=c(-10, 10)){
+                           interval=c(-10, 10),
+                           finite_sample = FALSE){
 
   .check.car_scheme(car_scheme, car_strata_cols)
 
@@ -69,6 +71,14 @@ robincar_covhr <- function(df,
   result <- adjust(model, data)
   result$original_df <- df
 
+  ## Finite sample correction
+  if (finite_sample){
+
+    result <- apply_correction(result)
+
+  }
+
+
   return(result)
 
 }

R/robincar-tte.R

@@ -22,6 +22,7 @@
 #' @param sparse_remove Remove sparse car_strata from calculation
 #' @param return_influence Return the influence function for the covariate-adjusted log-rank test for each individual
 #' @param id_col Name of the column in df indicating individual ids
+#' @param finite_sample Perform a finite-sample size correction
 #'
 #' @export
 #'
@@ -31,7 +32,7 @@ robincar_tte <- function(df,
                          adj_method,
                          car_strata_cols=NULL, covariate_cols=NULL,
                          p_trt=0.5, ref_arm=NULL, sparse_remove=TRUE,
-                         car_scheme="simple", return_influence=FALSE, id_col=NULL){
+                         car_scheme="simple", return_influence=FALSE, id_col=NULL, finite_sample = FALSE){
 
   .check.car_scheme(car_scheme, car_strata_cols)
   data <- .make.data(
@@ -61,6 +62,13 @@ robincar_tte <- function(df,
   result <- adjust(model, data, return_influence)
   result$original_df <- df
 
+  ## Finite sample correction
+  if (finite_sample){
+
+    result <- apply_correction(result)
+
+  }
+
   return(result)
 
 }

R/utils-results.R

@@ -16,3 +16,51 @@ format_results <- function(labels, estimates, varcov, label_name="treat", ...){
   )
   return(result)
 }
+
+
+
+## function to apply finite sample correction
+apply_correction <- function(result){
+
+  if (!(result$settings$method %in% c("CL", "CSL"))) return(result)  
+
+  n <- result$data$n
+  n_0 <- n * result$data$pie[1,1]
+  n_1 <- n * result$data$pie[2,1]
+
+
+  ### work out how many terms come from the covariate columns...
+  form <- paste0(result$data$response_col, " ~ ", paste(c("1", result$data$covariate_cols), collapse = "+"))
+  fit_lin <- lm(as.formula(form), data = result$original_df)
+  n_terms_cov <- length(coef(fit_lin)) - 1
+
+  ## add on the covariate strata terms (if needed)
+  n_terms_strata <- length(result$data$joint_strata_levels) - 1
+
+  n_terms <- ifelse(result$settings$method == "CL" && result$settings$car_scheme != "simple",
+                    n_terms_cov + n_terms_strata,
+                    n_terms_cov)
+
+  ## formula from Krotka & Magirr
+  correction <- (1 + n_terms * (((n_0 / n)^2) / (n_1 - n_terms - 2) + ((n_1 / n)^2) / (n_0 - n_terms - 2))) *
+    (n / (n - ((n_0 / n_1) + (n_1 / n_0) - 1) * n_terms - ((n_0 / n_1) + (n_1 / n_0))))
+
+  correction <- unname(correction)
+
+  ## what to do depends on estimation vs testing
+  if (attr(result$settings, "class")[1] == "CovHR"){
+
+    result$result$se_theta_CL <- result$result$se_theta_CL * sqrt(correction)
+
+
+  }
+  else {
+
+    result$result$se <- result$result$se * sqrt(correction)
+    result$result$statistic <- result$result$U / result$result$se
+
+  }
+
+  return(result)
+
+}

tests/testthat/test-finite-sample.R

@@ -0,0 +1,108 @@
+library(dplyr)
+
+test_that("test finite sample correction", {
+
+  set.seed(10)
+  n <- 100
+
+  # Simulate data
+  data.simu <- data_gen2(
+    n=n,
+    theta=0,
+    randomization="permuted_block",
+    p_trt=0.5,
+    case="case1"
+  )
+
+  data.simu$strata_2 <- factor(rbinom(100, 1, 0.5))
+  data.simu$strata_3 <- factor(rbinom(100, 1, 0.5))
+
+  cl_1_F <- robincar_logrank(
+    df=data.simu,
+    treat_col="I1",
+    response_col="t",
+    event_col="delta",
+    car_strata_cols=c("strata_2", "strata_3"),
+    covariate_cols=c("model_w3"),
+    car_scheme="permuted-block",
+    adj_method="CL",
+    ref_arm=0,
+    finite_sample = FALSE
+  )
+
+  cl_1_T <- robincar_logrank(
+    df=data.simu,
+    treat_col="I1",
+    response_col="t",
+    event_col="delta",
+    car_strata_cols=c("strata_2", "strata_3"),
+    covariate_cols=c("model_w3"),
+    car_scheme="permuted-block",
+    adj_method="CL",
+    ref_arm=0,
+    finite_sample = TRUE
+  )
+
+  expect_equal(round(cl_1_F$result$se, 6), 0.045337)
+  expect_equal(round(cl_1_T$result$se, 6), 0.047828)
+
+  cl_2_F <- robincar_logrank(
+    df=data.simu,
+    treat_col="I1",
+    response_col="t",
+    event_col="delta",
+    car_strata_cols=NULL,
+    covariate_cols=c("model_w3", "strata_2", "strata_3"),
+    car_scheme="simple",
+    adj_method="CL",
+    ref_arm=0,
+    finite_sample = FALSE
+  )
+
+  cl_2_T <- robincar_logrank(
+    df=data.simu,
+    treat_col="I1",
+    response_col="t",
+    event_col="delta",
+    car_strata_cols=NULL,
+    covariate_cols=c("model_w3", "strata_2", "strata_3"),
+    car_scheme="simple",
+    adj_method="CL",
+    ref_arm=0,
+    finite_sample = TRUE
+  )
+
+  expect_equal(round(cl_2_F$result$se, 6), 0.045343)
+  expect_equal(round(cl_2_T$result$se, 6), 0.047302)
+
+
+  csl_1_F <- robincar_logrank(
+    df=data.simu,
+    treat_col="I1",
+    response_col="t",
+    event_col="delta",
+    car_strata_cols=c("strata_2", "strata_3"),
+    covariate_cols=c("model_w3"),
+    car_scheme="permuted-block",
+    adj_method="CSL",
+    ref_arm=0,
+    finite_sample = FALSE
+  )
+
+  csl_1_T <- robincar_logrank(
+    df=data.simu,
+    treat_col="I1",
+    response_col="t",
+    event_col="delta",
+    car_strata_cols=c("strata_2", "strata_3"),
+    covariate_cols=c("model_w3"),
+    car_scheme="permuted-block",
+    adj_method="CSL",
+    ref_arm=0,
+    finite_sample = TRUE
+  )
+
+  expect_equal(round(csl_1_F$result$se, 6), 0.043468)
+  expect_equal(round(csl_1_T$result$se, 6), 0.044374)
+
+})