antecedent_validate/
bootstrap_refute.rs1#![allow(clippy::cast_possible_truncation, clippy::cast_precision_loss, clippy::cast_sign_loss)]
14
15use std::sync::Arc;
16
17use antecedent_core::ExecutionContext;
18use antecedent_data::TableView;
19use antecedent_estimate::{EstimationWorkspace, LinearAdjustmentAte};
20use antecedent_kernels::unbiased_index;
21
22use crate::common::{
23 RefutationProblem, RefutationReport, complete_case_rows, linear_estimator_no_bootstrap,
24 refit_effect, with_resampled_rows,
25};
26use crate::error::ValidationError;
27
28#[derive(Clone, Debug)]
34pub struct BootstrapRefute {
35 pub replicates: u32,
37 pub ci_level: f64,
39 pub estimator: LinearAdjustmentAte,
41}
42
43impl Default for BootstrapRefute {
44 fn default() -> Self {
45 Self::new()
46 }
47}
48
49impl BootstrapRefute {
50 #[must_use]
52 pub fn new() -> Self {
53 Self { replicates: 200, ci_level: 0.95, estimator: linear_estimator_no_bootstrap() }
54 }
55
56 pub fn refute(
62 &self,
63 problem: &RefutationProblem<'_>,
64 workspace: &mut EstimationWorkspace,
65 ctx: &ExecutionContext,
66 ) -> Result<RefutationReport, ValidationError> {
67 if self.replicates < 2 {
68 return Err(ValidationError::NotApplicable {
69 message: "bootstrap CI coverage requires replicates >= 2",
70 });
71 }
72 if !(self.ci_level > 0.0 && self.ci_level < 1.0) {
73 return Err(ValidationError::NotApplicable {
74 message: "bootstrap CI coverage requires ci_level in (0, 1)",
75 });
76 }
77 let n = problem.data.row_count();
78 let mut resample_ids = vec![problem.treatment(), problem.outcome()];
79 if problem.temporal.is_none() {
81 resample_ids.extend_from_slice(&problem.estimand.adjustment_set);
82 } else {
83 for v in problem.data.schema().variables() {
85 let id = v.id;
86 if id != problem.treatment() && id != problem.outcome() {
87 resample_ids.push(id);
88 }
89 }
90 }
91 let (keep, valid) = complete_case_rows(problem.data, &resample_ids)?;
94 if valid.len() < 2 {
95 return Err(ValidationError::NotApplicable {
96 message: "bootstrap CI coverage requires at least 2 complete-case rows",
97 });
98 }
99 let mut rng = ctx.rng.stream(0xA7E0_0009_0000_u64);
100 let mut row_idx = vec![0usize; n];
101 let mut ates = Vec::with_capacity(self.replicates as usize);
102 for _ in 0..self.replicates {
103 for slot in &mut row_idx {
104 *slot = valid[unbiased_index(&mut rng, valid.len())];
105 }
106 let data = with_resampled_rows(problem.data, &resample_ids, &row_idx, &keep)?;
107 let est = refit_effect(
108 problem,
109 &data,
110 problem.estimand,
111 &[],
112 &self.estimator,
113 workspace,
114 ctx,
115 )?;
116 ates.push(est.ate);
117 }
118 ates.sort_by(|a, b| a.partial_cmp(b).unwrap_or(std::cmp::Ordering::Equal));
119 let m = ates.len();
120 let lo_frac = (1.0 - self.ci_level) / 2.0;
121 let hi_frac = 1.0 - lo_frac;
122 let lo_idx = ((lo_frac * (m - 1) as f64).round() as usize).min(m - 1);
123 let hi_idx = ((hi_frac * (m - 1) as f64).round() as usize).min(m - 1);
124 let lo = ates[lo_idx];
125 let hi = ates[hi_idx];
126 let mean_ate = ates.iter().sum::<f64>() / m as f64;
127 let width = hi - lo;
128 let passed = problem.original.ate >= lo && problem.original.ate <= hi;
129 Ok(RefutationReport {
130 refuter: Arc::from("bootstrap.ci_coverage"),
131 original_ate: problem.original.ate,
132 refuted_ate: mean_ate,
133 comparison: width,
134 informative: true,
135 passed,
136 failure_condition: if passed {
137 None
138 } else {
139 Some(Arc::from(format!(
140 "original ATE {} outside {}% bootstrap CI [{lo}, {hi}] \
141 (coverage check of the point estimate, not a placebo falsification)",
142 problem.original.ate,
143 self.ci_level * 100.0
144 )))
145 },
146 replicates: self.replicates,
147 })
148 }
149}