Counterfactual Evaluation in Ads: IPS, SNIPS, and Doubly Robust

You have a new ranking model. You want to know if it’s better than the one in production before you ship it. The honest answer is: you should run an A/B test. But A/B tests take weeks, they require splitting real traffic, and they expose users to an unvalidated model. For a recommender system serving hundreds of millions of impressions per day, that cost is real.

Counterfactual evaluation offers a shortcut: use the logs you already have to estimate how a new policy would have performed. Those logs come from your logging policy — the production policy that chose which ad to show, and therefore the only data we have to learn from, even when the policy we want to evaluate would have made different choices. The math is elegant, the engineering is subtle, and the assumptions are easy to violate in ways that will mislead you completely.

This article covers four estimators — Direct Method (DM), Inverse Propensity Scoring (IPS), Self-Normalizing IPS (SNIPS), and Doubly Robust (DR) — with verified implementations and the numbers that show exactly when each one breaks.

Disclaimer: The opinions expressed in this article are my own and do not represent the views of Google. This content is based solely on publicly available information.

The Big Picture

Before diving into formulas and code, here is the full counterfactual evaluation pipeline. Production generates logs with four fields per impression; an estimator turns those logs into a single number — the estimated policy value of a candidate policy you have not yet deployed. Each path through below diagram makes a different bet: DM bets on your reward model being accurate; IPS bets on your propensities being correct; DR hedges both bets simultaneously. The rest of this article unpacks what those bets mean and when each one pays off.

The Setup: Logged Bandit Data

Before examining the four estimators, it helps to define the data structure they all operate on and the exact quantity each one is trying to estimate. Every recommendation system logs its decisions. For each impression you recorded:

• Context x: user features, query, session
• Action a: the item or ad shown (chosen by the logging policy π₀)
• Reward r: whether the user clicked, converted, or engaged
• Propensity π₀(a|x): the probability that the logging policy chose this action

The goal is to estimate the expected reward of a target policy π_new — one you have not deployed yet — using only this logged data.

The ground truth:

V(π_new) = E_{x~D}[ ∑_a π_new(a|x) · r(x, a) ]

In our verified simulation: a 3-action bandit where the logging policy strongly prefers action 0 (80% probability), but action 1 is actually better (true reward ≈ 0.60 vs 0.20). The target policy is uniform (33% each). Ground truth V(π_new) = 0.2996.

The mismatch between the logging policy and the target policy is the core challenge. The logs are dominated by observations from action 0. Any estimator that naively averages logged rewards will conclude “action 0 is fine” — missing the fact that action 1 is three times better but was shown only 10% of the time. This is selection bias, and every estimator below attempts to correct for it in a different way.

Simulation Setup

The simulate_bandit function generates a realistic biased logging scenario. Action 0 gets 80% of the traffic even though action 1 is three times better — exactly what happens when a production system exploits a known-good item and rarely explores alternatives. The ground truth is computed by averaging rewards under the uniform policy across all samples, giving V = 0.2996. Every estimator below is judged against this known value.

from __future__ import annotations
from dataclasses import dataclass
import numpy as np

@dataclass
class BanditData:
    contexts: np.ndarray       # shape (n, context_dim)
    actions: np.ndarray        # shape (n,) — chosen action index
    rewards: np.ndarray        # shape (n,) — observed 0/1 reward
    propensities: np.ndarray   # shape (n,) — π₀(a|x) for chosen action

