PITS: Power of an ITS — CDSS/CFR worked example

The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

PITS package

2026-06-23

Motivating example

A hospital is planning to deploy a clinical decision support system (CDSS) intended to reduce its case fatality rate (CFR). The research team wants to answer, prospectively:

“How many months of post-implementation follow-up do we need to reliably detect a meaningful reduction in CFR?”

Step 1 — Load pre-intervention data

data("example_cfr_data")
head(example_cfr_data)
#>   time outcome
#> 1    1   12.79
#> 2    2   14.58
#> 3    3   16.44
#> 4    4   18.41
#> 5    5   14.47
#> 6    6   13.21

plot(example_cfr_data$time, example_cfr_data$outcome,
     type = "o", pch = 16, col = "steelblue", lwd = 1.2,
     xlab = "Month", ylab = "CFR (%)",
     main = "Pre-intervention case fatality rate",
     las = 1, bty = "l")
abline(lm(outcome ~ time, data = example_cfr_data),
       col = "firebrick", lwd = 2, lty = 2)
legend("topright", legend = c("Observed", "Linear trend"),
       col = c("steelblue", "firebrick"), lty = 1,
       pch = c(16, NA), bty = "n", cex = 0.85)
grid(NULL, NULL, lwd = 0.6, col = rgb(0, 0, 0, 0.1))

Monthly CFR over 24 pre-intervention months.

Step 2 — Estimate nuisance parameters

params <- estimate_its_params(example_cfr_data)
#> 
#> ============================================================
#>   PITS - ITS Parameter Estimation
#> ============================================================ 
#> 
#>   Estimation method:      GLS-ML
#>   Observations (n_pre):   24
#>   Baseline:               14.5853
#>   Sigma (residual SD):    1.4876  (10.2% of baseline)
#>   Rho (AR1):              0.3832
#>   Pre-trend (per unit):   0.0039
#> 
#>  ------------------------------------------------------------ 
#>   Copy these into calculate_power():
#> 
#>     n_pre        = 24
#>     baseline     = 14.5853
#>     sigma        = 1.4876
#>     rho          = 0.3832
#>     level_change = ???  # YOUR clinical hypothesis
#>     slope_change = 0    # set if testing slope
#> 
#> ============================================================

{r diag-plots, fig.width = 8, fig.height = 6, fig.cap = "Diagnostic plots for pre-intervention residuals."} diagnose_params(example_cfr_data)

The residual ACF confirms moderate positive autocorrelation (rho ≈ 0.37), typical for monthly hospital data. Ignoring this would overstate power.

Step 3 — Specify the clinical hypothesis

level_change is not estimated from data — it is a clinical judgement: the minimum reduction the team would consider practice-changing.

level_change <- -3.0   # 3 pp reduction in CFR

Step 4 — Calculate power for a single design

Note on simulation size. To keep this vignette fast to build, the examples below use a small number of Monte Carlo replications (n_sim = 100–200). For a reportable analysis, use n_sim = 1000 or more for stable power estimates.

result_24 <- calculate_power(
  n_pre = params$n_pre, n_post = 24,
  baseline = params$baseline, level_change = level_change,
  sigma = params$sigma, rho = params$rho,
  test = "level", n_sim = 200, seed = 123
)
print(result_24)
#> 
#> ============================================================
#>   PITS - ITS Power Estimate
#> ============================================================ 
#> 
#>   Design:          Single site
#>   Test:            level
#>   n_pre:           24
#>   n_post:          24
#>   Baseline:        14.585
#>   Level change:    -3.000
#>   Sigma:           1.488
#>   Rho:             0.383
#>   Alpha:           0.050
#>   N simulations:   200
#>   N converged:     200 (100.0%)
#> 
#> ------------------------------------------------------------
#>   POWER = 77.5%  -  Borderline (60-79%)
#> ------------------------------------------------------------

24 months is borderline. The sweep below identifies the minimum adequate duration.

Step 5 — Optimise with a design sweep

sweep <- power_sweep(
  sweep_post   = c(12, 18, 24, 30, 36, 42, 48),
  n_pre        = params$n_pre,
  baseline     = params$baseline,
  level_change = level_change,
  sigma        = params$sigma,
  rho          = params$rho,
  test         = "level",
  n_sim        = 200, seed = 123
)
#> 
#> =======================================================
#>   PITS - Design optimisation sweep
#> ======================================================= 
#>   n_post       Power  Interpretation
#> -------------------------------------------------------
#>   12           70.5%  Borderline (60-79%)
#>   18           73.5%  Borderline (60-79%)
#>   24           77.5%  Borderline (60-79%)
#>   30           80.5%  Adequate (>= 80%)  <-- adequate
#>   36           87.5%  Adequate (>= 80%)  <-- adequate
#>   42           87.0%  Adequate (>= 80%)  <-- adequate
#>   48           91.5%  Adequate (>= 80%)  <-- adequate
#> =======================================================

plot_power_curve(
  sweep,
  main = "Power vs follow-up duration — CDSS/CFR example",
  xlab = "Post-implementation months (n_post)",
  ylab = "Estimated power (%)"
)

