Skip to content

add semi-synthetic causal simulator with observed treatment - #173

Open
kirilklein wants to merge 20 commits into
mainfrom
improved-simulation
Open

add semi-synthetic causal simulator with observed treatment#173
kirilklein wants to merge 20 commits into
mainfrom
improved-simulation

Conversation

@kirilklein

@kirilklein kirilklein commented Apr 3, 2026

Copy link
Copy Markdown
Owner

Introduces a new simulation DGP where treatment assignment is kept from real data and only the outcome is simulated via hand-crafted oracle features extracted from patient histories. Features are split into baseline risk (r_B) and longitudinal (r_L) groups, with configurable coefficients and treatment effects defined in YAML.

New files:

  • config_semisynthetic.py: dataclass configs + factory function
  • oracle_features.py: 10 feature extractors (7 baseline, 3 longitudinal)
  • semisynthetic_simulator.py: simulator class with same output format
  • simulate_semisynthetic.py: pipeline entry point
  • calibrate_semisynthetic.py: standalone calibration diagnostics + plots
  • simulate_semisynthetic.yaml: default config
  • tests for feature extraction and simulator (17 tests)

Summary by CodeRabbit

  • New Features

    • Adds a semi-synthetic causal simulation pipeline with configurable outcome models, deterministic seeds, calibration and simulation CLIs, Azure component entrypoints, and outputs (per-outcome CSVs, ITEs, counterfactuals, event records, and diagnostic figures).
    • Oracle-style patient feature extraction (temporal windows, motifs, standardization).
    • Multi-run config generator for deterministic replicate studies.
  • Documentation

    • New README and docs for local/Azure runs, calibration, features, and multi-run workflows.
  • Tests

    • Unit and end-to-end tests for feature extraction and the simulator.

Introduces a new simulation DGP where treatment assignment is kept from
real data and only the outcome is simulated via hand-crafted oracle
features extracted from patient histories. Features are split into
baseline risk (r_B) and longitudinal (r_L) groups, with configurable
coefficients and treatment effects defined in YAML.

New files:
- config_semisynthetic.py: dataclass configs + factory function
- oracle_features.py: 10 feature extractors (7 baseline, 3 longitudinal)
- semisynthetic_simulator.py: simulator class with same output format
- simulate_semisynthetic.py: pipeline entry point
- calibrate_semisynthetic.py: standalone calibration diagnostics + plots
- simulate_semisynthetic.yaml: default config
- tests for feature extraction and simulator (17 tests)

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
@coderabbitai

coderabbitai Bot commented Apr 3, 2026

Copy link
Copy Markdown
Contributor

Note

Reviews paused

It looks like this branch is under active development. To avoid overwhelming you with review comments due to an influx of new commits, CodeRabbit has automatically paused this review. You can configure this behavior by changing the reviews.auto_review.auto_pause_after_reviewed_commits setting.

Use the following commands to manage reviews:

  • @coderabbitai resume to resume automatic reviews.
  • @coderabbitai review to trigger a single review.

Use the checkboxes below for quick actions:

  • ▶️ Resume reviews
  • 🔍 Trigger review
📝 Walkthrough

Walkthrough

Adds a semi-synthetic causal simulation subsystem: new YAML config, typed simulation configs, oracle feature extractor, a two-pass SemiSyntheticCausalSimulator producing factual/counterfactual outcomes and ITEs, CLI entrypoints and Azure components for simulate/calibrate, unit tests, docs, and multi-run config generator.

Changes