def simulate_bandit(
    n_samples: int = 10_000,
    n_actions: int = 3,
    logging_bias: float = 0.80,
    seed: int = 42,
) -> tuple[BanditData, np.ndarray, float]:
    """
    Returns logged data, true reward table (n_samples x n_actions), and ground truth V.
    Logging policy prefers action 0 with probability logging_bias.
    Action 1 has the highest true reward (~0.60); action 0 is mediocre (~0.20).
    Target policy is uniform: 1 / n_actions for every action.
    """
    rng = np.random.default_rng(seed)
    contexts = rng.standard_normal((n_samples, 8)).astype(np.float32)

    # True reward table — action 1 dominates (~0.60), action 0 mediocre (~0.20).
    base = np.array([0.2, 0.6] + [0.1] * (n_actions - 2))
    noise = rng.standard_normal((n_samples, n_actions)) * 0.05
    true_rewards = np.clip(base + noise, 0.0, 1.0).astype(np.float32)

    # Logging policy: probability `logging_bias` on action 0, rest split evenly.
    logging_probs = np.full(n_actions, (1 - logging_bias) / (n_actions - 1))
    logging_probs[0] = logging_bias

    actions = rng.choice(n_actions, size=n_samples, p=logging_probs)
    propensities = logging_probs[actions].astype(np.float32)
    rewards = rng.binomial(
        1, true_rewards[np.arange(n_samples), actions]
    ).astype(np.float32)

    data = BanditData(contexts=contexts, actions=actions,
                      rewards=rewards, propensities=propensities)

    # Ground truth under uniform target policy.
    target_probs = np.full(n_actions, 1.0 / n_actions)
    ground_truth = float((true_rewards * target_probs).sum(axis=1).mean())

    return data, true_rewards, ground_truth

The Four Estimators

With the simulation environment established, we can now implement and compare the four estimators — each making a different trade-off between bias, variance, and robustness to assumption violations.

Direct Method (DM)

You have observations only for actions the logging policy chose. If action 1 was shown just 10% of the time, you have very few (context, action 1, reward) examples. DM sidesteps this by fitting a reward model on all logged data, then using that model to predict rewards for every action — including the ones rarely shown.

Train a reward model r̂(x, a) on the logged data, then pretend it is the oracle and ask: “what would the expected reward be if the target policy chose actions according to π_new?” No importance weighting needed; just evaluate the model across all actions.

When data is scarce, importance-weighted estimators have enormous variance (one bad weight can dominate). DM avoids this by using model predictions rather than re-weighting observations, giving low variance at the cost of trusting the model. Train a reward model r̂(x, a) on the logged data, then use it to predict the target policy’s reward directly:

V_DM = E_x[ ∑_a π_new(a|x) · r̂(x, a) ]

# BanditData and simulate_bandit defined in the Simulation Setup block above (BanditData provides the .actions and .propensities arrays used here)

def importance_weights(
    data: "BanditData",
    target_probs: np.ndarray,
    clip: float = np.inf,
) -> np.ndarray:
    w = target_probs[data.actions] / np.maximum(data.propensities, 1e-9)
    return np.minimum(w, clip)

def ips(
    data: "BanditData",
    target_probs: np.ndarray,
    clip: float = np.inf,
) -> float:
    """IPS: unbiased under correct propensities; high variance at low propensity."""
    w = importance_weights(data, target_probs, clip=clip)
    return float((w * data.rewards).mean())

• Advantage: No variance inflation — there are no importance weights. Works well with large datasets and accurate reward models.
• Fatal flaw: If r̂ is biased, V_DM inherits that bias with no correction mechanism. In our simulation, a model that underestimates the good action’s reward by 40% produces a DM estimate 26.7% away from ground truth — and more data does not help, because the model error is structural.

Inverse Propensity Scoring (IPS)

The logged data is a biased sample — actions the logging policy liked are over-represented. IPS de-biases the sample by re-weighting each observation to reflect how likely the target policy would have been to make the same choice.

Think of importance weighting as rebalancing a survey. If you surveyed mostly action-0 users but want to know what a uniform-policy world looks like, you down-weight action-0 observations (they are over-sampled) and up-weight action-1 observations (they are under-sampled). The weight w_i = π_new(a_i|x_i) / π_0(a_i|x_i) is exactly this rebalancing factor.

DM requires an accurate model; IPS requires accurate propensities. If your propensities are logged correctly (which is achievable with discipline — see the engineering section below), IPS is mathematically guaranteed to be unbiased regardless of how wrong your reward model is. Weight each logged reward by how much more (or less) likely the target policy was to choose that action compared to the logging policy:

V_IPS = E_i[ (π_new(a_i|x_i) / π_0(a_i|x_i)) · r_i ]

The ratio w_i = π_new / π_0 is the importance weight. When the logging policy rarely showed an action that the target policy likes, the weight is large — amplifying the signal from those rare observations.

