chestR workflow: global Cox to local biomarker maps

Richard Jackson

2026-09-03

Overview

chestR re-fits a Cox model at each point on a biomarker grid using kernel weights. Local coefficients show how treatment (and other) effects vary across biomarker space. The usual workflow is:

  1. Fit a global Cox model.
  2. Call chestr() on biomarker values.
  3. Visualise a coefficient with plot() / plot.chestr().
  4. Optionally run chestr_test() for a permutation global test. This vignette uses a small simulated survival dataset so it builds without external files.

Simulate data and fit a global Cox model

library(survival)
library(chestR)

set.seed(20250806)
n <- 120
dat <- data.frame(
  trt = rbinom(n, 1, 0.5),
  cov = rnorm(n),
  biom1 = rnorm(n),
  biom2 = rnorm(n)
)

# Mild treatment-effect modification near the origin of biomarker space
lp <- 0.2 * dat$cov +
  log(0.7) * dat$trt *
  exp(-0.5 * (dat$biom1^2 + dat$biom2^2))
dat$time <- rexp(n, rate = exp(lp))
dat$status <- as.integer(dat$time < 4)
dat$time[dat$status == 0L] <- 4

base <- coxph(Surv(time, status) ~ trt + cov, data = dat)
summary(base)$coefficients
#>           coef exp(coef)  se(coef)         z   Pr(>|z|)
#> trt -0.2064674 0.8134528 0.1879316 -1.098630 0.27192939
#> cov  0.2656864 1.3043260 0.1096953  2.422039 0.01543368

Local estimates with chestr()

biom <- dat[, c("biom1", "biom2")]
cr <- chestr(
  base,
  biom,
  grid.size = 8,
  method = "legacy",
  kern.adj = 2,
  min_events_per_df = 5
)
#> Warning: Skipped 42 of 64 grid point(s) with fewer than 5 Kish effective events
#> per model df (2 df; require >= 10 Kish ESS events). Coefficients set to NA.

cr
#> chestr object
#>   biomarkers : biom1, biom2 
#>   method     : legacy 
#>   data stored: yes 
#>   grid points:64 (reliable 22)
#> 
#> Estimates (head):
#>        biom1     biom2       trt       cov   trt.se  cov.se      eff.n
#> 1 -2.8548316 -1.787448        NA        NA       NA      NA 0.01556175
#> 2 -2.0889085 -1.787448        NA        NA       NA      NA 0.26021997
#> 3 -1.3229854 -1.787448        NA        NA       NA      NA 1.94667894
#> 4 -0.5570622 -1.787448        NA        NA       NA      NA 5.45277603
#> 5  0.2088609 -1.787448 0.7109821 0.3431944 1.033327 0.57304 5.21807707
#> 6  0.9747840 -1.787448        NA        NA       NA      NA 2.11586918
#>        eff.e  eff.event ess_events events_per_df reliable
#> 1 0.01556175 0.01556175   3.565665      1.782833    FALSE
#> 2 0.26021981 0.26021981   2.751248      1.375624    FALSE
#> 3 1.94666793 1.94666793   6.639199      3.319600    FALSE
#> 4 5.45263171 5.45263171   9.941330      4.970665    FALSE
#> 5 5.21784860 5.21784860  12.162757      6.081378     TRUE
#> 6 2.11584143 2.11584143   5.967228      2.983614    FALSE
head(cr$estimates[, c("biom1", "biom2", "trt", "ess_events",
                      "events_per_df", "reliable")])
#>        biom1     biom2       trt ess_events events_per_df reliable
#> 1 -2.8548316 -1.787448        NA   3.565665      1.782833    FALSE
#> 2 -2.0889085 -1.787448        NA   2.751248      1.375624    FALSE
#> 3 -1.3229854 -1.787448        NA   6.639199      3.319600    FALSE
#> 4 -0.5570622 -1.787448        NA   9.941330      4.970665    FALSE
#> 5  0.2088609 -1.787448 0.7109821  12.162757      6.081378     TRUE
#> 6  0.9747840 -1.787448        NA   5.967228      2.983614    FALSE

The object has class "chestr" and stores:

Plot local treatment effects

plot(cr, trt.param = "trt", col.scale = "obs", reliable_only = TRUE)

plot.chestr() takes everything it needs from cr (base and biom are not passed again). By default only reliable grid points are drawn.

Permutation global test

Pointwise local estimates from chestr() are descriptive. Formal evidence that the local treatment surface is not flat can be assessed with chestr_test(), which permutes treatment labels and recomputes T_L2 and T_MAX.

cr <- chestr(
  base, biom,
  grid.size = 8, method = "legacy", kern.adj = 2,
  min_events_per_df = 5,
  treat_term = "trt"
)
# data, treat_term, and fit settings come from cr:
tst <- chestr_test(cr, B = 99, seed = 1)
tst

Session info

sessionInfo()
#> R version 4.4.2 (2024-10-31)
#> Platform: aarch64-apple-darwin20
#> Running under: macOS 26.6.2
#> 
#> Matrix products: default
#> BLAS:   /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib 
#> LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
#> 
#> locale:
#> [1] C/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
#> 
#> time zone: Europe/London
#> tzcode source: internal
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] chestR_0.1.0   survival_3.8-3
#> 
#> loaded via a namespace (and not attached):
#>  [1] vctrs_0.7.2        cli_3.6.5          knitr_1.50         rlang_1.1.7       
#>  [5] xfun_0.54          generics_0.1.4     S7_0.2.1           jsonlite_2.0.0    
#>  [9] labeling_0.4.3     glue_1.8.0         htmltools_0.5.8.1  sass_0.4.10       
#> [13] scales_1.4.0       rmarkdown_2.30     grid_4.4.2         tibble_3.3.0      
#> [17] evaluate_1.0.5     jquerylib_0.1.4    fastmap_1.2.0      yaml_2.3.10       
#> [21] lifecycle_1.0.5    compiler_4.4.2     dplyr_1.1.4        RColorBrewer_1.1-3
#> [25] pkgconfig_2.0.3    rstudioapi_0.17.1  farver_2.1.2       lattice_0.22-7    
#> [29] digest_0.6.38      R6_2.6.1           tidyselect_1.2.1   pillar_1.11.1     
#> [33] splines_4.4.2      magrittr_2.0.4     bslib_0.9.0        Matrix_1.7-4      
#> [37] withr_3.0.2        tools_4.4.2        gtable_0.3.6       ggplot2_4.0.2     
#> [41] cachem_1.1.0