Targeted Maximum Likelihood Estimation (TMLE)

Targeted maximum likelihood estimation (TMLE)

is a general framework (van der Laan & Rubin, 2006) for estimating a low-dimensional causal target parameter — most commonly the average treatment effect (ATE), E[Y^1] − E[Y^0] — from observational data with high-dimensional nuisance parameters.

It is a plug-in estimator: it produces an estimate of the target by substituting fitted components into the parameter's mapping, but unlike a naive plug-in it adds a targeting (fluctuation) step that updates the initial fit specifically to reduce bias in the target parameter and to make the final estimate solve the efficient influence-curve (EIC) estimating equation.

Solving the EIC equation is what delivers two properties at once: (1) double robustness — the ATE estimate is consistent if either the outcome model or the propensity (treatment) model is correctly specified, not necessarily both; and (2) semiparametric efficiency — when both nuisance models are consistent at fast-enough rates, TMLE achieves the lowest possible asymptotic variance, and its variance is estimated directly from the empirical variance of the influence curve, giving valid Wald confidence intervals.

The algorithm, concretely (binary treatment A, outcome Y, covariates W)

1. Fit an initial outcome regression Qbar0(A,W) = E[Y | A, W] and obtain the predicted outcomes under treatment Qbar0(1,W) and control Qbar0(0,W).
2. Fit the propensity score g(W) = P(A=1 | W).
3. Construct the clever covariate H(A,W) = A/g(W) − (1−A)/(1−g(W)) (the score of the EIC for the ATE).

(4) Run the fluctuation: a no-intercept logistic regression of Y on the clever covariate with the initial prediction Qbar0 as an offset (on the logit scale), yielding a single fitted coefficient epsilon. (5) Update: Qbar*(A,W) = expit( logit(Qbar0(A,W)) + epsilon · H(A,W) ); compute the updated potential-outcome predictions Qbar*(1,W) and Qbar*(0,W).

(6) The TMLE of the ATE is the mean difference (1/n) Σ [Qbar*(1,Wi) − Qbar*(0,Wi)]. The variance is the sample variance of the estimated influence curve, IC_i = H(Ai,Wi)·(Yi − Qbar*(Ai,Wi)) + (Qbar*(1,Wi) − Qbar*(0,Wi)) − ATE, divided by n. Because step (1) and step (2) can each be fit with Super Learner (cross-validated ensemble stacking of flexible learners — GBM, random forests, splines, elastic net), TMLE allows machine learning for the nuisance functions while the targeting step restores valid, EIC-based inference that plain ML plug-in estimates do not have.

Core conceptual distinctions

TMLE sits in the doubly-robust family with AIPW (augmented inverse-probability weighting) and is related to g-computation (the parametric outcome-regression standardization). (1) vs g-computation: g-computation is the un-targeted plug-in — fit Qbar, standardize — and is consistent only if the outcome model is correct, with no double robustness and no influence-curve-based variance. TMLE is g-computation plus a targeting update that buys double robustness and valid inference. (2) vs AIPW: both are doubly robust and asymptotically efficient and both use the EIC.

AIPW is a one-step additive bias correction (it adds the mean of the influence function to the plug-in), so its estimate of a probability/risk can fall outside [0,1]; TMLE is a substitution estimator whose targeting step keeps the updated outcome predictions within the model's natural bounds, which is more stable in finite samples and under near-positivity violations. (3) TMLE estimates a marginal causal contrast (the population ATE, or ATT/risk ratio/odds ratio via the appropriate clever covariate and mapping), not a conditional regression coefficient.

Pros, cons, and trade-offs

• vs `propensity-score-methods-psm-iptw` (IPTW alone): TMLE is doubly robust — it tolerates misspecification of either the treatment or outcome model — and is efficient, whereas IPTW is consistent only if the propensity model is correct and is sensitive to extreme weights/positivity violations. Prefer TMLE when you want robustness to one model being wrong and you can credibly estimate both nuisance functions; prefer plain IPTW for transparency or when the outcome model is hard to specify and you trust the PS.
• vs AIPW (one-step doubly-robust): Both are doubly robust and efficient; TMLE's substitution/bounding makes it more stable when predicted probabilities are near 0/1 and respects the parameter space, while AIPW is simpler to implement and its theory is more transparent. Prefer TMLE under near-positivity or bounded outcomes; AIPW when simplicity and a closed-form correction are preferred. They are close cousins, not opposites.
• vs `g-estimation-structural-nested-models`: g-estimation of structural nested models targets effects in the presence of time-varying confounding affected by prior treatment and is natural for continuous/longitudinal exposures; TMLE (and its longitudinal extension, LTMLE) targets marginal contrasts and integrates ML nuisance estimation. Prefer g-estimation / SNMs for effect modification by time-varying covariates and structural-nested questions; prefer (L)TMLE for marginal parameters with flexible ML nuisance fits.
• vs `predictive-and-causal-ml-models-rwe` (plain ML plug-in): A bare machine-learning outcome model plugged into a standardization gives a biased estimate with invalid inference (the bias of ML for prediction does not vanish for the target parameter). TMLE's targeting step is exactly what corrects that bias and restores EIC-based CIs — this is the reason TMLE exists.