# BanditData and simulate_bandit defined in the Simulation Setup block above (BanditData provides the .actions and .propensities arrays used here)

def importance_weights(
    data: "BanditData",
    target_probs: np.ndarray,
    clip: float = np.inf,
) -> np.ndarray:
    w = target_probs[data.actions] / np.maximum(data.propensities, 1e-9)
    return np.minimum(w, clip)

def ips(
    data: "BanditData",
    target_probs: np.ndarray,
    clip: float = np.inf,
) -> float:
    """IPS: unbiased under correct propensities; high variance at low propensity."""
    w = importance_weights(data, target_probs, clip=clip)
    return float((w * data.rewards).mean())

• Advantage: Unbiased — if propensities are logged correctly, IPS converges to the true policy value regardless of how different π_new is from π₀.
• Fatal flaw: Variance explodes when propensities are small. If the logging policy almost never showed action 1 (propensity = 0.01) but the target policy shows it 33% of the time, the importance weight is 33. A few lucky observations with high weight dominate the estimate. The weight distribution histogram in Figure 1 makes this visible.

SNIPS (Self-Normalizing IPS)

IPS can produce estimates outside the valid reward range (e.g., negative click rates) when a few extreme weights dominate. The weights should sum to n in expectation, but with finite samples they sum to something different — and this discrepancy injects extra variance.

Instead of dividing by the number of samples (which IPS does), divide by the actual sum of weights. This self-normalization ensures the weights always average to 1.0, preventing any single observation from contributing more than its share. It is the difference between taking a weighted average (SNIPS) and a weighted sum divided by n (IPS).

In practice, at logging bias = 0.90 or higher, IPS standard deviation roughly doubles relative to SNIPS. The normalization costs a small asymptotic bias — the kind that disappears as n grows — while buying meaningful variance reduction at every finite sample size. Divide by the sum of weights instead of the count of samples:

V_SNIPS = ∑ w_i · r_i / ∑ w_i

# BanditData and simulate_bandit defined in the Simulation Setup block above; importance_weights defined in the IPS block above

def snips(
    data: "BanditData",
    target_probs: np.ndarray,
    clip: float = np.inf,
) -> float:
    """SNIPS: lower variance than IPS via self-normalization; small asymptotic bias."""
    w = importance_weights(data, target_probs, clip=clip)
    return float((w * data.rewards).sum() / max(w.sum(), 1e-9))

• Advantage: Consistently lower variance than IPS — especially at high bias. The normalisation prevents a single massive weight from dominating.
• Trade-off: Introduces a small asymptotic bias (it converges to the correct answer, but more slowly in some regimes). In practice, this bias-variance trade-off is almost always worth it.

Doubly Robust (DR)

You have both a reward model and propensity logs, but neither is perfect. DM fails when the model is wrong; IPS fails when propensities are wrong. DR is designed to be correct when either component is correct — halving the probability of catastrophic failure.

DR uses the DM estimate as a baseline and then corrects it using IPS — but only on the residual (the gap between the model’s prediction and the actual observed reward). If the model is perfect, the residual is zero and DR = DM. If the model is imperfect but propensities are correct, the IPS correction picks up the model’s error on logged actions. You only need one of the two components to be trustworthy for the full estimator to be approximately unbiased.

In production, you rarely know which of your two inputs (model or propensities) is more reliable. DR provides insurance: it degrades gracefully when one component fails, whereas IPS and DM each fail completely when their single assumption is violated.

Combine both a reward model and importance weighting:

V_DR = E_x[ ∑_a π_new(a|x) · r̂(x,a) ] + E_i[ w_i · (r_i − r̂(x_i, a_i)) ]

The first term is the DM estimate. The second is an IPS-weighted correction for the model’s error on logged actions.

# BanditData and simulate_bandit defined in the Simulation Setup block above; importance_weights defined in the IPS block above; reward_model constructed in the Direct Method block above