Cohort / File(s) Summary
Config & entrypoint YAML
corebehrt/configs/causal/simulate_semisynthetic.yaml
New simulation YAML: logging, input/output paths, seed/debug/min_num_codes/exposure_code, feature code-prefixes and time windows, and two outcome definitions (OUTCOME, OUTCOME_NULL).
Typed config model
corebehrt/modules/simulation/config_semisynthetic.py
New dataclasses and create_semisynthetic_config to parse YAML into typed simulation config (paths, features, outcomes, treatment-effect specs).
Simulation core
corebehrt/modules/simulation/semisynthetic_simulator.py
New SemiSyntheticCausalSimulator: two-pass workflow (global feature stats, per-shard simulate), exposure derivation, pre-index filtering, min-code filtering, feature standardization, parametric logistic outcome generation with interactions/noise, constant/heterogeneous treatment effects, returns counterfactuals/ite/event DataFrames and finalization (stats/plots).
Oracle features
corebehrt/modules/simulation/oracle_features.py
New extract_oracle_features and standardize_features: compute counts/burdens/age/recency/burst/motif features using code-prefix filtering and time windows with imputation/warnings.
CLI scripts
corebehrt/main_causal/simulate_semisynthetic.py, corebehrt/main_causal/calibrate_semisynthetic.py
New main_simulate and main_calibrate entrypoints: simulate writes per-outcome CSVs and optionally finalizes; calibrate aggregates shard diagnostics (feature SMDs, P0/P1 stats, ATE/ATT/ATC/RR), saves plots.
Azure components
corebehrt/azure/components/simulate_semisynthetic.py, corebehrt/azure/components/calibrate_semisynthetic.py
Azure wrappers declaring INPUTS/OUTPUTS and wiring CLI functions via job.run_main.
Tests
tests/test_modules/test_simulation/test_oracle_features.py, tests/test_modules/test_simulation/test_semisynthetic.py
Unit tests for oracle feature extraction (features, standardization, motifs) and end-to-end simulator behavior (output format, reproducibility, probability ranges, treatment-effect behavior, exposure derivation).
Docs & utilities
experiments/semisynthetic_simulation/*, experiments/semisynthetic_simulation/generate_configs.py
New experiment README and docs (local/Azure/calibration/features/multiple-runs) and generate_configs.py to produce per-run YAMLs for multi-seed experiments.

Sequence Diagram(s)

sequenceDiagram
    participant User
    participant CLI as simulate_semisynthetic.py
    participant ShardLoader as ShardLoader
    participant Simulator as SemiSyntheticCausalSimulator
    participant Features as extract_oracle_features()
    participant Storage as OutcomesDir

    User->>CLI: run simulate CLI (config)
    CLI->>ShardLoader: instantiate with cfg.paths
    loop per shard
        ShardLoader->>CLI: provide shard DataFrame
        CLI->>Simulator: compute_global_feature_stats / simulate_dataset
        Simulator->>Simulator: determine exposure, filter pre-index events
        Simulator->>Features: compute oracle features (history, pids, index_dates)
        Features-->>Simulator: return features DataFrame
        Simulator->>Simulator: compute P0/P1, ITEs, sample outcomes, build DataFrames
        Simulator-->>CLI: return outcome DataFrames
    end
    CLI->>Storage: write per-outcome CSVs and plots
Loading
sequenceDiagram
    participant User
    participant CLI as calibrate_semisynthetic.py
    participant ShardLoader as ShardLoader
    participant Simulator as SemiSyntheticCausalSimulator
    participant Diagnostics as Plotter
    participant Storage as FiguresDir

    User->>CLI: run calibrate CLI (config)
    CLI->>ShardLoader: iterate shards
    loop per shard
        ShardLoader->>Simulator: extract_features_and_probabilities(shard)
        Simulator-->>CLI: features, exposure, P0, P1, ITE
    end
    CLI->>Diagnostics: compute SMDs, P0/P1 stats, ATE/ATT/ATC, RR
    Diagnostics-->>Storage: save histograms, ITE plots, love-plot
Loading

Estimated code review effort

🎯 4 (Complex) | ⏱️ ~60 minutes

Possibly related PRs

  • Separate Causal #34 — Related to reorganized causal constants and CausalDirectoryPreparer used by the new simulator and CLI.
  • Merge/march 23 #38 — Related to Azure job-run infrastructure (job.run_main) invoked by the new Azure component wrappers.

Poem

🐰 I nibble configs beneath Synthea light,

I count motifs and hop from shard to shard,
I sketch P0 and P1 with a carrot's delight,
I leave tiny ITEs where outcomes are starred,
Hooray — simulated carrots brighten the yard! 🥕✨

🚥 Pre-merge checks | ✅ 2 | ❌ 1

❌ Failed checks (1 warning)

Check name Status Explanation Resolution
Docstring Coverage ⚠️ Warning Docstring coverage is 38.46% which is insufficient. The required threshold is 80.00%. Write docstrings for the functions missing them to satisfy the coverage threshold.
✅ Passed checks (2 passed)
Check name Status Explanation
Description Check ✅ Passed Check skipped - CodeRabbit’s high-level summary is enabled.
Title check ✅ Passed The title 'add semi-synthetic causal simulator with observed treatment' accurately describes the main change: adding a new simulator that uses observed treatment assignments from real data alongside simulated outcomes.

✏️ Tip: You can configure your own custom pre-merge checks in the settings.

✨ Finishing Touches
🧪 Generate unit tests (beta)
  • Create PR with unit tests
  • Commit unit tests in branch improved-simulation

Thanks for using CodeRabbit! It's free for OSS, and your support helps us grow. If you like it, consider giving us a shout-out.

❤️ Share

Comment @coderabbitai help to get the list of available commands.

The updated simulation description drops the r_B/r_L feature grouping.
The outcome model is now eta^(0) = beta_0 + f(r_i) with a single flat
feature vector. Merges baseline_coefficients and longitudinal_coefficients
into a single coefficients dict across config, simulator, YAML, and tests.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
@kirilklein
kirilklein requested a review from Montgomeryyyy April 3, 2026 09:30

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 6

🧹 Nitpick comments (1)
corebehrt/main_causal/calibrate_semisynthetic.py (1)

50-72: Pull the shared preprocessing/math into a public API.

This script reaches into five _... methods on SemiSyntheticCausalSimulator. That hidden coupling means a refactor of the simulator internals can break calibration even when the public simulator behavior is unchanged. A shared public helper/module would keep both entrypoints aligned.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/main_causal/calibrate_semisynthetic.py` around lines 50 - 72, The
calibration loop in calibrate_semisynthetic.py depends on five private methods
on SemiSyntheticCausalSimulator (_extract_treatment_and_index_dates,
_filter_to_pre_index, _apply_min_num_codes, _compute_eta_0, _compute_tau),
coupling this script to simulator internals; add a public API on
SemiSyntheticCausalSimulator (e.g., a method named prepare_calibration_batch or
compute_preindex_features) that encapsulates the preprocessing and math
currently performed by those five methods and returns (features_df, pids,
is_exposed, index_dates, eta_0, tau) or separate logical pieces, then update
calibrate_semisynthetic.py to call that new public method instead of the private
helpers so the calibration code uses the simulator's stable public interface.
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Inline comments:
In `@corebehrt/main_causal/calibrate_semisynthetic.py`:
- Around line 70-74: The calibration currently computes p0/p1 as expit(eta_0)
and expit(eta_0 + tau) but the simulator adds a noise term before the sigmoid;
update the calibration to match the simulator by incorporating the same noise
model (use the simulator or sim_config noise_scale used in
corebehrt/modules/simulation/semisynthetic_simulator.py when computing
probabilities) so p0 = expit(eta_0 + noise) and p1 = expit(eta_0 + tau + noise),
or alternatively explicitly set noise_scale=0 for this script; reference eta_0,
tau, expit, simulator._compute_eta_0, simulator._compute_tau and the
simulator/sim_config noise_scale to locate where to change.

In `@corebehrt/main_causal/simulate_semisynthetic.py`:
- Around line 43-55: The loop buffers all shard DataFrames into
simulated_outcomes causing unbounded memory growth; modify the shard processing
loop (the shard_loader() iteration and use of simulator.simulate_dataset) to
stream each non-empty DataFrame to disk immediately instead of appending to
simulated_outcomes: for each key returned by simulate_dataset, open/join the
target file in outcomes_dir (use append mode) and write the DataFrame, writing
the header only for the first write for that key (track a simple set of written
keys), and remove the defaultdict(list)/post-concat code so you no longer
accumulate df_list for keys like simulated_outcomes[k]; ensure behavior for
multiple shards per key remains identical (rows appended in same CSV).

In `@corebehrt/modules/simulation/oracle_features.py`:
- Around line 164-175: In _compute_age, handle the case where no BIRTH_CODE rows
exist: compute ages for patients with DOB as you already do, then compute
mean_age from the non-null entries of age_series (or dob_per_patient-derived
ages); if that mean is NaN (i.e., dob_per_patient is empty), set mean_age to a
sensible fallback constant (e.g., 50.0) before calling
age_series.fillna(mean_age) so the returned Series contains no nulls; reference
symbols: _compute_age, history_df, BIRTH_CODE, dob_per_patient, age_series,
pids, index_dates.
- Around line 253-257: The _standardize function uses DataFrame.std() with
default ddof=1 which yields NaN for single-row DataFrames and the current
stds.replace(0, 1) does not fix NaNs; update _standardize to compute stds with
ddof=0 (features_df.std(ddof=0)) and then replace any remaining NaN or zero
values with 1 via fillna(1) (or equivalent) before dividing so single-patient
shards and zero-variance columns produce valid standardized features.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 145-154: The filtering uses a non-strict cutoff (df =
df[df[TIMESTAMP_COL] <= patient_index]) which lets events at the index timestamp
through; change this to a strict pre-index cutoff by using < instead of <= when
comparing df[TIMESTAMP_COL] to patient_index (refer to the variables
patient_index, TIMESTAMP_COL and index_dates in semisynthetic_simulator.py) and,
if same-timestamp inclusion was intentional, update the docstring to state that
events at the index time are considered pre-index instead of changing the
operator.
- Around line 303-339: The per-shard simulate_dataset flow currently writes
simulation_stats, theoretical_max_roc_auc and figures inside
semisynthetic_simulator.SemisyntheticSimulator (calls to
_calculate_and_save_simulation_stats, _calculate_theoretical_roc_auc,
plot_probability_distributions, plot_true_effects_vs_risk_differences), causing
overwrites when simulate_dataset() is run for each shard; move all stats/plot
generation out of the per-shard path by removing those calls from
simulate_dataset and instead implement a dataset-level finalizer (e.g.,
SemisyntheticSimulator.finalize_simulation_outputs or
generate_dataset_level_reports) that accepts the aggregated ite_records,
cf_records, all_probas, pids and is_exposed and performs the calls to
_calculate_and_save_simulation_stats, _calculate_theoretical_roc_auc,
plot_probability_distributions and plot_true_effects_vs_risk_differences; update
simulate_dataset to return the raw records (ite_records, cf_records, all_probas,
pids, is_exposed) so the caller (simulate_semisynthetic.py) can aggregate across
shards and invoke the new finalizer exactly once after aggregation.

---

Nitpick comments:
In `@corebehrt/main_causal/calibrate_semisynthetic.py`:
- Around line 50-72: The calibration loop in calibrate_semisynthetic.py depends
on five private methods on SemiSyntheticCausalSimulator
(_extract_treatment_and_index_dates, _filter_to_pre_index, _apply_min_num_codes,
_compute_eta_0, _compute_tau), coupling this script to simulator internals; add
a public API on SemiSyntheticCausalSimulator (e.g., a method named
prepare_calibration_batch or compute_preindex_features) that encapsulates the
preprocessing and math currently performed by those five methods and returns
(features_df, pids, is_exposed, index_dates, eta_0, tau) or separate logical
pieces, then update calibrate_semisynthetic.py to call that new public method
instead of the private helpers so the calibration code uses the simulator's
stable public interface.
🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

  • Push a commit to this branch (recommended)
  • Create a new PR with the fixes

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro

Run ID: 13722394-db16-48a2-ba35-5de4e654efed

📥 Commits

Reviewing files that changed from the base of the PR and between 545e880 and 1246ea5.

📒 Files selected for processing (9)
  • corebehrt/configs/causal/simulate_semisynthetic.yaml
  • corebehrt/main_causal/calibrate_semisynthetic.py
  • corebehrt/main_causal/simulate_semisynthetic.py
  • corebehrt/modules/simulation/config_semisynthetic.py
  • corebehrt/modules/simulation/oracle_features.py
  • corebehrt/modules/simulation/semisynthetic_simulator.py
  • tests/test_modules/test_simulation/__init__.py
  • tests/test_modules/test_simulation/test_oracle_features.py
  • tests/test_modules/test_simulation/test_semisynthetic.py

Comment thread corebehrt/main_causal/calibrate_semisynthetic.py Outdated
Comment thread corebehrt/main_causal/simulate_semisynthetic.py
Comment thread corebehrt/modules/simulation/oracle_features.py Outdated
Comment thread corebehrt/modules/simulation/oracle_features.py Outdated
Comment thread corebehrt/modules/simulation/semisynthetic_simulator.py
Comment thread corebehrt/modules/simulation/semisynthetic_simulator.py
Bug fixes from CodeRabbit review:
- Handle all-missing-DOB with default age fallback (65 years)
- Fix single-patient standardization: use ddof=0 to avoid NaN
- Move stats/plots to finalize() method called after shard aggregation
- Add clarifying comment re intentional noise omission in calibration

New files:
- Azure components for simulate and calibrate
- experiments/semisynthetic_simulation/ docs and config generator

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 3

♻️ Duplicate comments (3)
corebehrt/modules/simulation/semisynthetic_simulator.py (2)

145-149: ⚠️ Potential issue | 🟠 Major

Use a strict pre-index cutoff.

<= includes same-timestamp events even though the docstring says "before". If treatment is assigned at index time, diagnoses/procedures recorded at that timestamp leak contemporaneous information into the oracle features. Use < here, or explicitly document that index-time events are considered pre-index.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 145 -
149, The cutoff currently uses a non-strict comparison (df =
df[df[TIMESTAMP_COL] <= patient_index]) which includes events at the exact index
timestamp; change it to a strict pre-index filter by using df =
df[df[TIMESTAMP_COL] < patient_index] (keep the same PID_COL, TIMESTAMP_COL and
patient_index variables) so events occurring at index time are excluded from the
oracle features.

297-335: ⚠️ Potential issue | 🟠 Major

Generate stats and plots once after shard aggregation.

simulate_dataset() is shard-scoped, but this block rewrites simulation_stats.csv, theoretical_max_roc_auc.csv, and all figures on every call. The event/counterfactual outputs can be aggregated across shards while the diagnostics only describe the last shard processed. Move these side effects to a dataset-level finalize step.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 297 -
335, The diagnostics and plotting block in simulate_dataset (calls to
_calculate_and_save_simulation_stats, _calculate_theoretical_roc_auc,
plot_probability_distributions, and plot_true_effects_vs_risk_differences) must
be moved out of the shard-scoped simulate_dataset and executed once after shard
aggregation; remove the block that writes simulation_stats.csv,
theoretical_max_roc_auc.csv, and creates figs from simulate_dataset, and instead
create or use a dataset-level finalize method (e.g., finalize_dataset or
similar) that runs after all shards are merged and receives aggregated pids,
is_exposed, cf_records, all_probas, ite_records and config; call
_calculate_and_save_simulation_stats(pids, is_exposed, cf_records, output_dir),
_calculate_theoretical_roc_auc(cf_records, is_exposed, output_dir), and the
plotting helpers (plot_probability_distributions and
plot_true_effects_vs_risk_differences) there, ensuring true_effects_config is
built from self.config.outcomes as shown so diagnostics reflect the full dataset
not the last shard.
corebehrt/main_causal/calibrate_semisynthetic.py (1)

70-74: ⚠️ Potential issue | 🟠 Major

Calibration is using a different probability model than simulation.

SemiSyntheticCausalSimulator._simulate_outcomes() adds noise_scale before expit, but this script reports expit(eta_0) and expit(eta_0 + tau) directly. Whenever noise_scale > 0, the prevalence and causal diagnostics here won't match the probabilities that generated the sampled outcomes. Reuse the simulator's noise model here, or explicitly require noise_scale == 0 for calibration runs.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/main_causal/calibrate_semisynthetic.py` around lines 70 - 74,
Calibration currently computes p0 = expit(eta_0) and p1 = expit(eta_0 + tau) but
SemiSyntheticCausalSimulator._simulate_outcomes() applies the simulator's
noise_scale before the expit, so prevalences will mismatch when noise_scale > 0;
fix by either (A) reusing the simulator's noise model when computing
probabilities here (i.e., apply the same noise-generation logic used in
_simulate_outcomes to eta_0 and eta_0 + tau before calling expit, using
simulator.noise_scale or sim_config.noise_scale and the same RNG), or (B)
enforce simulator.noise_scale == 0 for calibration runs by adding an explicit
assertion/guard that raises if noise_scale > 0 so the current p0/p1 are valid.
Ensure you reference SemiSyntheticCausalSimulator._simulate_outcomes,
_compute_eta_0, _compute_tau, and the noise_scale setting when making the
change.
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Inline comments:
In `@corebehrt/main_causal/calibrate_semisynthetic.py`:
- Around line 64-67: The current code calls extract_oracle_features (which
standardizes per shard), appends features_df to all_features, then computes
global diagnostics, which mixes shard-local z-scores; instead either (A) disable
standardization inside extract_oracle_features for calibration mode (add a flag
like standardize=False or use sim_config.calibrate_no_std) and return raw
features for concatenation, or (B) collect raw unstandardized features from each
shard (use the raw output from extract_oracle_features or a new function),
concatenate into all_features, compute a single global scaler (mean/std) and
apply that global transform before computing quantiles/SMDs and the love plot;
update any calls between lines 85-110 that depend on per-shard standardization
to use the global-scaled features or the no-standardization path. Ensure
references: extract_oracle_features, features_df, all_features, and
sim_config.features are updated accordingly.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 254-263: The loop silently ignores misspelled feature names and
invalid treatment effect modes, altering the DGP; update validation in
semisynthetic_simulator.py to (1) check outcome_model.coefficients keys and each
interaction["features"] and heterogeneous_coefficients references against
features_df.columns and raise a clear ValueError listing missing feature names,
(2) validate treatment_effect.mode explicitly (allow only "constant" and any
documented modes) and raise on unknown values instead of falling through to
heterogeneous logic, and (3) apply the same presence/validation logic used here
to the analogous block around heterogeneous_coefficients (lines ~274-281) so
typos or bad config fail fast with informative errors mentioning the offending
symbol names (e.g., outcome_model.interactions, heterogeneous_coefficients,
treatment_effect.mode).
- Around line 84-86: The per-shard call to extract_oracle_features(history_df,
pids, index_dates, self.config.features) allows FeatureConfig.standardize=True
to z-score each shard independently, causing shard-dependent logits; fix by
computing a single global mean/std for the entire cohort once (e.g., call
extract_oracle_features over the full cohort or add a helper to compute global
scaler from FeatureConfig) and then apply that same scaler inside the shard loop
(or set FeatureConfig.standardize=False for per-shard calls and explicitly
normalize each shard using the precomputed global mean/std); update calls to
extract_oracle_features (and/or FeatureConfig) accordingly so features_df for
every shard is standardized with the same global statistics.

---

Duplicate comments:
In `@corebehrt/main_causal/calibrate_semisynthetic.py`:
- Around line 70-74: Calibration currently computes p0 = expit(eta_0) and p1 =
expit(eta_0 + tau) but SemiSyntheticCausalSimulator._simulate_outcomes() applies
the simulator's noise_scale before the expit, so prevalences will mismatch when
noise_scale > 0; fix by either (A) reusing the simulator's noise model when
computing probabilities here (i.e., apply the same noise-generation logic used
in _simulate_outcomes to eta_0 and eta_0 + tau before calling expit, using
simulator.noise_scale or sim_config.noise_scale and the same RNG), or (B)
enforce simulator.noise_scale == 0 for calibration runs by adding an explicit
assertion/guard that raises if noise_scale > 0 so the current p0/p1 are valid.
Ensure you reference SemiSyntheticCausalSimulator._simulate_outcomes,
_compute_eta_0, _compute_tau, and the noise_scale setting when making the
change.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 145-149: The cutoff currently uses a non-strict comparison (df =
df[df[TIMESTAMP_COL] <= patient_index]) which includes events at the exact index
timestamp; change it to a strict pre-index filter by using df =
df[df[TIMESTAMP_COL] < patient_index] (keep the same PID_COL, TIMESTAMP_COL and
patient_index variables) so events occurring at index time are excluded from the
oracle features.
- Around line 297-335: The diagnostics and plotting block in simulate_dataset
(calls to _calculate_and_save_simulation_stats, _calculate_theoretical_roc_auc,
plot_probability_distributions, and plot_true_effects_vs_risk_differences) must
be moved out of the shard-scoped simulate_dataset and executed once after shard
aggregation; remove the block that writes simulation_stats.csv,
theoretical_max_roc_auc.csv, and creates figs from simulate_dataset, and instead
create or use a dataset-level finalize method (e.g., finalize_dataset or
similar) that runs after all shards are merged and receives aggregated pids,
is_exposed, cf_records, all_probas, ite_records and config; call
_calculate_and_save_simulation_stats(pids, is_exposed, cf_records, output_dir),
_calculate_theoretical_roc_auc(cf_records, is_exposed, output_dir), and the
plotting helpers (plot_probability_distributions and
plot_true_effects_vs_risk_differences) there, ensuring true_effects_config is
built from self.config.outcomes as shown so diagnostics reflect the full dataset
not the last shard.
🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

  • Push a commit to this branch (recommended)
  • Create a new PR with the fixes

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro

Run ID: b8da52c3-246d-4ac9-9bc3-fef0784b9c8c

📥 Commits

Reviewing files that changed from the base of the PR and between 1246ea5 and bcbbd41.

📒 Files selected for processing (7)
  • corebehrt/configs/causal/simulate_semisynthetic.yaml
  • corebehrt/main_causal/calibrate_semisynthetic.py
  • corebehrt/modules/simulation/config_semisynthetic.py
  • corebehrt/modules/simulation/oracle_features.py
  • corebehrt/modules/simulation/semisynthetic_simulator.py
  • tests/test_modules/test_simulation/test_oracle_features.py
  • tests/test_modules/test_simulation/test_semisynthetic.py
✅ Files skipped from review due to trivial changes (2)
  • corebehrt/configs/causal/simulate_semisynthetic.yaml
  • tests/test_modules/test_simulation/test_semisynthetic.py
🚧 Files skipped from review as they are similar to previous changes (2)
  • corebehrt/modules/simulation/oracle_features.py
  • tests/test_modules/test_simulation/test_oracle_features.py

Comment thread corebehrt/main_causal/calibrate_semisynthetic.py Outdated
Comment thread corebehrt/modules/simulation/semisynthetic_simulator.py Outdated
Comment thread corebehrt/modules/simulation/semisynthetic_simulator.py
- Vectorize _compute_age, _compute_event_recency, and
  _compute_sequence_motif_count using pandas merge/reindex instead
  of Python loops over patients
- Use None sentinel instead of "" for prefix in
  _filter_by_prefix_and_window
- Add public extract_features_and_probabilities() on simulator so
  calibration script no longer calls private methods

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
@kirilklein

Copy link
Copy Markdown
Owner Author

@Montgomeryyyy for context please read the description of simulation in the overleaf

Features were z-scored per shard independently, making simulated
probabilities depend on which shard a patient lands in. Now the
simulator does two passes: pass 1 computes global mean/std across
all shards, pass 2 standardizes using those global stats before
computing logits.

- extract_oracle_features now returns raw (unstandardized) features
- New standardize_features(df, means, stds) function for explicit control
- compute_global_feature_stats() iterates shards once to collect stats
- Entry points (simulate + calibrate) call pass 1 then pass 2
- Removed standardize flag from FeatureConfig (always standardize globally)

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 2

♻️ Duplicate comments (3)
corebehrt/modules/simulation/semisynthetic_simulator.py (3)

184-188: ⚠️ Potential issue | 🟠 Major

Use a strict pre-index cutoff.

Line 187 still includes same-timestamp events even though the docstring says "before" the index date. If treatment is assigned at index time, those events leak contemporaneous information into the oracle features.

🛠️ Suggested fix
-        df = df[df[TIMESTAMP_COL] <= patient_index]
+        df = df[df[TIMESTAMP_COL] < patient_index]
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 184 -
188, The current filtering uses df = df[df[TIMESTAMP_COL] <= patient_index]
which allows events exactly at the index date to leak; change the cutoff to
strict pre-index by using a strict less-than comparison (use df[TIMESTAMP_COL] <
patient_index) while keeping the same patient_index construction via
index_dates.reindex(df[PID_COL]).values and the prior PID_COL membership
filtering so alignment/broadcasting remains correct.

291-300: ⚠️ Potential issue | 🟠 Major

Fail fast on bad feature names and unsupported treatment modes.

Lines 292 and 317 silently drop misspelled features, and anything other than "constant" currently falls into the heterogeneous branch. A YAML typo will quietly change the DGP instead of failing.

🛠️ Suggested fix
+        missing = sorted(set(outcome_model.coefficients) - set(features_df.columns))
+        if missing:
+            raise ValueError(
+                f"Unknown outcome_model.coefficients features: {missing}"
+            )
         for name, coeff in outcome_model.coefficients.items():
-            if name in features_df.columns:
-                eta += coeff * features_df[name].values
+            eta += coeff * features_df[name].values
@@
         if treatment_effect.mode == "constant":
             return np.full(n, treatment_effect.delta)
+        if treatment_effect.mode != "heterogeneous":
+            raise ValueError(
+                f"Unsupported treatment_effect.mode: {treatment_effect.mode}"
+            )

Please apply the same validation pattern to interaction feature pairs and to heterogeneous_coefficients.

Also applies to: 311-319

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 291 -
300, The code currently silently ignores misspelled or missing feature names and
treats any non-"constant" mode as heterogeneous; update validation so that when
building eta you check that every key in outcome_model.coefficients exists in
features_df.columns and raise a clear exception listing missing names (same for
outcome_model.interactions: verify both feature names from
interaction["features"] exist and raise if not), and apply the same presence
checks to heterogeneous_coefficients; additionally validate the treatment mode
string explicitly (accept only "constant" or "heterogeneous" and raise on
anything else) so a YAML typo fails fast; reference outcome_model.coefficients,
outcome_model.interactions, heterogeneous_coefficients, features_df, and the
treatment mode handling when adding these checks.

84-86: ⚠️ Potential issue | 🟠 Major

Shard-local standardization still makes the DGP shard-dependent.

extract_oracle_features() standardizes inside the call, so invoking it once per shard here means identical raw histories can get different z-scores/logits depending on shard composition. Calibration then concatenates incompatible scales. Compute one cohort-level scaler and reuse it for every shard, or disable standardization in the shard loop.

Also applies to: 131-133

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 84 -
86, extract_oracle_features is being called per-shard which applies
standardization per shard and makes the DGP shard-dependent; instead, fit a
single cohort-level scaler once and reuse it for all shards. Concretely: call
extract_oracle_features (or a one-shot feature extraction variant) over the
entire cohort to obtain raw features and fit a scaler (or obtain the scaler if
extract_oracle_features can return it), store that scaler as cohort_scaler, then
in the shard loop replace per-shard calls that produce features_df (the calls
using history_df, pids, index_dates, self.config.features) with extraction that
does not standardize and apply cohort_scaler.transform (or pass the
cohort_scaler into extract_oracle_features if you add a parameter) so every
shard uses the same global z-score/logit scaling; apply the same change to the
other occurrence where features are extracted (the second
extract_oracle_features call around the later shard handling).
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Inline comments:
In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 519-539: The current code computes auc_if_all_treated and
auc_if_all_control by scoring p_treated and p_control against the factual label
y_true (from outcome_col); change this to score p_treated against the simulated
exposed label and p_control against the simulated control label (the CF record
columns named like SIMULATED_OUTCOME_EXPOSED_* and SIMULATED_OUTCOME_CONTROL_*
in cf_records) when calling roc_auc_score so auc_treated and auc_control reflect
the correct counterfactual targets; keep auc_factual_dgp computed against
outcome_col as before and still guard with the existing len(np.unique(...))
check.

In `@experiments/semisynthetic_simulation/docs/multiple_runs.md`:
- Around line 18-24: Update the two unannotated code fences that show the ASCII
tree structures so they declare the language as plain text: change the opening
``` to ```text for the snippet that begins with "generated_configs/" and for the
snippet that begins with "outputs/causal/semisynthetic_study/runs/". This will
satisfy markdownlint MD040 without changing the displayed content.

---

Duplicate comments:
In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 184-188: The current filtering uses df = df[df[TIMESTAMP_COL] <=
patient_index] which allows events exactly at the index date to leak; change the
cutoff to strict pre-index by using a strict less-than comparison (use
df[TIMESTAMP_COL] < patient_index) while keeping the same patient_index
construction via index_dates.reindex(df[PID_COL]).values and the prior PID_COL
membership filtering so alignment/broadcasting remains correct.
- Around line 291-300: The code currently silently ignores misspelled or missing
feature names and treats any non-"constant" mode as heterogeneous; update
validation so that when building eta you check that every key in
outcome_model.coefficients exists in features_df.columns and raise a clear
exception listing missing names (same for outcome_model.interactions: verify
both feature names from interaction["features"] exist and raise if not), and
apply the same presence checks to heterogeneous_coefficients; additionally
validate the treatment mode string explicitly (accept only "constant" or
"heterogeneous" and raise on anything else) so a YAML typo fails fast; reference
outcome_model.coefficients, outcome_model.interactions,
heterogeneous_coefficients, features_df, and the treatment mode handling when
adding these checks.
- Around line 84-86: extract_oracle_features is being called per-shard which
applies standardization per shard and makes the DGP shard-dependent; instead,
fit a single cohort-level scaler once and reuse it for all shards. Concretely:
call extract_oracle_features (or a one-shot feature extraction variant) over the
entire cohort to obtain raw features and fit a scaler (or obtain the scaler if
extract_oracle_features can return it), store that scaler as cohort_scaler, then
in the shard loop replace per-shard calls that produce features_df (the calls
using history_df, pids, index_dates, self.config.features) with extraction that
does not standardize and apply cohort_scaler.transform (or pass the
cohort_scaler into extract_oracle_features if you add a parameter) so every
shard uses the same global z-score/logit scaling; apply the same change to the
other occurrence where features are extracted (the second
extract_oracle_features call around the later shard handling).
🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

  • Push a commit to this branch (recommended)
  • Create a new PR with the fixes

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro

Run ID: 55b6c13f-e849-4a40-8f0c-1cd88060db3d

📥 Commits

Reviewing files that changed from the base of the PR and between bcbbd41 and 187bda6.

📒 Files selected for processing (13)
  • corebehrt/azure/components/calibrate_semisynthetic.py
  • corebehrt/azure/components/simulate_semisynthetic.py
  • corebehrt/main_causal/calibrate_semisynthetic.py
  • corebehrt/main_causal/simulate_semisynthetic.py
  • corebehrt/modules/simulation/oracle_features.py
  • corebehrt/modules/simulation/semisynthetic_simulator.py
  • experiments/semisynthetic_simulation/README.md
  • experiments/semisynthetic_simulation/docs/azure.md
  • experiments/semisynthetic_simulation/docs/calibration.md
  • experiments/semisynthetic_simulation/docs/features.md
  • experiments/semisynthetic_simulation/docs/local.md
  • experiments/semisynthetic_simulation/docs/multiple_runs.md
  • experiments/semisynthetic_simulation/generate_configs.py
✅ Files skipped from review due to trivial changes (6)
  • experiments/semisynthetic_simulation/README.md
  • experiments/semisynthetic_simulation/docs/features.md
  • corebehrt/azure/components/simulate_semisynthetic.py
  • experiments/semisynthetic_simulation/docs/local.md
  • experiments/semisynthetic_simulation/docs/calibration.md
  • experiments/semisynthetic_simulation/docs/azure.md
🚧 Files skipped from review as they are similar to previous changes (2)
  • corebehrt/modules/simulation/oracle_features.py
  • corebehrt/main_causal/simulate_semisynthetic.py

Comment thread corebehrt/modules/simulation/semisynthetic_simulator.py
Comment on lines +18 to +24
```
generated_configs/
├── my_scenario_run_01.yaml # seed=43, outcomes → runs/run_01/my_scenario/
├── my_scenario_run_02.yaml # seed=44, outcomes → runs/run_02/my_scenario/
├── ...
└── my_scenario_run_50.yaml # seed=92, outcomes → runs/run_50/my_scenario/
```

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🟡 Minor

Add a language to the tree-style code fences.

Lines 18 and 53 currently trip markdownlint's MD040 warning. Marking these snippets as plain text keeps the doc lint-clean.

🛠️ Suggested fix
-```
+```text
 generated_configs/
 ├── my_scenario_run_01.yaml    # seed=43, outcomes → runs/run_01/my_scenario/
 ├── my_scenario_run_02.yaml    # seed=44, outcomes → runs/run_02/my_scenario/
 ├── ...
 └── my_scenario_run_50.yaml    # seed=92, outcomes → runs/run_50/my_scenario/

@@
- +text
outputs/causal/semisynthetic_study/runs/
├── run_01/my_scenario/
│ ├── counterfactuals.csv
│ ├── ite.csv
│ └── ...
├── run_02/my_scenario/
│ └── ...
└── run_50/my_scenario/
└── ...

Also applies to: 53-63

🧰 Tools
🪛 markdownlint-cli2 (0.22.0)

[warning] 18-18: Fenced code blocks should have a language specified

(MD040, fenced-code-language)

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@experiments/semisynthetic_simulation/docs/multiple_runs.md` around lines 18 -
24, Update the two unannotated code fences that show the ASCII tree structures
so they declare the language as plain text: change the opening ``` to ```text
for the snippet that begins with "generated_configs/" and for the snippet that
begins with "outputs/causal/semisynthetic_study/runs/". This will satisfy
markdownlint MD040 without changing the displayed content.

kirilklein and others added 2 commits April 3, 2026 12:26
auc_if_all_treated now uses Y(1) as labels (not factual Y), and
auc_if_all_control uses Y(0). Previously both scored against the
factual outcome which mixes treated and untreated labels.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
The auc_if_all_treated/control metrics intentionally score P(Y(1)) and
P(Y(0)) against the factual observed outcome, matching the existing
RealisticCausalSimulator. This measures how well counterfactual risk
models discriminate observed outcomes — a more informative diagnostic
than scoring against counterfactual labels (which would be near-perfect
since Y(a) ~ Bernoulli(P(Y(a)))).

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 3

♻️ Duplicate comments (2)
corebehrt/modules/simulation/semisynthetic_simulator.py (2)

224-230: ⚠️ Potential issue | 🟠 Major

Use a strict pre-index cutoff here.

The docstring says "before" the index date, but <= still keeps same-timestamp non-exposure events. That leaks contemporaneous information into the oracle features whenever other codes land exactly at index time.

🔒 Suggested fix
-        df = df[df[TIMESTAMP_COL] <= patient_index]
+        df = df[df[TIMESTAMP_COL] < patient_index]
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 224 -
230, In _filter_to_pre_index ensure the cutoff is strict "before" the index
date: in semisynthetic_simulator.py inside method _filter_to_pre_index replace
the current comparison df[TIMESTAMP_COL] <= patient_index with a strict
less-than df[TIMESTAMP_COL] < patient_index so events with timestamp equal to
the index date are excluded; keep the existing patient_index reindexing logic
and DataFrame filtering otherwise unchanged.

334-343: ⚠️ Potential issue | 🟠 Major

Fail fast on bad outcome-model config instead of silently changing the DGP.

Misspelled feature names in coefficients, interactions, or heterogeneous_coefficients are ignored, and any mode other than "constant" currently falls through to the heterogeneous branch. A YAML typo should raise here, not change the simulator behavior quietly.

🛡️ Suggested validation sketch
     def _compute_eta_0(
         self,
         features_df: pd.DataFrame,
         outcome_model: OutcomeModelConfig,
     ) -> np.ndarray:
         """Compute the baseline log-odds: beta_0 + f(r_i) + interactions."""
         n = len(features_df)
         eta = np.full(n, outcome_model.beta_0)
+        valid_features = set(features_df.columns)
+        missing = set(outcome_model.coefficients) - valid_features
+        missing.update(
+            feat
+            for interaction in outcome_model.interactions
+            for feat in interaction.get("features", [])
+            if feat not in valid_features
+        )
+        if missing:
+            raise ValueError(f"Unknown outcome-model features: {sorted(missing)}")
@@
     def _compute_tau(
         self,
         features_df: pd.DataFrame,
         treatment_effect: TreatmentEffectConfig,
     ) -> np.ndarray:
         """Compute individual treatment effects on the logit scale."""
         n = len(features_df)
         if treatment_effect.mode == "constant":
             return np.full(n, treatment_effect.delta)
+        if treatment_effect.mode != "heterogeneous":
+            raise ValueError(
+                f"Unsupported treatment_effect.mode: {treatment_effect.mode}"
+            )
+        missing = (
+            set(treatment_effect.heterogeneous_coefficients) - set(features_df.columns)
+        )
+        if missing:
+            raise ValueError(
+                f"Unknown heterogeneous effect features: {sorted(missing)}"
+            )
 
         # heterogeneous mode
         tau = np.full(n, treatment_effect.delta_0)

Also applies to: 354-361

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 334 -
343, The simulator currently ignores misspelled feature names and silently
treats unknown outcome modes as heterogeneous; update validation in the code
that builds eta (use symbols outcome_model.coefficients,
outcome_model.interactions, outcome_model.heterogeneous_coefficients,
features_df, and eta) to fail fast: before computing eta, verify that every
feature name referenced in outcome_model.coefficients and in each
interaction["features"] exists in features_df.columns and raise a clear
ValueError if any are missing, and validate outcome_model.mode against allowed
values (e.g., "constant" and the documented alternatives) raising an error on
unknown modes instead of falling through to the heterogeneous branch; apply the
same presence checks/validation for the heterogeneous_coefficients handling
referenced later (addresses the similar block around lines 354-361).
🧹 Nitpick comments (1)
tests/test_modules/test_simulation/test_semisynthetic.py (1)

114-123: Please cover the new pass-1 standardization workflow in this suite.

These tests call simulate_dataset() directly, so they never exercise compute_global_feature_stats() or assert shard-invariant z-scoring. Add at least one multi-shard test that runs pass 1 first and checks the standardized features or deterministic P0/P1 outputs from extract_features_and_probabilities() are unchanged by shard boundaries.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@tests/test_modules/test_simulation/test_semisynthetic.py` around lines 114 -
123, The current tests call SemiSyntheticCausalSimulator.simulate_dataset()
directly and never exercise the pass-1 standardization path
(compute_global_feature_stats), so add a new multi-shard unit test that: (1)
creates two or more shards via _make_test_shard() (or variants) and a simulator
via SemiSyntheticCausalSimulator(self.config); (2) runs
simulator.compute_global_feature_stats(shards) to perform pass‑1 global
z‑scoring; (3) for each shard call
simulator.extract_features_and_probabilities(shard) (or simulate_dataset(shard)
if that returns the same feature/prob outputs) and capture the standardized
features and deterministic P0/P1; and (4) assert that features and P0/P1 are
identical whether you process shards separately after pass‑1 or if you process
the same data as a single combined shard, ensuring shard boundaries do not
change standardization or deterministic outputs. Ensure you reference
compute_global_feature_stats, extract_features_and_probabilities,
simulate_dataset and use _make_test_shard() to build test inputs.
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Inline comments:
In `@corebehrt/modules/simulation/oracle_features.py`:
- Around line 198-205: In _compute_event_recency make the same finite-fallback
guard used in _compute_age: after computing mean_recency check if it is finite
(not NaN/inf) and if not, replace it with the same fallback constant/value used
by _compute_age before calling recency.fillna(mean_recency); ensure you
reference the symbols last_event, idx_aligned, recency and mean_recency and
apply the fallback so recency.fillna(...) never receives NaN.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 543-605: The file fails CI formatting; run the code formatter
(ruff format) on the modified file(s) so the _calculate_theoretical_roc_auc
function and surrounding code match the project's style rules. Specifically, run
`ruff format` (or your project's formatting command) targeting
corebehrt/modules/simulation/semisynthetic_simulator.py (or the whole repo),
stage the updated file, and re-run tests to ensure the formatting job passes; no
logic changes are required—only formatting adjustments.
- Around line 70-97: The current compute_global_feature_stats in
semisynthetic_simulator.py materializes all per-patient feature DataFrames into
all_features and pd.concat which can blow memory; instead compute and accumulate
per-feature counts, sums and sums-of-squares per shard and reduce them to global
aggregates. Modify compute_global_feature_stats to, for each shard returned by
shard_loader(): call _extract_treatment_and_index_dates, _filter_to_pre_index,
_apply_min_num_codes and extract_oracle_features as now, but do not append the
full features_df—instead compute per-column count (non-null), sum, and sumsq
(features_df.count(), features_df.sum(), (features_df**2).sum()) and add those
into running totals (e.g., global_counts, global_sums, global_sumsq); after the
loop compute global means as global_sums / global_counts and variances as
(global_sumsq / global_counts) - mean**2, then set _global_stds =
sqrt(max(variance, 0)) with zeros replaced by 1 and fillna(1) as before;
preserve the existing empty-cohort warning and the logger.info message (use len
aggregated count) and ensure same behavior for ddof=0 by using the population
variance formula.

---

Duplicate comments:
In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 224-230: In _filter_to_pre_index ensure the cutoff is strict
"before" the index date: in semisynthetic_simulator.py inside method
_filter_to_pre_index replace the current comparison df[TIMESTAMP_COL] <=
patient_index with a strict less-than df[TIMESTAMP_COL] < patient_index so
events with timestamp equal to the index date are excluded; keep the existing
patient_index reindexing logic and DataFrame filtering otherwise unchanged.
- Around line 334-343: The simulator currently ignores misspelled feature names
and silently treats unknown outcome modes as heterogeneous; update validation in
the code that builds eta (use symbols outcome_model.coefficients,
outcome_model.interactions, outcome_model.heterogeneous_coefficients,
features_df, and eta) to fail fast: before computing eta, verify that every
feature name referenced in outcome_model.coefficients and in each
interaction["features"] exists in features_df.columns and raise a clear
ValueError if any are missing, and validate outcome_model.mode against allowed
values (e.g., "constant" and the documented alternatives) raising an error on
unknown modes instead of falling through to the heterogeneous branch; apply the
same presence checks/validation for the heterogeneous_coefficients handling
referenced later (addresses the similar block around lines 354-361).

---

Nitpick comments:
In `@tests/test_modules/test_simulation/test_semisynthetic.py`:
- Around line 114-123: The current tests call
SemiSyntheticCausalSimulator.simulate_dataset() directly and never exercise the
pass-1 standardization path (compute_global_feature_stats), so add a new
multi-shard unit test that: (1) creates two or more shards via
_make_test_shard() (or variants) and a simulator via
SemiSyntheticCausalSimulator(self.config); (2) runs
simulator.compute_global_feature_stats(shards) to perform pass‑1 global
z‑scoring; (3) for each shard call
simulator.extract_features_and_probabilities(shard) (or simulate_dataset(shard)
if that returns the same feature/prob outputs) and capture the standardized
features and deterministic P0/P1; and (4) assert that features and P0/P1 are
identical whether you process shards separately after pass‑1 or if you process
the same data as a single combined shard, ensuring shard boundaries do not
change standardization or deterministic outputs. Ensure you reference
compute_global_feature_stats, extract_features_and_probabilities,
simulate_dataset and use _make_test_shard() to build test inputs.
🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

  • Push a commit to this branch (recommended)
  • Create a new PR with the fixes

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro

Run ID: 541c4efa-b975-44a6-8694-e9ccb4816705

📥 Commits

Reviewing files that changed from the base of the PR and between 187bda6 and 50dd4f5.

📒 Files selected for processing (8)
  • corebehrt/configs/causal/simulate_semisynthetic.yaml
  • corebehrt/main_causal/calibrate_semisynthetic.py
  • corebehrt/main_causal/simulate_semisynthetic.py
  • corebehrt/modules/simulation/config_semisynthetic.py
  • corebehrt/modules/simulation/oracle_features.py
  • corebehrt/modules/simulation/semisynthetic_simulator.py
  • tests/test_modules/test_simulation/test_oracle_features.py
  • tests/test_modules/test_simulation/test_semisynthetic.py
✅ Files skipped from review due to trivial changes (2)
  • corebehrt/configs/causal/simulate_semisynthetic.yaml
  • tests/test_modules/test_simulation/test_oracle_features.py
🚧 Files skipped from review as they are similar to previous changes (1)
  • corebehrt/main_causal/simulate_semisynthetic.py

Comment thread corebehrt/modules/simulation/oracle_features.py
Comment thread corebehrt/modules/simulation/semisynthetic_simulator.py
Comment on lines +543 to +605
def _calculate_theoretical_roc_auc(
self,
cf_records: Dict[str, np.ndarray],
is_exposed: np.ndarray,
output_dir: str,
) -> Dict[str, float]:
theoretical_aucs = {}
results_data = []

for outcome_name in self.config.outcomes:
outcome_col = f"{OUTCOME_COL}_{outcome_name}"
p_exposed_col = f"{SIMULATED_PROBAS_EXPOSED}_{outcome_name}"
p_control_col = f"{SIMULATED_PROBAS_CONTROL}_{outcome_name}"

if outcome_col not in cf_records:
continue
if p_exposed_col not in cf_records or p_control_col not in cf_records:
continue

y_true = cf_records[outcome_col]
y1_col = f"{SIMULATED_OUTCOME_EXPOSED}_{outcome_name}"
y0_col = f"{SIMULATED_OUTCOME_CONTROL}_{outcome_name}"
p_treated = cf_records[p_exposed_col]
p_control = cf_records[p_control_col]
y_prob_factual = np.where(is_exposed, p_treated, p_control)

if len(np.unique(y_true)) > 1:
auc_factual = roc_auc_score(y_true, y_prob_factual)
# Score counterfactual probabilities against their matching labels
y1 = cf_records[y1_col]
y0 = cf_records[y0_col]
auc_treated = roc_auc_score(y1, p_treated) if len(np.unique(y1)) > 1 else np.nan
auc_control = roc_auc_score(y0, p_control) if len(np.unique(y0)) > 1 else np.nan
theoretical_aucs[outcome_name] = auc_factual

results_data.append(
{
"outcome": outcome_name,
"auc_factual_dgp": auc_factual,
"auc_if_all_treated": auc_treated,
"auc_if_all_control": auc_control,
"n_positive": int(np.sum(y_true)),
"n_total": len(y_true),
"prevalence": np.mean(y_true),
}
)
logger.info(
f"{outcome_name}: Theoretical max ROC AUC = {auc_factual:.4f}"
)
else:
logger.warning(
f"Cannot calculate ROC AUC for {outcome_name}: only one class present."
)
theoretical_aucs[outcome_name] = np.nan

if results_data:
results_df = pd.DataFrame(results_data)
os.makedirs(output_dir, exist_ok=True)
results_path = join(output_dir, "theoretical_max_roc_auc.csv")
results_df.to_csv(results_path, index=False)
logger.info(f"Theoretical ROC AUC results saved to {results_path}")

return theoretical_aucs

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🟡 Minor

Run ruff format on this file before merge.

The format job is currently red on corebehrt/modules/simulation/semisynthetic_simulator.py, so CI will stay failing until this file is reformatted.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 543 -
605, The file fails CI formatting; run the code formatter (ruff format) on the
modified file(s) so the _calculate_theoretical_roc_auc function and surrounding
code match the project's style rules. Specifically, run `ruff format` (or your
project's formatting command) targeting
corebehrt/modules/simulation/semisynthetic_simulator.py (or the whole repo),
stage the updated file, and re-run tests to ensure the formatting job passes; no
logic changes are required—only formatting adjustments.

kirilklein and others added 2 commits April 3, 2026 13:10
Same pattern as the DOB/age fallback: if all patients in a shard
have no events, mean_recency is NaN. Fall back to 365 days.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

@coderabbitai coderabbitai Bot left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actionable comments posted: 1

♻️ Duplicate comments (1)
corebehrt/modules/simulation/semisynthetic_simulator.py (1)

224-241: ⚠️ Potential issue | 🟠 Major

Clarify or fix the pre-index cutoff semantics.

The docstring says "before each patient's index date" but line 230 uses <=, which includes events at the exact index timestamp. If treatment is assigned at index time, same-timestamp events could leak contemporaneous information into oracle features. Either change to strict < or update the docstring to clarify that "at or before" is intentional.

🔧 Option 1: Use strict pre-index cutoff
-        df = df[df[TIMESTAMP_COL] <= patient_index]
+        df = df[df[TIMESTAMP_COL] < patient_index]
📝 Option 2: Update docstring if inclusive is intentional
     def _filter_to_pre_index(
         self, df: pd.DataFrame, index_dates: pd.Series
     ) -> pd.DataFrame:
-        """Keep only events before each patient's index date, excluding special codes."""
+        """Keep only events at or before each patient's index date, excluding special codes."""
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 224 -
241, The _filter_to_pre_index method claims to keep events "before each
patient's index date" but currently uses a non-strict comparison
(df[TIMESTAMP_COL] <= patient_index); change this to a strict cutoff
(df[TIMESTAMP_COL] < patient_index) to avoid including events that occur exactly
at index time, and update the method docstring to state "Keep only events
strictly before each patient's index date, excluding special codes." Reference
symbols to edit: _filter_to_pre_index, TIMESTAMP_COL, patient_index, and the
docstring at the top of the method.
🧹 Nitpick comments (2)
corebehrt/modules/simulation/semisynthetic_simulator.py (2)

562-591: Metric names may be misleading given factual-outcome scoring.

The commit notes indicate scoring auc_if_all_treated and auc_if_all_control against the factual outcome was intentional (to avoid near-perfect scores). However, the metric names suggest they measure how well each probability predicts the respective counterfactual outcome. Consider renaming to clarify semantics, e.g., auc_p_treated_vs_factual and auc_p_control_vs_factual.

🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 562 -
591, The metric keys and variable names auc_treated/auc_control are misleading
because they compute ROC AUC of p_treated and p_control against the factual
outcome; rename the output dict keys and variables to make semantics explicit
(e.g., auc_p_treated_vs_factual and auc_p_control_vs_factual) and update where
they are assigned and used (references to auc_treated, auc_control, and the dict
keys "auc_if_all_treated"/"auc_if_all_control") as well as the results_data
entries and any downstream code that reads results_data or theoretical_aucs to
use the new names consistently while leaving the factual AUC stored as
auc_factual/theoretical_aucs[outcome_name].

70-97: Consider streaming stats computation for large cohorts.

The current approach materializes all per-shard feature DataFrames in memory before computing statistics. For large cohorts, this could spike memory usage. A bounded-memory approach would accumulate per-feature counts, sums, and squared sums per shard instead.

♻️ Bounded-memory alternative
     def compute_global_feature_stats(self, shard_loader):
         """Pass 1: compute global mean/std across all shards for standardization."""
-        all_features = []
+        n_rows = 0
+        sum_ = None
+        sumsq = None
         for shard, _ in shard_loader():
             pids, is_exposed, index_dates = self._extract_treatment_and_index_dates(
                 shard
             )
             if len(pids) == 0:
                 continue
             history_df = self._filter_to_pre_index(shard, index_dates)
             history_df, pids, is_exposed, index_dates = self._apply_min_num_codes(
                 history_df, pids, is_exposed, index_dates
             )
             if len(pids) == 0:
                 continue
             features_df = extract_oracle_features(
                 history_df, pids, index_dates, self.config.features
             )
-            all_features.append(features_df)
+            shard_sum = features_df.sum().astype(float)
+            shard_sumsq = (features_df**2).sum().astype(float)
+            sum_ = shard_sum if sum_ is None else sum_.add(shard_sum, fill_value=0.0)
+            sumsq = shard_sumsq if sumsq is None else sumsq.add(shard_sumsq, fill_value=0.0)
+            n_rows += len(features_df)
 
-        if not all_features:
+        if n_rows == 0:
             logger.warning("No patients found during global stats computation")
             return
 
-        combined = pd.concat(all_features)
-        self._global_means = combined.mean()
-        self._global_stds = combined.std(ddof=0).replace(0, 1).fillna(1)
-        logger.info("Computed global feature stats from %d patients", len(combined))
+        self._global_means = sum_ / n_rows
+        variance = (sumsq / n_rows) - self._global_means.pow(2)
+        self._global_stds = np.sqrt(variance.clip(lower=0)).replace(0, 1).fillna(1)
+        logger.info("Computed global feature stats from %d patients", n_rows)
🤖 Prompt for AI Agents
Verify each finding against the current code and only fix it if needed.

In `@corebehrt/modules/simulation/semisynthetic_simulator.py` around lines 70 -
97, compute_global_feature_stats currently appends all per-shard features to
all_features and materializes them, which can OOM for large cohorts; change it
to stream aggregations by iterating shard_loader and for each shard calling
_extract_treatment_and_index_dates, _filter_to_pre_index, _apply_min_num_codes
and extract_oracle_features (same helpers) but instead of appending features_df,
update running per-feature accumulators: counts (non-NaN per feature), sums, and
sums of squares; after the loop compute self._global_means = sums / counts and
self._global_stds = sqrt(sumsq / counts - mean**2) (handle numeric stability,
ddof=0) and replace zeros/NaNs with 1 as before; ensure you skip shards with no
pids and treat NaNs in features_df by ignoring them in the per-feature
count/sum/sumsq updates.
🤖 Prompt for all review comments with AI agents
Verify each finding against the current code and only fix it if needed.

Inline comments:
In `@corebehrt/modules/simulation/oracle_features.py`:
- Around line 1-11: This file fails the linter formatting check; run the
project's formatter (ruff format) on
corebehrt/modules/simulation/oracle_features.py (or the repo root) to fix
import/order/whitespace issues so the top-level docstring, imports (numpy,
pandas, constants, FeatureConfig) and logger =
logging.getLogger("oracle_features") line conform to the project's formatting
rules and then re-run CI.

---

Duplicate comments:
In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 224-241: The _filter_to_pre_index method claims to keep events
"before each patient's index date" but currently uses a non-strict comparison
(df[TIMESTAMP_COL] <= patient_index); change this to a strict cutoff
(df[TIMESTAMP_COL] < patient_index) to avoid including events that occur exactly
at index time, and update the method docstring to state "Keep only events
strictly before each patient's index date, excluding special codes." Reference
symbols to edit: _filter_to_pre_index, TIMESTAMP_COL, patient_index, and the
docstring at the top of the method.

---

Nitpick comments:
In `@corebehrt/modules/simulation/semisynthetic_simulator.py`:
- Around line 562-591: The metric keys and variable names
auc_treated/auc_control are misleading because they compute ROC AUC of p_treated
and p_control against the factual outcome; rename the output dict keys and
variables to make semantics explicit (e.g., auc_p_treated_vs_factual and
auc_p_control_vs_factual) and update where they are assigned and used
(references to auc_treated, auc_control, and the dict keys
"auc_if_all_treated"/"auc_if_all_control") as well as the results_data entries
and any downstream code that reads results_data or theoretical_aucs to use the
new names consistently while leaving the factual AUC stored as
auc_factual/theoretical_aucs[outcome_name].
- Around line 70-97: compute_global_feature_stats currently appends all
per-shard features to all_features and materializes them, which can OOM for
large cohorts; change it to stream aggregations by iterating shard_loader and
for each shard calling _extract_treatment_and_index_dates, _filter_to_pre_index,
_apply_min_num_codes and extract_oracle_features (same helpers) but instead of
appending features_df, update running per-feature accumulators: counts (non-NaN
per feature), sums, and sums of squares; after the loop compute
self._global_means = sums / counts and self._global_stds = sqrt(sumsq / counts -
mean**2) (handle numeric stability, ddof=0) and replace zeros/NaNs with 1 as
before; ensure you skip shards with no pids and treat NaNs in features_df by
ignoring them in the per-feature count/sum/sumsq updates.
🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

  • Push a commit to this branch (recommended)
  • Create a new PR with the fixes

ℹ️ Review info
⚙️ Run configuration

Configuration used: defaults

Review profile: CHILL

Plan: Pro

Run ID: 795cb1df-21bf-46f5-851f-b076b927376e

📥 Commits

Reviewing files that changed from the base of the PR and between 50dd4f5 and 351c22b.

📒 Files selected for processing (2)
  • corebehrt/modules/simulation/oracle_features.py
  • corebehrt/modules/simulation/semisynthetic_simulator.py

Comment thread corebehrt/modules/simulation/oracle_features.py
kirilklein and others added 11 commits June 26, 2026 16:34
Adds a study runner that mirrors the resampling study but uses the
semi-synthetic simulator: each outer run is an independent simulation
(own seed) followed by K inner reshuffle fits (finetune -> calibrate ->
estimate), reusing the existing causal mains unchanged.

- simulate_semisynthetic: optional config-gated patient sampling
  (SampledShardLoader), default off so existing behaviour is unchanged
- config_semisynthetic: optional `cohort` path for saving sampled pids
- run_study.py: lean N x K runner with templated base configs
- run_semisynthetic_study azure component + job CLI registration
- submit_runs.sh: one command submits N parallel outer-run jobs
- base_configs + job template + docs/study.md (smoke test + full study)

Smoke-tested locally: sampling + null/medium outcomes recover true
effects (NULL RD 0.000, MEDIUM delta=0.5 -> RD 0.063) and emit the
exposure/outcome/counterfactual files the downstream steps consume.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
A coefficient or interaction naming a feature the extractor doesn't
produce was silently ignored (`if name in features_df.columns`), so a
typo in the outcome config would silently zero-weight that term and
change the true data-generating process with no warning — dangerous for
a simulation study where the config defines ground truth.

- add ORACLE_FEATURE_NAMES (canonical feature set) in oracle_features
- validate all outcome coefficients/interactions/heterogeneous terms at
  simulator init, raising ValueError listing unknown names
- tests: drift guard (extractor output == ORACLE_FEATURE_NAMES) + raises
  on unknown coefficient/interaction names

Addresses CodeRabbit PR #173 finding (semisynthetic_simulator.py:343).

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
The gist badge update fails with 401 (GIST_SECRET expired/unauthorized),
which turned the Test/Docstring Coverage checks red. Every other step in
these jobs is already continue-on-error; the badge step was missing it.
Add continue-on-error to the badge step so a failed badge update no
longer fails the check (the badge itself still updates once the gist
token is refreshed).

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
The simulator required an assigned_index_date column in the MEDS shards,
which only the pre-processed example data has. Real MEDS is raw and has
no such column (KeyError on TRACE).

Index dates are an analysis artifact (cohort selection output), so add
paths.index_dates: when set, load per-patient index dates from a cohort
dir / index_dates.csv (subject_id, time) and use them; otherwise fall
back to the assigned_index_date column. Treatment is still taken from the
data (presence of the exposure code).

Validated locally against example data via a constructed index_dates.csv:
recovers the true effect (ATE 0.062 vs true 0.063), all outputs produced.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Reads every estimate_results.csv under a study dir (each already carries
the point estimate, CI, and appended true_effect), groups by
model x method x outcome, and reports bias, empirical SD, mean estimated
SE, SE-calibration (SD_emp/mean_SE), and 95% CI coverage.

- single estimate per group (Phase 1): bias + coverage meaningful,
  SD/calibration undefined (NaN)
- many estimates (bootstrap refits / outcome redraws): full table

Validated on synthetic estimate_results fixtures.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Reworks the semi-synthetic study per the design discussion:

- Run on a FIXED, shared pre-existing cohort (e.g. a diabetes cohort).
  Drop per-run select_cohort and per-run sampling; the cohort's
  index_dates.csv + cohort_config.yaml define the fixed design, treatment
  is the real exposure. Add `cohort` as a component/runner input.
- noise_scale=0 and three effect sizes (null/medium/large) -> a single
  fixed true effect per scenario.
- Inner refits (-k) are the SE mechanism: k=1 single fit (estimate's
  internal-bootstrap CI), k>1 each refit trains on a bootstrap resample
  (reshuffle + top-level bootstrap), summarizer combines the K.
- Add the CatBoost baseline path alongside BERT (method + baseline).
- Update component inputs, job template (cohort path), submit_runs.sh
  (drop sampling, default k=1, forward extra runner flags), and docs.

Config generation, outcome-file match, component import, bash syntax and
lint all validated locally; downstream steps need Azure (pretrain +
CausalEstimate) so the Phase-1 run is the first end-to-end check.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
On a real cohort there is no generic EXPOSURE code in raw MEDS, so the
simulator found 0 exposed and never wrote exposure.csv -> prepare then
failed its (bogus) exposure pre-check.

Treatment is observed and already recorded by select_cohort_full in
cohort/exposures.csv (same file prepare reads). So:
- simulator: load exposed pids from cohort/exposures.csv when index_dates
  points at a cohort dir; exposure_code is now only a fallback.
- prepare config: drop the exposures/exposure keys; prepare reads exposure
  from the cohort (load_cohort_data), so the sim-output pre-check was wrong.

Validated locally with a bogus exposure_code: exposure still resolves from
the cohort, exposure.csv is written, true effect recovered (ATE 0.064 vs 0.063).

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
- resample training folds with replacement per refit; validation folds stay
  fixed so every patient keeps one out-of-fold score
- preserve duplicate pids through dataset resampling for bert and baseline
- collect per-refit patient-bootstrap draws in the estimate step and
  aggregate the k refit estimates in the study summary
- fail fast when estimator n_bootstrap < 1: the previously pushed runner
  set 0 for k>1 and crashed only after all refits had trained, and no
  unit test exercises that wiring

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant