| Item | Value |
|---|---|
| Source file count | 62 (31 production + 31 test) |
| Models ported | 16 / 16 (incl. cross-section DID) |
| Tests | 241 / 241 on all 4 backends (native, wasm-gc, wasm, js) |
| Warnings | 0 (under moon test --deny-warn) |
| Python cross-checks | 16 / 16 PASS |
| License | Apache-2.0 |
$ moon test --deny-warn
Total tests: 241, passed: 241, failed: 0.
$ moon run cmd/main
=== MoonBit DML PLR (partialling out) ===
true theta_0 = 1
estimated theta = 0.9763281577675552
standard error = 0.08740490230937094
...
$ python validate_irm_with_python.py
======================================================================
PASS |mb - handrolled_nrep5| (theta) = 5.24e-02 < max(MODEL_TOL=0.1, 2.0*handrolled_n5_se) = 1.91e-01let data = @dml.DoubleMLData::new(x, y, d) // x : Matrix, y / d : Array[Double]
let fitted = @dml.DoubleMLPLR::new(data, n_folds=2, n_rep=1, seed=3141).fit()
let coef = fitted.coef() // Double
let se = fitted.se() // Double
let (lo, hi) = fitted.confint()| Driver | Model | DGP | True θ |
|---|---|---|---|
| cmd/main | DoubleMLPLR (and 5 others) | Simple partially linear, n=500, p=5 | 1.0 |
| cmd/datasets | DoubleMLPLR + DoubleMLIRM | Synthetic 401(k)-style, n=5000, 9 controls | 1.5 |
| cmd/did_binary | DoubleMLDIDBinary | 2-period panel DID, 400 units | 1.0 |
| cmd/did_cs | DoubleMLDIDCS | Staggered CS-DID, 4 cohorts × 4 periods | 1.0 |
| cmd/did_multi | DoubleMLDIDMulti | Top-level multi-period DID + aggregation | 1.0 |
| cmd/did_cross_section | DoubleMLDIDCrossSection | Sant'Anna-Zhao 2020 cross-section DID, 500 units | 1.0 |
$ moon run cmd/main
$ moon run cmd/datasets
$ moon run cmd/did_binary
$ moon run cmd/did_cs
$ moon run cmd/did_multi
$ moon run cmd/did_cross_section| Model | Score | Description |
|---|---|---|
| DoubleMLPLR | partialling-out | partially linear regression |
| DoubleMLIRM | ATE | interactive regression model |
| DoubleMLPLIV | partialling-out | partially linear IV regression |
| DoubleMLIIVM | LATE | interactive IV model |
| DoubleMLDID | observational | difference-in-differences |
| DoubleMLSSM | MAR | sample selection (missing-at-random) |
| DoubleMLAPO | APO | average potential outcome |
| DoubleMLAPOS | APOS | average potential outcome (share) |
| DoubleMLPQ | PQ | potential quantile |
| DoubleMLQTE | PQ | quantile treatment effect |
| DoubleMLLPQ | LPQ | local potential quantile (compliers) |
| DoubleMLCVAR | CVaR | conditional value-at-risk |
| DoubleMLRDD | (sharp / fuzzy) | regression discontinuity |
| DoubleMLBLP | BLP | best linear predictor |
| DoubleMLPolicyTree | (depth-N) | policy tree on weighted-variance-reduction gain |
moon test --target native --deny-warn # 241/241
moon test --target wasm-gc --deny-warn # 241/241
moon test --target wasm --deny-warn # 241/241
moon test --target js --deny-warn # 241/241$ for s in validate_*_with_python.py; do echo "=== $s ==="; python $s | tail -1; done
=== validate_blp_policy_with_python.py === BLP/PolicyTree reference checks passed
=== validate_bootstrap_with_python.py === Multipliers match: PASS
=== validate_did_with_python.py === PASS |mb - handrolled_nrep5| (theta) = 1.29e-02 ...
=== validate_did_binary_with_python.py === Reference: run `moon run cmd/did_binary` for the MoonBit output.
=== validate_did_cross_section_with_python.py === Cross-section DID reference: PASS
=== validate_did_cs_with_python.py === Reference: run `moon run cmd/did_cs` for the MoonBit output.
=== validate_gain_statistics_with_python.py === Gain statistics match: PASS
=== validate_iivm_with_python.py === PASS |mb - handrolled_nrep5| (theta) = 8.58e-03 ...
=== validate_irm_with_python.py === PASS |mb - handrolled_nrep5| (theta) = 5.24e-02 ...
=== validate_padjust_with_python.py === Romano-Wolf reference matches: PASS
=== validate_pava_with_python.py === PAVA cross-check passed
=== validate_pliv_with_python.py === PASS |mb - handrolled_nrep5| (theta) = 1.41e-01 ...
=== validate_quantile_with_python.py === reference checks passed
=== validate_rdd_with_python.py === RDD reference checks passed
=== validate_ssm_with_python.py === SSM reference checks passed
=== validate_with_python.py === Sanity check: true theta = 1.0 is inside every confidence interval.// Deterministic chacha8 RNG keyed by an integer seed.
// Replaces the 3-line `seed_to_bytes -> Bytes::from_array -> chacha8`
// boilerplate that used to live in every test file.
let rng = @dml.chacha8_rng(3141)
// Stratum labels for stratified K-fold partitioning.
// Pairs with `stratified_kfold` in `resampling.mbt` to balance
// (G, T) cells across folds when used inside `DoubleMLDID`.
let strata : Array[Int] = []
for i = 0; i < n; i = i + 1 {
strata = strata + [g_indicator[i] + 2 * t_indicator[i]]
}
// Propensity-score processing. Default clips to `[1e-2, 1 - 1e-2]`.
let psp = @dml.PSProcessor::new() // defaults
let psp = @dml.PSProcessor::new(config=@dml.PSProcessorConfig::new(clipping_threshold=0.05))
let out = psp.adjust_ps(ps_array, treatment_array)
// v0.14.0+: isotonic (PAVA) calibration of the propensity
// scores. The calibrated output is a step function from
// PAVA on `(ps, treatment)`, then clipped to
// `[clipping_threshold, 1 - clipping_threshold]`. Combine
// with `cv_calibration=true` for K-fold cross-validated
// predictions (matches upstream `cross_val_predict(cv=5)`).
let cfg = @dml.PSProcessorConfig::new(
calibration_method="isotonic", // v0.14.0+: PAVA fit
)
let psp_iso = @dml.PSProcessor::new(config=cfg)
let out_iso = psp_iso.adjust_ps(ps_array, treatment_array)
// 5-fold CV calibration with the deterministic kfold split.
let cfg_cv = @dml.PSProcessorConfig::new(
calibration_method="isotonic",
cv_calibration=true,
)
let psp_cv = @dml.PSProcessor::new(config=cfg_cv)
let out_cv = psp_cv.adjust_ps(ps_array, treatment_array)
// v0.15.0+: multiplier bootstrap for joint confidence
// intervals on `DoubleMLDIDMulti`. Draws `n_rep_boot` weight
// vectors from the chosen multiplier distribution ("normal"
// / "Bayes" / "wild") and computes per-cell t-statistics.
// Joint CIs use the empirical 95th percentile of the
// max-abs-t distribution as the critical value; pointwise
// CIs use 1.96. `joint=true` CIs are wider (more
// conservative).
let fitted = @dml.DoubleMLDIDMulti::new(data, n_folds=2, seed=3141).fit()
let booted = fitted.bootstrap(method_name="normal", n_rep_boot=500, seed=2024)
let ci_pw = booted.confint(joint=false) // Wald-style (1.96 * se)
let ci_joint = booted.confint(joint=true) // bootstrap critical value
// v0.16.0+: multiple-testing p-value adjustment.
// "romano-wolf" (default) requires the bootstrap; "holm",
// "bonferroni", "bh", "by" don't. Returns an Array[Double]
// of length n_combinations.
let pv_rw = booted.p_adjust(method_name="romano-wolf")
let pv_holm = booted.p_adjust(method_name="holm")
let pv_bonf = booted.p_adjust(method_name="bonferroni")
// v0.18.0+: FDR-controlling adjustments (Benjamini-Hochberg
// and Benjamini-Yekutieli). Both consume only the unadjusted
// p-values, so they don't require `bootstrap()`.
let pv_bh = fitted.p_adjust(method_name="bh")
let pv_by = fitted.p_adjust(method_name="by")
// v0.24.0+: two-stage FDR (more powerful when some
// hypotheses are non-null). "tsbh" / "fdr_tsbh" and
// "tsby" / "fdr_tsbky" are statsmodels-compatible
// aliases.
let pv_tsbh = fitted.p_adjust(method_name="tsbh")
let pv_tsby = fitted.p_adjust(method_name="tsby")
// v0.17.0+: gain statistics for sensitivity parameter
// benchmarks. Pass two `GainStatsSource` (one for the
// "long" model with all confounders, one for the
// "short" model with benchmark confounders excluded);
// returns `cf_y / cf_d / rho / delta_theta` per
// coefficient. Use as the upper bound on the
// sensitivity parameters in `sensitivity_analysis`.
let src_long = @dml.GainStatsSource::new(
var_y_residuals_long, nu2_long, all_coef_long, n_rep, var_y,
)
let src_short = @dml.GainStatsSource::new(
var_y_residuals_short, nu2_short, all_coef_short, n_rep, var_y,
)
let gs = @dml.gain_statistics(src_long, src_short)
// v0.19.0+: `from_blp(blp)` auto-populates a
// `GainStatsSource` from a fitted `DoubleMLBLP`.
// `var_y_residuals = rss / n_obs`,
// `nu2[k] = var_y_residuals / (n * se[k]^2)`,
// `all_coef = blp.coef()`, `var_y = blp.var_y()`.
let src = @dml.GainStatsSource::from_blp(blp_fitted)pub suberror BootstrapMethodError {
UnknownMethod(String)
}pub suberror BracketSignError {
UpperSignFailed
}pub suberror CalibrationFittingError {
IncompleteCVPartition
}pub suberror ClusterDataError {
MissingUnit(Int)
}pub suberror DIDDataError {
NonBinaryTreatment(Int)
}pub suberror EmptyArrayErrorpub suberror InvalidCalibrationError {
UnknownMethod(String)
}pub suberror PSConfigError {
InconsistentCVCalibration
}pub suberror VarEstClusterError {
JTooSmall(Double, Double, Int)
}type DidCsDatatype DidMultiDatapub struct DoubleMLAPO {
data : DoubleMLData
treatment_level : Double
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
g_hat : Array[Double]
m_hat : Array[Double]
coef : Double
se : Double
fitted : Bool
} derive(Debug)fn DoubleMLAPO::new(data : DoubleMLData, treatment_level? : Double, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLAPOpub struct DoubleMLAPOS {
data : DoubleMLData
treatment_levels : Array[Double]
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
coefs : Array[Double]
ses : Array[Double]
fitted : Bool
} derive(Debug)fn DoubleMLAPOS::causal_contrast(self : DoubleMLAPOS, reference_levels : Array[Double]) -> Array[Array[Double]]fn DoubleMLAPOS::new(data : DoubleMLData, treatment_levels : Array[Double], n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLAPOSfn DoubleMLBLP::confint_joint(self : DoubleMLBLP, contrast : Matrix, level? : Double) -> Array[(Double, Double)] Y = expit(D * theta_0 + r_0(X)), Y in {0, 1}pub struct DoubleMLCVAR {
data : DoubleMLData
treatment : Double
quantile : Double
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
normalize_ipw : Bool
g_hat : Array[Double]
m_hat : Array[Double]
coef : Double
se : Double
fitted : Bool
} derive(Debug)fn DoubleMLCVAR::new(data : DoubleMLData, treatment? : Double, quantile? : Double, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double, normalize_ipw? : Bool) -> DoubleMLCVARpub struct DoubleMLDID {
data : DoubleMLDIDData
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
ps_processor : PSProcessor
score : String
in_sample_normalization : Bool
strata : Array[Int]
g0_hat : Array[Double]
g1_hat : Array[Double]
m_hat : Array[Double]
coef : Double
se : Double
psi_a : Array[Double]
psi_b : Array[Double]
fitted : Bool
} derive(Debug) Y_post = g_0(0, X) + D * theta + U_post, E[U | D, X] = 0
Y_pre = g_0(0, X) + U_pre, E[U_pre | X] = 0 dY = Y_post - Y_pre = g_1(1, X) - g_0(0, X) + D * theta + (U_post - U_pre) g0(X) = E[Y | D = 0, X] (trained on D = 0)
g1(X) = E[Y | D = 1, X] (trained on D = 1)
m(X) = P(D = 1 | X) (trained on all obs; observational
only, clipped to [eps, 1 - eps]) resid_d0 = Y - g0
p_hat = mean(D)
weight_psi_a = D / p_hat
weight_resid = (D - m) / (p_hat * (1 - m))
psi_b = (D - m)/(p_hat (1 - m)) * (Y - g0) [the g1 term cancels in ATT]
psi_a = -D / p_hat
psi(theta) = theta * psi_a + psi_b theta_hat = -mean(psi_b) / mean(psi_a)
J = mean(psi_a)
gamma = mean(psi(theta_hat)^2)
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2).fn DoubleMLDID::fit(self : DoubleMLDID, ml_g? : LinearRegression, ml_m? : LinearRegression) -> DoubleMLDIDfn DoubleMLDID::new(data : DoubleMLDIDData, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double, ps_processor? : PSProcessor, score? : String, in_sample_normalization? : Bool, strata? : Array[Int]) -> DoubleMLDIDpub struct DoubleMLDIDBinary {
data : DoubleMLDIDBinaryData
g_value : Int
t_value_pre : Int
t_value_eval : Int
control_group : String
anticipation_periods : Int
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
ps_processor : PSProcessor
score : String
in_sample_normalization : Bool
eval_idx : Array[Int]
inner : DoubleMLDID
fitted : Bool
} derive(Debug)fn DoubleMLDIDBinary::fit(self : DoubleMLDIDBinary, ml_g? : LinearRegression, ml_m? : LinearRegression) -> DoubleMLDIDBinaryfn DoubleMLDIDBinary::new(data : DoubleMLDIDBinaryData, g_value : Int, t_value_pre : Int, t_value_eval : Int, control_group? : String, anticipation_periods? : Int, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double, ps_processor? : PSProcessor, score? : String, in_sample_normalization? : Bool) -> DoubleMLDIDBinaryfn DoubleMLDIDBinaryData::new(x : Matrix, y : Array[Double], d : Array[Double], t : Array[Int], g : Array[Int], id : Array[Int]) -> DoubleMLDIDBinaryDatapub struct DoubleMLDIDCS {
data : DoubleMLDIDCSData
control_group : String
anticipation_periods : Int
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
ps_processor : PSProcessor
in_sample_normalization : Bool
coef_matrix : Array[Double]
se_matrix : Array[Double]
psi_matrix : Array[Double]
n_groups : Int
n_periods : Int
fitted : Bool
} derive(Debug)fn DoubleMLDIDCS::fit(self : DoubleMLDIDCS, ml_g? : LinearRegression, ml_m? : LinearRegression) -> DoubleMLDIDCSfn DoubleMLDIDCS::new(data : DoubleMLDIDCSData, control_group? : String, anticipation_periods? : Int, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double, ps_processor? : PSProcessor, in_sample_normalization? : Bool) -> DoubleMLDIDCSpub struct DoubleMLDIDCSBinary {
data : DoubleMLDIDCSData
g_value : Int
t_value_pre : Int
t_value_eval : Int
control_group : String
anticipation_periods : Int
n_folds : Int
n_rep : Int
seed : Int
propensity_clip : Double
ps_processor : PSProcessor
score : String
in_sample_normalization : Bool
coef : Double
se : Double
psi_a : Array[Double]
psi_b : Array[Double]
g_d0_t0 : Array[Double]
g_d0_t1 : Array[Double]
g_d1_t0 : Array[Double]
g_d1_t1 : Array[Double]
m_hat : Array[Double]
n_obs_subset : Int
n_g_subset : Int
n_c_subset : Int
fitted : Bool
} derive(Debug)fn DoubleMLDIDCSBinary::fit(self : DoubleMLDIDCSBinary, ml_g? : LinearRegression, ml_m? : LinearRegression) -> DoubleMLDIDCSBinaryfn DoubleMLDIDCSBinary::new(data : DoubleMLDIDCSData, g_value : Int, t_value_pre : Int, t_value_eval : Int, control_group? : String, anticipation_periods? : Int, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double, ps_processor? : PSProcessor, score? : String, in_sample_normalization? : Bool) -> DoubleMLDIDCSBinaryfn DoubleMLDIDCSData::new(x : Matrix, y : Array[Double], d : Array[Double], t : Array[Int], id : Array[Int], g : Array[Int]) -> DoubleMLDIDCSDatapub struct DoubleMLDIDCrossSection {
data : DoubleMLDIDCrossSectionData
n_folds : Int
n_rep : Int
seed : Int
score : String
in_sample_normalization : Bool
propensity_clip : Double
ps_processor : PSProcessor
coef : Double
se : Double
psi_a : Array[Double]
psi_b : Array[Double]
g_d0_t0 : Array[Double]
g_d0_t1 : Array[Double]
g_d1_t0 : Array[Double]
g_d1_t1 : Array[Double]
m_hat : Array[Double]
fitted : Bool
boot_t_stat : Array[Double]
boot_method : String
n_rep_boot : Int
boot_seed : Int
} derive(Debug)fn DoubleMLDIDCrossSection::bootstrap(self : DoubleMLDIDCrossSection, method_name? : String, n_rep_boot? : Int, seed? : Int) -> DoubleMLDIDCrossSectionfn DoubleMLDIDCrossSection::confint(self : DoubleMLDIDCrossSection, joint? : Bool, level? : Double) -> (Double, Double)fn DoubleMLDIDCrossSection::new(data : DoubleMLDIDCrossSectionData, n_folds? : Int, n_rep? : Int, seed? : Int, score? : String, in_sample_normalization? : Bool, propensity_clip? : Double, ps_processor? : PSProcessor) -> DoubleMLDIDCrossSectionfn DoubleMLDIDCrossSectionData::new(x : Matrix, y : Array[Double], d : Array[Double], t : Array[Int], name? : String) -> DoubleMLDIDCrossSectionDatafn DoubleMLDIDData::new(x : Matrix, y : Array[Double], d : Array[Double]) -> DoubleMLDIDData raise DIDDataErrorpub struct DoubleMLDIDMulti {
data : DoubleMLDIDCSData
gt_combinations : Array[(Int, Int, Int)]
control_group : String
anticipation_periods : Int
n_folds : Int
n_rep : Int
seed : Int
ps_processor : PSProcessor
in_sample_normalization : Bool
group_sizes : Array[Int]
inner : DoubleMLDIDCS
boot_t_stat : Array[Double]
boot_method : String
n_rep_boot : Int
boot_seed : Int
fitted : Bool
} derive(Debug)fn DoubleMLDIDMulti::bootstrap(self : DoubleMLDIDMulti, method_name? : String, n_rep_boot? : Int, seed? : Int) -> DoubleMLDIDMultifn DoubleMLDIDMulti::confint(self : DoubleMLDIDMulti, joint? : Bool, level? : Double) -> Array[(Double, Double)]fn DoubleMLDIDMulti::new(data : DoubleMLDIDCSData, gt_combinations? : Array[(Int, Int, Int)], gt_combinations_keyword? : String, control_group? : String, anticipation_periods? : Int, n_folds? : Int, n_rep? : Int, seed? : Int, ps_processor? : PSProcessor, in_sample_normalization? : Bool) -> DoubleMLDIDMultifn DoubleMLData::new(x : Matrix, y : Array[Double], d : Array[Double], cluster_vars? : Array[Int]) -> DoubleMLData Y = theta * D + g_0(D, X) + U, E[U | D, X] = 0
D = m_0(X, Z) + V, E[V | X, Z] = 0 g0(X) = E[Y | Z = 0, X] (trained only on Z = 0)
g1(X) = E[Y | Z = 1, X] (trained only on Z = 1)
m(X) = E[Z | X] (trained on all obs, then
clipped to [eps, 1 - eps])
r0(X) = E[D | Z = 0, X] (trained only on Z = 0)
r1(X) = E[D | Z = 1, X] (trained only on Z = 1) u_hat0 = Y - g0, u_hat1 = Y - g1
w_hat0 = D - r0, w_hat1 = D - r1 psi_b = (g1 - g0) + Z u_hat1 / m - (1 - Z) u_hat0 / (1 - m)
psi_a = -(r1 - r0) - Z w_hat1 / m + (1 - Z) w_hat0 / (1 - m)
psi(theta) = theta * psi_a + psi_b theta_hat = -mean(psi_b) / mean(psi_a)
J = mean(psi_a)
gamma = mean(psi(theta_hat)^2)
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2).fn DoubleMLIIVM::fit(self : DoubleMLIIVM, ml_g? : LinearRegression, ml_m? : LinearRegression, ml_r? : LinearRegression, max_attempts? : Int) -> DoubleMLIIVMfn DoubleMLIIVM::new(data : DoubleMLIIVMData, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLIIVMfn DoubleMLIIVMData::new(x : Matrix, y : Array[Double], d : Array[Double], z : Array[Double], cluster_vars? : Array[Int]) -> DoubleMLIIVMData Y = g_0(D, X) + U, E[U | D, X] = 0
D = m_0(X) + V, E[V | X] = 0 g0(X) = E[Y | D=0, X]
g1(X) = E[Y | D=1, X]
m(X) = P(D=1 | X) (the propensity score)
u0 = Y - g0(X)
u1 = Y - g1(X)
psi_b = (g1 - g0) + (D u1 / m - (1 - D) u0 / (1 - m))
psi_a = -1
psi(theta) = theta * psi_a + psi_b theta_hat = -mean(psi_b) / mean(psi_a) = mean(psi_b) J = mean(psi_a) = -1
gamma = mean(psi(theta_hat)^2)
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2)fn DoubleMLIRM::fit(self : DoubleMLIRM, ml_g? : LinearRegression, ml_m? : LinearRegression, max_attempts? : Int) -> DoubleMLIRMfn DoubleMLIRM::new(data : DoubleMLData, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLIRMpub struct DoubleMLLPLR {
data : DoubleMLBinaryData
score : String
n_folds : Int
n_folds_inner : Int
n_rep : Int
seed : Int
coef_ : Double
se_ : Double
r_hat : Array[Double]
m_hat : Array[Double]
a_hat : Array[Double]
fitted : Bool
} derive(Debug)fn DoubleMLLPLR::new(data : DoubleMLBinaryData, score? : String, n_folds? : Int, n_folds_inner? : Int, n_rep? : Int, seed? : Int) -> DoubleMLLPLRpub struct DoubleMLLPQ {
data : DoubleMLLPQData
treatment : Double
quantile : Double
n_folds : Int
seed : Int
propensity_clip : Double
coef : Double
se : Double
fitted : Bool
predictions_g0 : Array[Double]
predictions_g1 : Array[Double]
predictions_m : Array[Double]
} derive(Debug)fn DoubleMLLPQ::new(data : DoubleMLLPQData, treatment? : Double, quantile? : Double, n_folds? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLLPQfn DoubleMLLPQData::new(x : Matrix, y : Array[Double], d : Array[Double], z : Array[Double]) -> DoubleMLLPQDatapub struct DoubleMLPLIV {
data : DoubleMLPLIVData
n_folds : Int
n_rep : Int
seed : Int
l_hat : Array[Double]
r_hat : Array[Double]
m_hat : Array[Double]
coef : Double
se : Double
fitted : Bool
} derive(Debug) Y = D * theta_0 + g_0(X) + zeta, E[zeta | D, X] = 0
D = m_0(X) + V, E[V | X] = 0
Z = ell_0(X) + xi, E[xi | X] = 0,
Cov(Z, V) != 0 (relevance) l_hat = E_hat[Y | X]
r_hat = E_hat[D | X]
m_hat = E_hat[Z | X]
u_hat = Y - l_hat
w_hat = D - r_hat
v_hat = Z - m_hat
psi_a = -w_hat * v_hat
psi_b = v_hat * u_hat
psi(theta) = theta * psi_a + psi_b theta_hat = -mean(psi_b) / mean(psi_a)
= mean(v_hat * u_hat) / mean(w_hat * v_hat)
J = mean(psi_a)
gamma = mean(psi(theta_hat)^2)
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2).fn DoubleMLPLIV::fit(self : DoubleMLPLIV, learner? : LinearRegression, max_attempts? : Int) -> DoubleMLPLIVfn DoubleMLPLIV::new(data : DoubleMLPLIVData, n_folds? : Int, n_rep? : Int, seed? : Int) -> DoubleMLPLIVfn DoubleMLPLIVData::new(x : Matrix, y : Array[Double], d : Array[Double], z : Array[Double], cluster_vars? : Array[Int]) -> DoubleMLPLIVDatapub struct DoubleMLPLPR {
panel : DoubleMLPanelData
approach : String
score : String
n_folds : Int
n_rep : Int
seed : Int
coef_ : Double
se_ : Double
l_hat : Array[Double]
m_hat : Array[Double]
fitted : Bool
} derive(Debug)fn DoubleMLPLPR::fit(self : DoubleMLPLPR, learner? : LinearRegression, max_attempts? : Int) -> DoubleMLPLPRfn DoubleMLPLPR::new(panel : DoubleMLPanelData, approach? : String, score? : String, n_folds? : Int, n_rep? : Int, seed? : Int) -> DoubleMLPLPRpub struct DoubleMLPLR {
data : DoubleMLData
n_folds : Int
n_rep : Int
seed : Int
l_hat : Array[Double]
m_hat : Array[Double]
coef : Double
se : Double
fitted : Bool
} derive(Debug) Y = D * theta_0 + g_0(X) + zeta, E[zeta | D, X] = 0
D = m_0(X) + V, E[V | X] = 0 psi_a(theta) = -(D - m_hat)^2,
psi_b(theta) = (D - m_hat) * (Y - l_hat),
psi(theta) = theta * psi_a + psi_b theta_hat = -mean(psi_b) / mean(psi_a)
= mean((D - m_hat)(Y - l_hat)) / mean((D - m_hat)^2). J = mean(psi_a) # expected derivative of psi w.r.t. theta
gamma = mean(psi(theta_hat)^2)
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2).fn DoubleMLPLR::fit(self : DoubleMLPLR, learner? : LinearRegression, max_attempts? : Int) -> DoubleMLPLRpub struct DoubleMLPQ {
data : DoubleMLData
treatment : Double
quantile : Double
n_folds : Int
seed : Int
propensity_clip : Double
coef : Double
se : Double
fitted : Bool
} derive(Debug)fn DoubleMLPQ::new(data : DoubleMLData, treatment? : Double, quantile? : Double, n_folds? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLPQ Y_it = D_it * theta_0 + g_0(X_it) + alpha_i + zeta_itfn DoubleMLPanelData::new(x : Matrix, y : Array[Double], d : Array[Double], t : Array[Int], id : Array[Int]) -> DoubleMLPanelDatapub struct DoubleMLPolicyTree {
features : Matrix
orth_signal : Array[Double]
depth : Int
root : PolicyTreeNode
split_feature : Int
split_value : Double
left_treatment : Int
right_treatment : Int
fitted : Bool
} derive(Debug)fn DoubleMLPolicyTree::new(features : Matrix, orth_signal : Array[Double], depth? : Int) -> DoubleMLPolicyTreepub struct DoubleMLQTE {
data : DoubleMLData
quantiles : Array[Double]
n_folds : Int
seed : Int
propensity_clip : Double
coefs : Array[Double]
ses : Array[Double]
} derive(Debug)fn DoubleMLQTE::new(data : DoubleMLData, quantiles? : Array[Double], n_folds? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLQTEpub struct DoubleMLRDD {
data : DoubleMLRDDData
cutoff : Double
bandwidth : Double
fuzzy : Bool
cov_type : String
coef : Double
se : Double
n_local : Int
fitted : Bool
} derive(Debug)fn DoubleMLRDD::new(data : DoubleMLRDDData, cutoff? : Double, bandwidth? : Double, fuzzy? : Bool, cov_type? : String) -> DoubleMLRDDfn DoubleMLRDDData::new(x : Matrix, y : Array[Double], d : Array[Double], score : Array[Double]) -> DoubleMLRDDData Y = theta * D + X @ beta * D + U, E[U | D, X] = 0
S = 1{D + gamma Z + X @ beta + V > 0}, E[V | X, D] = 0 g_d1(X) = E[Y | D = 1, S = 1, X] (trained on D=1 ∧ S=1,
features = X only —
Bug #1 fix: previously
`pi_hat` was appended as
an extra feature, which
leaked the test-fold pi
into the training fold)
g_d0(X) = E[Y | D = 0, S = 1, X] (trained on D=0 ∧ S=1,
features = X only)
m(X) = P(D = 1 | X) (trained on all obs,
clipped to [eps, 1 - eps])
pi(X, D) = P(S = 1 | D, X) (trained on (X, D),
clipped to [eps, 1 - eps]) psi_a = -1
psi_b1 = (D == 1) * S * (Y - g_d1) / (m * pi) + g_d1
psi_b0 = (D == 0) * S * (Y - g_d0) / ((1 - m) * pi) + g_d0
psi_b = psi_b1 - psi_b0 theta_hat = -mean(psi_b) / mean(psi_a) = mean(psi_b)
J = mean(psi_a) = -1
gamma = mean(psi(theta_hat)^2)
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2).fn DoubleMLSSM::fit(self : DoubleMLSSM, ml_g? : LinearRegression, ml_m? : LinearRegression, ml_pi? : LinearRegression) -> DoubleMLSSMfn DoubleMLSSM::new(data : DoubleMLSSMData, n_folds? : Int, n_rep? : Int, seed? : Int, propensity_clip? : Double) -> DoubleMLSSMfn DoubleMLSSMData::new(x : Matrix, y : Array[Double], d : Array[Double], s : Array[Double]) -> DoubleMLSSMDatafn GainStatsSource::from_blp_cv_repeated(blp : DoubleMLBLP, n_folds? : Int, n_repeats? : Int, seed? : Int) -> GainStatsSourcefn GainStatsSource::new(var_y_residuals : Array[Double], nu2 : Array[Double], all_coef : Array[Double], n_rep : Int, var_y : Double) -> GainStatsSourcetype IivmDatatype IrmConfoundedDatatype IrmDatatype IrmDiscreteDatatype IrmHeterogeneousDatafn LinearRegression::fit(self : LinearRegression, x : Matrix, y : Array[Double]) -> LinearRegressionfn LinearRegression::fit_weighted(self : LinearRegression, x : Matrix, y : Array[Double], w : Array[Double]) -> LinearRegressionfn LinearRegression::sandwich_se(self : LinearRegression, x : Matrix, y : Array[Double]) -> Array[Double]fn LinearRegression::sandwich_se_weighted(self : LinearRegression, x : Matrix, y : Array[Double], w : Array[Double]) -> Array[Double] cov(beta_hat) = (X'WX)^{-1} (X' diag(w · e^2) X) (X'WX)^{-1}fn LogisticRegression::fit(self : LogisticRegression, x : Matrix, y : Array[Double], max_iter? : Int, tol? : Double, ridge? : Double) -> LogisticRegressionfn LogisticRegression::predict_class(self : LogisticRegression, x : Matrix, threshold? : Double) -> Array[Double]type LplrDatafn PSProcessor::adjust_ps(self : PSProcessor, ps : Array[Double], treatment : Array[Double], cv? : Array[(Array[Int], Array[Int])]?) -> Array[Double]pub(all) struct PSProcessorConfig {
clipping_threshold : Double
extreme_threshold : Double
calibration_method : String
cv_calibration : Bool
} derive(Debug)fn PSProcessorConfig::new(clipping_threshold? : Double, extreme_threshold? : Double, calibration_method? : String, cv_calibration? : Bool) -> PSProcessorConfig raise PSConfigErrortype PlivClusterDatatype PlivDatatype PlprDatatype PlrCcddhnr2018type PlrConfoundedDatapub enum PolicyTreeNode {
Leaf(Int)
Split(Int, Double, PolicyTreeNode, PolicyTreeNode)
} derive(Debug)type RddSimpleDatatype SsmData theta_hat = median(theta_1, ..., theta_R)
ub_r = theta_r + 1.96 * se_r for each rep r
ub_hat = median(ub_1, ..., ub_R)
se_hat = (ub_hat - theta_hat) / 1.96 theta_hat = median([c]) = c
ub_hat = median([c + 1.96 * s]) = c + 1.96 * s
se_hat = (c + 1.96 * s - c) / 1.96 = s#callsite(autofill(loc))
fn check(condition : Bool, loc~ : SourceLoc) -> Unit raise PreconditionErrorfn draw_bootstrap_weights(method_name : String, n_rep_boot : Int, n_obs : Int, seed : Int) -> Array[Double] raise BootstrapMethodError theta = - sum_k w_k * sum_{i in k} psi_b
/ sum_k w_k * sum_{i in k} psi_a, w_k = 1/|I_k|.fn expit(x : Double) -> Doublefn g_cross_fit_calls() -> Int f_hat(x) = (1 / (n * h * sqrt(2*pi))) * sum_i exp(-(x - y_i)^2 / (2 * h^2)) f_hat(x) = (1 / (h * sqrt(2*pi))) * sum_i w[i] * exp(-(x - y[i])^2 / (2*h^2)) d/dtheta mean(psi_ipw) = (1/n) * sum_i w_i * delta(y_i - theta)
≈ (1/n) * f_hat_weighted(theta) s = 0.0
for x in arr: s = s + x sum = 0.0
c = 0.0 // compensation for low-order bits lost in the next add
for x in arr:
y = x - c // align x with the running sum's precision
t = sum + y // primary add
c = (t - sum) - y // low-order bits of `t` that did not fit
sum = t
return sumfn logit(p : Double, eps? : Double) -> Doublefn lpq_score_ipw(data : DoubleMLLPQData, treated : Array[Double], m : Array[Double], comp : Double, theta : Double, q : Double, sign : Double) -> Array[Double]fn make_confounded_plr_data(n_obs : Int, dim_x : Int, theta : Double, seed : Int) -> PlrConfoundedDatafn make_irm_confounded_data(n_obs : Int, dim_x : Int, theta : Double, seed : Int) -> IrmConfoundedDatafn make_irm_discrete_treatments(n_obs : Int, dim_x : Int, theta : Double, seed : Int) -> IrmDiscreteDatafn make_pliv_multiway_cluster(n_obs : Int, dim_x : Int, theta : Double, seed : Int) -> PlivClusterDatafn norm_cdf(x : Double) -> Doublefn norm_ppf(p : Double) -> Doublefn norm_sf(x : Double) -> Double#callsite(autofill(loc))
fn require(condition : Bool, loc~ : SourceLoc) -> Unit raise PreconditionErrorfn reset_g_cross_fit_count() -> Unitfn solve_pq(data : DoubleMLData, treatment : Double, q : Double, n_folds : Int, seed : Int, clip : Double) -> (Double, Array[Double], Double) raise BracketSignError theta_hat = -mean(psi_b) / mean(psi_a)
J = mean(psi_a)
gamma = mean(psi(theta_hat)^2) where psi(theta) = theta * psi_a + psi_b
sigma2 = gamma / (J^2 * n)
se = sqrt(sigma2) gamma += S_g^2 / |I_k|, S_g = sum of psi over unit g,
J += (sum of psi_deriv over the fold's rows) / |I_k|,
sigma2 = (gamma / npc) / (N_units * (J / npc)^2).Install
Download zip