def doubly_robust(
    data: "BanditData",
    target_probs: np.ndarray,
    reward_model: np.ndarray,
    clip: float = np.inf,
) -> float:
    """DR: unbiased if reward model OR propensities are correct (doubly robust)."""
    w      = importance_weights(data, target_probs, clip=clip)
    dm     = (reward_model * target_probs).sum(axis=1)
    r_hat  = reward_model[np.arange(len(data.actions)), data.actions]
    resid  = w * (data.rewards - r_hat)
    return float((dm + resid).mean())

If the reward model is perfect, the residual term is zero and DR = DM. If the propensities are correct, the residual term is an unbiased correction of the model’s error. DR is unbiased if either component is correct — this is the “doubly robust” property.

The Numbers

The four estimators have very different theoretical properties — now we verify those properties against a simulation where the ground truth is known. From a verified simulation run (n=10,000, 3 actions, logging bias=80%):

With ground truth of0.2996, IPS, SNIPS, and DR all come within 1% because the propensities are correct in this simulation. DM fails because the reward model underestimates action 1 — and there’s nothing to correct it.

The 26.7% error on DM is not a quirk of this example. It is a direct consequence of the model underestimating action 1’s reward by 40% (multiplied by a factor of 0.6 in the simulation). The error is proportional to the model’s bias on the action the target policy prefers — exactly the action that was least represented in training data, which is the direction a biased reward model almost always fails.

The Variance Explosion

The numbers above assume a fixed logging bias of 80%. The next question is how estimator reliability degrades as that bias increases — the regime where most production systems operate because logging policies are rarely close to uniform.

The left panel of Figure 1 tells the story through the weight distribution shape. At balanced logging (bias=0.50, green), importance weights cluster tightly around 1.0 — the target policy and logging policy are similar enough that re-weighting barely changes anything. As the logging policy concentrates on action 0 (bias=0.95, amber), a secondary spike appears near the maximum weight: the rare observations of action 1 and action 2 receive disproportionate weight. At extreme bias (bias=0.99, red), the distribution has a fat right tail extending past 10×, meaning a handful of observations can swing the estimate by 10 reward points.

The right panel translates this into variance numbers. IPS standard deviation (50 trials, n=500) rises from 0.024 at balanced logging to 0.190 at bias=0.99 — roughly 8× worse. SNIPS consistently tracks below IPS, but the gap between them also widens: at bias=0.99, SNIPS is about 84% of IPS variance rather than 92% at low bias. Self-normalization helps most precisely where the problem is worst. The max weight axis (amber, right scale) shows why: as bias goes from 0.50 to 0.99, the theoretical maximum importance weight jumps from 1.3× to 66.7×. This is why clipping matters. Setting clip=M bounds the worst-case weight at M, trading a small bias for a large variance reduction. The optimal clip threshold is data-dependent; a common heuristic is clip = sqrt(n).

Figure 1: IPS variance explodes as the logging policy concentrates on a single action.

Estimator Robustness Across Failure Modes

The variance plot answers “what happens when the logging policy is biased?” This figure answers the orthogonal question: “what happens when your model or your propensities are wrong?” As shown in Figure 2, reading from left to right, the columns represent four diagnostic scenarios: both components correct (column 1), wrong model only (column 2), wrong propensity only (column 3), and both wrong (column 4). The heights directly quantify how much each estimator is hurt by each failure.

Figure 2: Mean absolute error of each estimator (40 trials, n=3,000) across four model × propensity correctness combinations.

Column 1 (correct model, correct propensity) is the baseline: all estimators perform well, with DR and IPS at roughly 0.009–0.011. DM’s slightly higher error here is normal; model predictions are still imperfect even when well-calibrated, and DM has no correction mechanism.

The doubly-robust property is starkest in column 2 (wrong model, correct propensity). DM’s error spikes to 0.071 — the model underestimates action 1 by 70%, and DM propagates that error directly. But DR stays around 0.011 — nearly as accurate as IPS alone — because the IPS-weighted residual correction repairs the model’s error on logged actions. DR’s IPS component sees the discrepancy between model predictions and actual rewards and corrects for it.

Column 3 (correct model, wrong propensity) shows DR’s other protection. Here IPS and SNIPS suffer because their weights are based on incorrect propensities. DR collapses to 0.009 — essentially matching DM performance — because the IPS residual term goes nearly to zero whenever the model is right, leaving DR ≈ DM regardless of propensity quality.

