ready set go

2026-06-03 11:28:14 -08:00
commit 097f88f85d
22 changed files with 8793 additions and 0 deletions
@@ -0,0 +1,18 @@
+# Simulate an ICU cohort for learning target trial emulation.
+# The simulation logic lives in R/simulate_icu_cohort.R.
+
+suppressPackageStartupMessages(library(readr))
+suppressPackageStartupMessages(library(skimr))
+
+source("R/simulate_icu_cohort.R")
+
+icu_data <- simulate_icu_cohort(n_patients = 1000, seed = 20260531)
+
+write_csv(
+  icu_data,
+  file = "data/icu_septic_shock_base_r.csv"
+)
+
+cli::cli_alert_success("Wrote simulated cohort to {.file data/icu_septic_shock_base_r.csv}.")
+
+skim(icu_data)
@@ -0,0 +1,91 @@
+suppressPackageStartupMessages(library(dplyr))
+suppressPackageStartupMessages(library(gt))
+suppressPackageStartupMessages(library(gtsummary))
+suppressPackageStartupMessages(library(readr))
+suppressPackageStartupMessages(library(tidyr))
+
+cli::cli_h1("Naive analysis of the simulated ICU cohort")
+
+icu_data <- read_csv("data/icu_septic_shock_base_r.csv", show_col_types = FALSE)
+
+if (!dir.exists("outputs")) {
+  dir.create("outputs")
+}
+
+eligible_data <- icu_data |>
+  filter(eligible == 1)
+
+analysis_summary <- eligible_data |>
+  summarize(
+    n_eligible = n(),
+    n_early = sum(early_vasopressor == 1),
+    n_not_early = sum(early_vasopressor == 0),
+    risk_early = mean(death_28d[early_vasopressor == 1]),
+    risk_not_early = mean(death_28d[early_vasopressor == 0]),
+    risk_difference = risk_early - risk_not_early,
+    risk_ratio = risk_early / risk_not_early
+  )
+
+analysis_results <- analysis_summary |>
+  pivot_longer(
+    cols = everything(),
+    names_to = "measure",
+    values_to = "value"
+  )
+
+analysis_table <- analysis_results |>
+  gt() |>
+  tab_header(title = "Naive 28-Day Mortality Comparison") |>
+  cols_label(
+    measure = "Measure",
+    value = "Value"
+  ) |>
+  fmt_number(columns = value, decimals = 3)
+
+baseline_table <- eligible_data |>
+  mutate(
+    early_vasopressor = factor(
+      early_vasopressor,
+      levels = c(0, 1),
+      labels = c("No early vasopressor", "Early vasopressor")
+    )
+  ) |>
+  select(early_vasopressor, age, sex, sofa_score, lactate, map, death_28d) |>
+  tbl_summary(
+    by = early_vasopressor,
+    statistic = list(
+      all_continuous() ~ "{mean} ({sd})",
+      all_categorical() ~ "{n} ({p}%)"
+    ),
+    missing = "no"
+  ) |>
+  add_overall()
+
+baseline_console <- eligible_data |>
+  group_by(early_vasopressor) |>
+  summarize(
+    mean_age = mean(age),
+    mean_sofa_score = mean(sofa_score),
+    mean_lactate = mean(lactate),
+    mean_map = mean(map),
+    .groups = "drop"
+  )
+
+cli::cli_h2("Core causal contrast")
+print(analysis_results)
+
+gtsave(
+  data = analysis_table,
+  filename = "outputs/naive_mortality_comparison.html"
+)
+
+cli::cli_alert_success("Wrote {.file outputs/naive_mortality_comparison.html}.")
+
+cli::cli_h2("Baseline comparison by observed treatment")
+print(baseline_console)
+
+baseline_table |>
+  as_gt() |>
+  gtsave(filename = "outputs/baseline_by_treatment.html")
+
+cli::cli_alert_success("Wrote {.file outputs/baseline_by_treatment.html}.")