When NOT to use — and when it is actively misleading or dangerous

• Under serious positivity violations. If g(W) is near 0 or 1 for some covariate strata (no comparable treated/untreated units), the clever covariate explodes, the fluctuation is unstable, and TMLE — like any IPW-based method — extrapolates; truncating the PS hides rather than fixes the violation. Diagnose and report positivity before trusting the estimate.
• When neither nuisance model can be credibly estimated. Double robustness protects against one wrong model, not two; if both the outcome and treatment models are badly misspecified (and ML cannot rescue them because of unmeasured confounding or too little data), TMLE is biased like everything else. It does not manufacture identification.
• As a fix for unmeasured confounding. TMLE assumes conditional exchangeability (no unmeasured confounding) given W; it is an estimation method, not an identification one. Reporting a tight TMLE CI on a confounded contrast launders bias into false precision.
• Without cross-fitting when using aggressive ML. Highly adaptive learners without sample-splitting (cross-fitting / CV-TMLE) can overfit the nuisance functions and invalidate the asymptotics; use CV-TMLE when learners are complex.
• For a conditional/effect-modified parameter without the right target. TMLE estimates the parameter you specify; using a marginal-ATE TMLE and narrating it as a subgroup effect is a target-parameter error.

Data-source operational depth

• Claims (FFS vs MA): The high-dimensional covariate bank (diagnoses, prior fills, utilization) is ideal Super Learner input for both nuisance models, and high-dimensional propensity-score covariate selection can feed W. Build the contrast on FFS-observable person-time only — Medicare Advantage enrollees lack fee-for-service claims, so MA-only spans give differentially incomplete confounders and a positivity structure that is an artifact of missingness, not biology. Check the estimated PS distribution by arm for positivity before targeting.
• EHR: Labs, vitals, and NLP-derived severity sharpen both Qbar and g (richer confounding control, more plausible exchangeability), but informative-presence/visit-driven capture creates a selection process W does not capture; treat informative loss to follow-up with the longitudinal extension (LTMLE) and censoring weights rather than a single cross-sectional TMLE.
• Registry / linked: Adjudicated outcomes and severity tighten exchangeability and improve the outcome model; linked claims-EHR-vital-records is the strongest substrate for both nuisance fits, with linkage selection reported as a separate bias. The arithmetic of the targeting step is unchanged; what changes is the credibility of the no-unmeasured-confounding assumption.

Interpreting the output

Consider the diabetes hospitalization study above. TMLE using SuperLearner for both nuisance models reports an adjusted risk difference (ATE) of −0.06 (95% CI −0.10, −0.02) for 1-year hospitalization, with confidence intervals derived from the efficient influence curve.

Formal interpretation: The TMLE estimate of −0.06 is the average treatment effect (ATE) — the population-averaged difference in 1-year hospitalization probability if all 1,000 patients had received the new drug versus if none had. It is doubly robust: the estimate converges to the true ATE if either the outcome model (predicting hospitalization from covariates and treatment) or the propensity model (predicting treatment receipt from covariates) is correctly specified, even if the other is misspecified.

The influence-curve-based 95% CI is valid even when machine-learning models are used for both nuisance functions — standard parametric CI formulas do not apply when SuperLearner is used without the targeting step.

The central untestable assumption is no unmeasured confounding: every variable that jointly determines treatment assignment and hospitalization must be measured and included. The CI means that across repeated datasets from the same data-generating process, intervals constructed this way would contain the true ATE in 95% of replications.

Practical interpretation: Receiving the new diabetes drug is estimated to reduce 1-year hospitalization risk by 6 percentage points on average across the full patient population. TMLE's double-robustness property means that a moderately misspecified outcome model does not invalidate the estimate if the propensity model is adequate, and vice versa — providing a meaningful safeguard that neither IPTW alone nor outcome regression alone offers.