Power curve. Dashed red = 80% target; green = minimum adequate n_post.

Step 6 — Plot a simulated ITS realisation

plot_its_example(
  n_pre = params$n_pre, n_post = 30,
  baseline = params$baseline, level_change = level_change,
  sigma = params$sigma, rho = params$rho,
  seed = 7,
  xlab = "Month", ylab = "CFR (%)",
  main = "Simulated ITS — CDSS effect on case fatality rate"
)

Simulated ITS under the alternative hypothesis with fitted segmented regression.

Step 7 — Sensitivity analysis

How sensitive is the required follow-up to uncertainty in sigma and rho?

grid <- build_param_grid(
  n_post       = c(18, 24, 30, 36),
  level_change = level_change,
  sigma        = c(1.0, round(params$sigma, 1), 2.0),
  rho          = c(0.2, 0.4, 0.6),
  n_pre        = params$n_pre,
  baseline     = params$baseline
)
grid_results <- run_power_grid(grid, n_sim = 100, seed = 42, verbose = FALSE)

{r heatmap-sigma, fig.height = 5, fig.cap = "Power by n_post and sigma (rho = 0.4). Dashed line = first adequate column."} plot_power_heatmap( grid_results[grid_results$rho == 0.4, ], x = "n_post", y = "sigma", main = "Power heatmap: n_post × sigma (rho = 0.4)", xlab = "Post-implementation months", ylab = "Residual SD (sigma)" )

{r heatmap-rho, fig.height = 5, fig.cap = "Power by n_post and rho (sigma fixed at estimated value)."} sigma_val <- round(params$sigma, 1) plot_power_heatmap( grid_results[round(grid_results$sigma, 1) == sigma_val, ], x = "n_post", y = "rho", main = sprintf("Power heatmap: n_post × rho (sigma = %.1f)", sigma_val), xlab = "Post-implementation months", ylab = "Autocorrelation (rho)" )

Step 8 — Multi-site analysis

If the CDSS is deployed simultaneously across three hospitals, pooling evidence increases power substantially.

sites <- list(
  list(name = "Hospital A", n_pre = 24, n_post = 24,
       baseline = params$baseline, level_change = level_change,
       slope_change = 0, sigma = params$sigma, rho = params$rho),
  list(name = "Hospital B", n_pre = 24, n_post = 24,
       baseline = 17.5, level_change = level_change,
       slope_change = 0, sigma = 2.0, rho = 0.35),
  list(name = "Hospital C", n_pre = 24, n_post = 24,
       baseline = 13.0, level_change = level_change,
       slope_change = 0, sigma = 1.8, rho = 0.42)
)
result_multi <- calculate_power_multi(
  sites, test = "level", n_sim = 200, seed = 123
)
print(result_multi)
#> 
#> ============================================================
#>   PITS - ITS Power Estimate
#> ============================================================ 
#> 
#>   Design:          Multi-site (3 sites)
#>   Test:            level
#>   Alpha:           0.050
#>   N simulations:   200
#>   N converged:     200 (100.0%)
#> 
#> ------------------------------------------------------------
#>   POWER = 96.0%  -  Adequate (>= 80%)
#> ------------------------------------------------------------

Summary table

Power by follow-up duration. baseline = 14.6%, sigma = 1.49, rho = 0.38, level_change = -3 pp.
n_post (months)	Power (%)	Interpretation
12	70.5	Borderline (60-79%)
18	73.5	Borderline (60-79%)
24	77.5	Borderline (60-79%)
30	80.5	Adequate (>= 80%)
36	87.5	Adequate (>= 80%)
42	87.0	Adequate (>= 80%)
48	91.5	Adequate (>= 80%)

Conclusion: The study requires approximately 30 months of post-implementation follow-up to achieve ≥ 80% power.

One-call shortcut

result_quick <- estimate_and_calculate(
  data = example_cfr_data,
  level_change = -3, n_post = 30,
  n_sim = 200, seed = 123
)
#> 
#> ============================================================
#>   PITS - ITS Parameter Estimation
#> ============================================================ 
#> 
#>   Estimation method:      GLS-ML
#>   Observations (n_pre):   24
#>   Baseline:               14.5853
#>   Sigma (residual SD):    1.4876  (10.2% of baseline)
#>   Rho (AR1):              0.3832
#>   Pre-trend (per unit):   0.0039
#> 
#>  ------------------------------------------------------------ 
#>   Copy these into calculate_power():
#> 
#>     n_pre        = 24
#>     baseline     = 14.5853
#>     sigma        = 1.4876
#>     rho          = 0.3832
#>     level_change = ???  # YOUR clinical hypothesis
#>     slope_change = 0    # set if testing slope
#> 
#> ============================================================
cat(sprintf("Power: %.1f%%  (%s)\n",
            result_quick$power_pct, result_quick$interpretation))
#> Power: 86.0%  (Adequate (>= 80%))

References

Lopez Bernal J, et al. (2017). Int J Epidemiol 46:348–355.
Zhang F, et al. (2011). J Clin Epidemiol 64:1252–1261.
Wagner AK, et al. (2002). J Clin Pharm Ther 27:299–309.

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.