Good question -- this is a really subtle distinction, and I think it is mostly practical. If you read on a few pages, you'll see that both methods have the same estimand (what we're trying to predict/model) -- namely $P(Y^a = 1)$, which is obtained through (informally) "adjustment" for confounders $L$. Thus, in a perfect world, where we know all of the conditional/joint/marginal distributions in the above arithmetic with certainty, the distinction is moot. But mathematical equivalency need not imply practical equivalence.
Subjectively speaking, I like to think that there are three stages in the causal inference pipeline:
- Identify the estimand of interest (i.e., treatment effect)
- Verify that the estimand is identifiable (e.g., via backdoor criterion)
- Fit model(s) to estimate the estimand
So, Step 1 seems to be the same for IPW vs. standardization as you've pointed out, and Step 2 is orthogonal to your question (we'll assume it's true). So the distinction comes from Step 3: modeling. Let's dive in.
The following explanation is adapted from What If (Hernan & Robins, 2020), Chapter 13.4, with some tweaks based on how I think of this issue.
Background: what would our estimators actually look like? The difference is in how we would actually build estimators in each case (i.e., modeling choices). Consider inverse probability weighting (IPW). We generally need to fit some propensity model $f(A\mid L)$ (e.g., logistic regression) to predict treatment $A$ from covariates $L$. The $f(A\mid L)$ are then used as weights for fitting an outcome prediction model $g(Y \mid A)$ (e.g., a regression). So our final estimator is something like:
$$\hat{P}(Y^a = 1) = \frac{1}{N} \sum_{i}^N g(a_i) \left(\frac{\mathbf{1}[a_i = 1]}{f(\ell_i)} + \frac{\mathbf{1}[a_i = 0]}{1 - f(\ell_i)}\right)$$
with some abuse of notation (treating $g(Y \mid A)$ as a function of $a_i$ that outputs estimates of outcome $Y$ and likewise for $f$). For $Y \approxeq g(a_i)$ I believe this equivalent to the more standard form of writing such estimators: $\mathbb{E}\left[\frac{\mathbf{1}[A = a]Y}{f(A \mid L)}\right]$ (e.g., Horvitz-Thomposon).
Then, consider standardization. For each value of $L= \ell$ in our data, we would fit a model $g_{L=\ell}(Y \mid A)$, and output an estimate
$$\hat{P}(Y^a = 1) = \frac{1}{N} \sum_{i}^N g_{L=\ell_i}(a_i),$$
with some abuse of notation ($g_{L=\ell}$ here is treated as a function that takes in $a_i$ and outputs outcome estimates of $Y$). In other words, for each $(a_i, \ell_i)$, we pick the model $g_{L=\ell_i}$ that we fit for that particular strata of $L$, and use it to output an outcome estimate. We then average all estimates. Convince yourself that this is equivalent to (at the single-example level) $P(Y = 1 \mid A = a, L = \ell) \cdot P(L = \ell)$.
OK, what's the punchline? Let's list the models that we fit for each method.
- IPW: We fit an estimator $g$ to predict outcome $Y$ from treatment $A$, and another propensity model $f$ to predict treatment $A$ from covariates $L$.
- Standardization: We fit multiple outcome estimators $g$ to predict outcome $Y$ from treatment $A$ for each strata of $L$.
So, we had to make different modeling choices to execute each method in practice. This matters because it's doubtful that any model we choose will be perfectly specified, so we're likely to incur differing amounts of misspecification bias between the two approaches. But it might be a good robustness check to use both to check if they yield "similar" estimates.
As for the statistical properties of IPW vs. standardization-based estimators, and reasons to prefer either from that angle (e.g., asymptotic analysis) -- unfortunately I wouldn't know, and I suspect that has to do more with model choice.
So, the takeaway is that there can be more than one way to estimate a quantity of interest.
In practice, more modern approaches (i.e., doubly-robust estimation) are generally used as baselines for causal effect estimation to combat such misspecification bias, but this is orthogonal to your question. A good place to start here is the augmented inverse probability weighting estimator, where we essentially use the propensity score (i.e., $P(A \mid L)$ estimates) as an extra covariate -- multiple methods follow a similar template (Bang & Robins, 2005; Scharfstein, 1999 -- I'm sure I'm forgetting a few).This explanation simply illustrates that there is sometimes more than one way to construct an estimator for a quantity.