Decision diagram

flowchart TD
  Q0[Initial outcome regression<br/>Qbar0 A,W via Super Learner] --> H
  G[Propensity model<br/>g W = P A=1 given W] --> H[Clever covariate<br/>H = A/g - 1-A/1-g]
  H --> Fluc[Fluctuation: logistic Y on H<br/>offset = logit Qbar0 -> epsilon]
  Fluc --> Upd[Update Qbar* = expit logit Qbar0 + epsilon H]
  Upd --> ATE[ATE = mean Qbar* 1,W - Qbar* 0,W]
  ATE --> IC[Variance from efficient influence curve<br/>valid Wald 95% CI]

flowchart TD
  Goal{What do you need?} -->|Prediction only| ML[Plain ML model<br/>biased + invalid CI for causal effect]
  Goal -->|Marginal causal contrast| DR{Can you fit BOTH<br/>outcome and propensity models?}
  DR -->|Yes| Pos{Positivity OK?<br/>g W away from 0 and 1}
  Pos -->|Yes| TMLE[TMLE / AIPW<br/>doubly robust + efficient + IC inference]
  Pos -->|No near-positivity| Warn[Diagnose positivity first<br/>truncation hides not fixes it]
  DR -->|Outcome model trusted only| Gcomp[g-computation<br/>no double robustness]
  DR -->|Propensity trusted only| IPTW[IPTW<br/>consistent if PS correct]

Worked example

Scenario

Imagine a cohort study of 1,000 patients to estimate whether a new oral diabetes drug lowers one-year hospitalization risk compared to usual care. Patients who are sicker tend to receive the new drug more often (confounding by indication). We compare three approaches: plain logistic regression (one model, outcome only), and TMLE (two models combined with a targeting step).

Dataset

Simplified summary of what the analyst observes for each patient. There is no pharmacy fill timeline here; the raw inputs are patient background variables, treatment assignment, and a binary outcome.

Steps

Result

In a stylized version of this scenario, plain regression might estimate the new drug reduces hospitalization risk by about 4 percentage points, but if the regression model misspecified how HbA1c relates to hospitalization, that estimate is biased with no fallback. TMLE, using both models together, produces a doubly robust estimate: as long as either the outcome model or the propensity model captures the main confounding relationships reasonably well, the risk difference estimate converges to the true value. In practice, with machine learning for both models and the targeting step, TMLE often gives better-calibrated estimates and honest confidence intervals in the same high-dimensional claims or EHR setting where plain regression is most likely to be misspecified.

Trade-offs

Runnable example

import numpy as np
import pandas as pd
from zepid.causal.doublyrobust import TMLE

# Illustrative cohort: confounders W1, W2 affect both treatment A and outcome Y.
rng = np.random.default_rng(11)
n = 3000
w1 = rng.normal(0, 1, n)
w2 = rng.binomial(1, 0.45, n)
a  = rng.binomial(1, 1 / (1 + np.exp(-(-0.3 + 0.8*w1 + 0.6*w2))))   # treatment
# True ATE on the risk scale is ~ +0.08 (A raises outcome risk by ~8 points).
py = 1 / (1 + np.exp(-(-1.0 + 0.5*a + 0.9*w1 + 0.7*w2)))
y  = rng.binomial(1, py)
df = pd.DataFrame(dict(A=a, Y=y, W1=w1, W2=w2))

# Single-step TMLE: specify the outcome model, then the propensity model, then fit.
tmle = TMLE(df, exposure="A", outcome="Y")
tmle.outcome_model("A + W1 + W2", print_results=False)      # initial Qbar(A,W)
tmle.exposure_model("W1 + W2", print_results=False)         # propensity g(W)
tmle.fit()                                                  # targeting + IC inference

print(f"ATE (risk difference) = {tmle.risk_difference:.4f}")
print(f"95% CI = ({tmle.risk_difference_ci[0]:.4f}, {tmle.risk_difference_ci[1]:.4f})")
print(f"Risk ratio = {tmle.risk_ratio:.3f}")

library(tmle)        # install.packages("tmle"); pulls in SuperLearner
set.seed(11)

# Illustrative cohort: confounders W1, W2 affect both treatment A and outcome Y.
n  <- 3000
W1 <- rnorm(n); W2 <- rbinom(n, 1, 0.45)
A  <- rbinom(n, 1, plogis(-0.3 + 0.8*W1 + 0.6*W2))            # treatment
Y  <- rbinom(n, 1, plogis(-1.0 + 0.5*A + 0.9*W1 + 0.7*W2))    # binary outcome
W  <- data.frame(W1 = W1, W2 = W2)