Column 4 (both wrong) is where DR finally fails: it inherits both the model’s bias and the propensity-amplified residual error, ending at 0.041 — worse than any single-component estimator. This is the scenario where the doubly robust guarantee breaks down: you need at least one component to be trustworthy. The key practical insight: the worst thing you can do is use IPS with incorrect propensities. If you do not log propensities — or you log them incorrectly — every IPS-based estimator breaks. DM, for all its limitations, at least fails gracefully.

Sample Efficiency

The log-log scale makes convergence rates visible as slopes as shown in Figure 3. IPS, SNIPS, and DR all fall along straight lines with slope close to −0.5 — the O(1/√n) convergence rate you expect from unbiased estimators. At n=100, all three cluster around 0.055–0.058 mean absolute error. By n=10,000 they have fallen to 0.005–0.007. The gap between DR and IPS is consistent across every sample size: DR is roughly 25–30% more accurate, which means you can hit the same target accuracy with about half the data if you use DR instead of IPS. For a system that logs millions of impressions per day this difference is irrelevant; for a small experiment with 1,000 logged impressions, it is meaningful.

Figure 3: Mean absolute error vs. dataset size on a log-log scale for all four estimators.

DM tells a completely different story. It starts at a comparable error to the other estimators at n=100, but as sample size grows, error stops decreasing around 0.10 and plateaus. More data does not help because the reward model’s underestimation of action 1 is a structural bias — every additional sample is used to train a model that is systematically wrong in the same direction. The DM flatline is the visual proof that model bias, unlike variance, does not average out. It is the most important practical lesson in this figure: IPS, SNIPS, and DR all continue to get better with more data; DM does not.

Weight Clipping in Practice

The variance explosion shown in Figure 1 calls for a practical mitigation. Weight clipping is a bias-variance knob that bounds the influence of extreme importance weights, and selecting the right clip threshold determines how much bias you trade for variance control. Without clipping, a single observation with weight 50× dominates 50 normal observations. With hard clipping at M, you cap the influence of any single data point at the cost of introducing bias proportional to how much mass exceeds the clip threshold.

# BanditData, simulate_bandit, and ips defined in earlier blocks

import numpy as np

def clip_analysis(
    data: "BanditData",
    target_probs: np.ndarray,
    ground_truth: float,
    clip_values: list[float],
    n_trials: int = 50,
    n_samples: int = 500,
    seed: int = 0,
) -> dict[float, dict[str, float]]:
    """Measure IPS bias and std dev at each clip threshold."""
    rng = np.random.default_rng(seed)
    n = len(data.actions)
    results = {}
    for clip in clip_values:
        estimates = []
        for _ in range(n_trials):
            idx = rng.choice(n, size=n_samples, replace=False)
            sub = BanditData(
                contexts=data.contexts[idx],
                actions=data.actions[idx],
                rewards=data.rewards[idx],
                propensities=data.propensities[idx],
            )
            estimates.append(ips(sub, target_probs, clip=clip))
        arr = np.array(estimates)
        results[clip] = {
            "bias": float(np.mean(arr) - ground_truth),
            "std":  float(np.std(arr)),
            "rmse": float(np.sqrt(np.mean((arr - ground_truth) ** 2))),
        }
    return results

And here is rule of thumb for clip threshold selection:

At very large scale (n > 1M), variance is no longer the binding constraint and you can use higher or no clip. At small scale, variance kills you and aggressive clipping pays off even at some bias cost.

The Engineering You Can’t Skip: Logging Propensities

Every estimator above assumes that propensities are available and correct. In practice this assumption is the most commonly violated one — not through negligence, but because production systems are rarely designed with counterfactual evaluation in mind from the start. All of this assumes you have correct propensities. Most teams do not, for several reasons:

Epsilon-greedy (ε-greedy) exploration is easy to log wrong. A system using ε-greedy with ε=0.1 chooses a random action 10% of the time. If you log the greedy action’s score as the propensity (common mistake), every IPS estimate will be wrong.

