Simulation Diagnostics for Doubly Robust Estimators

Why Simulation Diagnostics Need a Design

AIPW and TMLE have useful robustness properties, but those properties do not guarantee identical finite-sample performance in every data set. Estimates can be variable when treatment probabilities are extreme, confidence intervals can miss their target more often than intended, and an estimator can fail to return a usable result. A simulation study makes these features visible under controlled conditions, because the data-generating process supplies a known causal target against which estimates can be checked.

The goal is not to declare one method universally superior. It is to learn how performance changes when one or both nuisance models are misspecified, when overlap weakens, and when the sample size changes. A credible comparison sets these conditions deliberately, uses the same simulated data for both methods, and reports failures as well as successful estimates.

Define the Target and Data-Generating Process

For a simple comparison, let \(A\) be a binary treatment, \(Y\) a binary outcome, and \(X\) a measured covariate. Write \(Y(1)\) and \(Y(0)\) for the potential outcomes under treatment and control. The average treatment effect (ATE) is \(\psi_0=E[Y(1)-Y(0)]\). A simulation designed to study AIPW and TMLE should state this target before generating data, then choose a process with a known \(\psi_0\).

For example, let \(X\sim\mathrm{Bernoulli}(0.5)\), and generate treatment according to \(P(A=1\mid X=0)=0.30\) and \(P(A=1\mid X=1)=0.70\). Generate potential outcomes with \(P(Y(0)=1\mid X)=0.15+0.10X\) and \(P(Y(1)=1\mid X)=0.35+0.10X\). Conditional on \(X\), treatment assignment is independent of the potential outcomes. The difference in the two conditional outcome probabilities is \(0.20\) at both values of \(X\), so the marginal ATE is also \(0.20\). Treatment probabilities lie between zero and one, and both treatment groups occur at each value of \(X\); thus positivity holds in this process.

To investigate overlap, keep the outcome process fixed but change the treatment probabilities, for instance to \(0.05\) and \(0.95\). Positivity still holds mathematically, but some covariate-treatment combinations are rare and inverse-probability weights can be large. To investigate model correctness, fit nuisance models that either represent the generating relationships or deliberately omit relevant structure. With binary \(X\), a treatment model that includes \(X\) can represent the specified propensity scores; a model omitting \(X\) generally cannot when the probabilities differ by \(X\). The outcome model should include the treatment and covariate structure needed to represent the conditional outcome means.

Choose Metrics Before Running the Simulation

For estimator \(j\), let \(\hat\psi_{jr}\) be its estimate in replicate \(r\), and let \(R_s\) be the number of successful estimates in a scenario. The empirical bias is the mean estimate minus the known target. The empirical variance describes replicate-to-replicate dispersion; the empirical standard deviation is its square root. These are distinct from the average standard error reported by the estimator.

Here \(\mathcal{S}\) is the set of successful replicates and \(\bar{\hat\psi}_j\) is their mean estimate. Report the number of attempts and successes as well as these conditional summaries. If failures are excluded from the bias calculation, the reported bias describes successful runs only; it should not conceal a high failure rate.

For a nominal 95% confidence interval, coverage is the proportion of intervals that contain \(\psi_0\). State how failures enter this calculation. One operational choice is to divide the number of covering intervals by all attempted replicates and count a failed replicate as not covered. Also report coverage among successful intervals, so readers can distinguish interval performance from computational failure. The failure rate is the number of failed attempts divided by the total number attempted. Define failure in advance, such as a nonfinite estimate, a nonconvergent targeting step, or inability to produce the prespecified interval; do not remove inconvenient results after seeing them.

Worked Performance Calculations

A simulation also has Monte Carlo error: repeating the entire simulation with a different random seed would produce slightly different performance summaries. For a coverage estimate \(\hat c\) based on \(R\) independent attempted replicates, a useful approximate Monte Carlo standard error is \(\sqrt{\hat c(1-\hat c)/R}\), when coverage is coded for every attempt according to the chosen failure rule. This follows from treating each replicate's coverage indicator as a Bernoulli variable with variance \(c(1-c)\); the variance of its sample mean is \(c(1-c)/R\). Substituting \(\hat c\) gives the estimated standard error.

Compare Methods Fairly and Interpret the Results

The strongest design compares AIPW and TMLE on the same replicate data sets, with the same target, nuisance-model information, and interval convention. This pairing reduces noise in their difference: an unusually difficult simulated data set affects both methods rather than one method being evaluated on a separate draw. Keep the data-generation code, fitting rules, interval formula, and failure definitions fixed across methods wherever the comparison permits.

The following hypothetical results illustrate how a report might separate conditions. They are not universal performance claims. Bias and standard deviation are shown in ATE units; coverage uses all 1000 attempts, with failures counted as not covered.

Practical Design and Common Pitfalls

A useful simulation plan specifies the estimand, data-generating mechanisms, nuisance-model correctness conditions, overlap settings, sample sizes, number of replicates, interval procedure, and failure rule before examining results. When varying one feature, hold others fixed if the purpose is to isolate that feature; factorial designs are also useful when interactions are of interest. In particular, weak overlap may be more consequential at small sample sizes, so reporting only a single sample size can miss an important pattern.

• Changing the target accidentally. Keep the target parameter fixed when changing overlap or sample size, unless a changed target is an explicit part of the design.
• Reporting only successful fits. A conditional bias or coverage summary cannot reveal how often the method failed. Show the attempted denominator and failure rate.
• Calling a finite-sample result a general guarantee. A simulation evaluates the specified mechanisms and fitting procedures; it cannot prove universal performance or identify which nuisance model is correct in an application.
• Comparing methods on different simulated samples. Unpaired runs add avoidable variation to method comparisons. Use shared replicate data and report paired differences when comparing summaries.
• Equating empirical spread with standard-error quality. Empirical standard deviation measures variation across replicates; average estimated standard error assesses the estimator's uncertainty estimates. Compare them, and assess coverage, rather than substituting one for another.