# Super Learner library for the Q (outcome) and g (propensity) fits, then targeting.
fit <- tmle(Y = Y, A = A, W = W,
            Q.SL.library = c("SL.glm", "SL.glmnet", "SL.mean"),
            g.SL.library = c("SL.glm", "SL.glmnet", "SL.mean"),
            family = "binomial")

print(fit$estimates$ATE)     # point estimate, IC-based variance, and 95% CI
cat(sprintf("ATE = %.4f (95%% CI %.4f, %.4f)\n",
            fit$estimates$ATE$psi,
            fit$estimates$ATE$CI[1], fit$estimates$ATE$CI[2]))

/* (1) Initial outcome model Qbar(A,W); predict under observed A, under A=1, and under A=0. */
data coh1; set work.coh; run;
proc logistic data=coh1 noprint;
  model Y(event='1') = A W1 W2;
  output out=q_obs p=q_aw;                 /* Qbar(A_i, W_i) */
run;
data c1; set work.coh; A=1; run;  data c0; set work.coh; A=0; run;
proc logistic data=coh1 noprint; model Y(event='1') = A W1 W2; score data=c1 out=q1(rename=(P_1=q1)); run;
proc logistic data=coh1 noprint; model Y(event='1') = A W1 W2; score data=c0 out=q0(rename=(P_1=q0)); run;

/* (2) Propensity model g(W) = P(A=1|W). */
proc logistic data=work.coh noprint;
  model A(event='1') = W1 W2;
  output out=g_out p=g;
run;

/* (3) Clever covariate H = A/g - (1-A)/(1-g); offset = logit(Qbar). */
data tgt;
  merge q_obs g_out q1(keep=q1) q0(keep=q0);
  H      = A/g - (1-A)/(1-g);
  offset = log(q_aw/(1-q_aw));             /* logit of the initial prediction */
run;

/* (4) Fluctuation: no-intercept logistic of Y on H with logit(Qbar) as offset -> epsilon. */
proc logistic data=tgt noprint;
  model Y(event='1') = H / noint offset=offset;
  ods output ParameterEstimates=eps_out;
run;
data _null_; set eps_out; if Variable='H' then call symputx('eps', Estimate); run;

/* (5) Update predictions, form ATE, and compute influence-curve variance. */
data fin; set tgt;
  H1 = 1/g;  H0 = -1/(1-g);
  qs1 = 1/(1+exp(-(log(q1/(1-q1)) + &eps*H1)));   /* Qbar*(1,W) */
  qs0 = 1/(1+exp(-(log(q0/(1-q0)) + &eps*H0)));   /* Qbar*(0,W) */
  qs_obs = 1/(1+exp(-(offset + &eps*H)));         /* Qbar*(A,W) */
  diff = qs1 - qs0;
run;
proc means data=fin noprint; var diff; output out=atem mean=ate; run;
data _null_; set atem; call symputx('ate', ate); run;
data icv; set fin;
  ic = H*(Y - qs_obs) + (qs1 - qs0) - &ate;       /* efficient influence curve */
run;
proc means data=icv noprint; var ic; output out=v var=icvar n=n; run;
data result; set v; ate=&ate; se=sqrt(icvar/n);
  lcl=ate-1.96*se; ucl=ate+1.96*se;
run;
proc print data=result noobs; var ate se lcl ucl; run;

Citations

1. [1]Schuler MS, Rose S. Targeted maximum likelihood estimation for causal inference in observational studies. American Journal of Epidemiology. 2017;185(1):65-73.
2. [2]van der Laan MJ, Rubin D. Targeted maximum likelihood learning. The International Journal of Biostatistics. 2006;2(1):Article 11.
3. [3]Funk MJ, Westreich D, Wiesen C, Sturmer T, Brookhart MA, Davidian M. Doubly robust estimation of causal effects. American Journal of Epidemiology. 2011;173(7):761-767.

1. [4]Luque‐Fernandez MA, Schomaker M, Rachet B, Schnitzer ME. Targeted maximum likelihood estimation for a binary treatment: a tutorial. Statistics in Medicine. 2018;37(16):2530-2546.

1. [5]Gruber S, van der Laan MJ. tmle: an R package for targeted maximum likelihood estimation. Journal of Statistical Software. 2012;51(13):1-35.