Stochastic ensembles. If your system samples from a distribution of models or uses beam search with temperature, the propensity requires marginalising over the ensemble — expensive and often approximated.

Non-logged propensities. Legacy systems often do not log any probabilistic information. If you want to retrofit counterfactual evaluation, you need to either reconstruct propensities from logs (fragile) or start an exploration experiment to generate clean propensity data.

The practical fix: When you build a new ranking system, log π₀(a|x) explicitly for every served impression, alongside the action and reward. It is one float per impression. The cost is negligible; the payoff is every counterfactual analysis you will ever want to run.

# What to write to your impression log

from dataclasses import dataclass, asdict
from typing import Optional
import json

@dataclass
class ImpressionLog:
    request_id: str
    user_id: str
    item_id: str
    position: int
    propensity: float          # π₀(action | context) — DO NOT OMIT
    reward: Optional[int]      # filled in asynchronously when click/conversion observed
    model_version: str

def log_impression(log: ImpressionLog, sink) -> None:
    """Write impression to log sink with propensity always included."""
    record = asdict(log)
    sink.write(json.dumps(record) + "\n")

When Counterfactual Evaluation Lies

Even with correct propensities and calibrated weights, there are structural conditions under which all four estimators will give confidently wrong answers. Understanding these failure modes is what separates offline evaluation you can trust from offline evaluation that misleads A/B test decisions.

Covariate shift. Your logs come from the distribution of contexts that π₀ actually served. If π_new would change what contexts appear — for example, by promoting a new category that attracts different users — the logged distribution does not represent what π_new would face in deployment.

Reward censoring. If your reward is “purchased within 7 days” and you evaluate on logs from last week, recent impressions have truncated reward windows. IPS will underestimate the value of recent actions.

Violation of overlap. IPS requires that every action with π_new(a|x) > 0 also has π₀(a|x) > 0. If your logging policy never showed certain items (zero propensity), no amount of importance weighting can estimate their reward. In practice, every large catalogue has a long tail of never-logged items.

Reward model extrapolation. DM predicts rewards for (context, action) pairs not in the logs. If the reward model extrapolates poorly to unseen combinations, DM will be systematically wrong.

The most dangerous failure mode is one you cannot detect from the data itself: when the evaluation looks plausible but the logged distribution is structurally misaligned with the deployment distribution. This is why offline evaluation is a complement to A/B testing, not a replacement. Use counterfactual estimates to filter candidate policies and prioritise what to A/B test — not to skip the test entirely.

Production Checklist

Before trusting any counterfactual estimate in a production decision:

1. Verify propensity coverage. What fraction of (context, action) pairs in the target policy have non-zero propensity in the logs? If it’s below 80%, the estimate is unreliable for the uncovered tail.
2. Run a bias check. Split the logs in half. Use the first half as “target policy” logs and evaluate on the second half with known ground truth. If your estimator can’t recover the known value within 5%, it won’t be accurate on the real target either.
3. Compare IPS, SNIPS, and DR. Large disagreement between them is a signal that assumptions are violated. Agreement doesn’t prove correctness, but disagreement proves something is wrong.
4. Check the effective sample size (ESS). ESS = (∑ w_i)² / ∑ w_i² — values below 10% of n indicate near-degenerate weighting.
5. Run sensitivity to clip threshold. If the estimate is unstable across clip=5/10/20, the high-weight observations are controlling the result.

Summary

Here are summary of estimator we discussed:

Three rules for production use:

1. Log propensities for every served impression — without them, IPS is meaningless.
2. Clip importance weights — unbounded IPS is rarely worth its variance.
3. Use DR as the default — it survives one component being wrong; IPS and DM do not.

Counterfactual evaluation does not replace A/B tests. It makes them cheaper by letting you screen out bad policies before they ever touch production traffic.

Thank you for reading my post, and I hope it was useful for you. If you enjoyed the article and would like to show your support, please consider taking the following actions:

👏 Give the story a round of applause (clap) to help it gain visibility.

📖 Follow me on Medium to access more of the content on my profile. Follow Now

🔔 Subscribe to the newsletter to not miss my latest posts: Subscribe Now

🛎 Connect with me on LinkedIn for updates.