Is DML a doubly robust estimator?

[james-belmonte]

I am curious if DML estimators like the one presented here are doubly robust, I read somewhere that they are, but it was on the stack overflow stats page and the poster didnt explain why, just asserted it, so I dont want to just blindly trust that.

Can you point to a paper that works this out, or a section of a paper that shows this?

[nmriabov]

DML estimators are doubly-robust. The mathematical construction of the proof of double robustness rests on the fact that all such estimators consist of two components (usually, Pr(treatment | ...) and E[outcome | ...]), and that the bias of the estimate of the ATE/CATE is the product of the bias of these two parts. Usually, the requirement is that this product goes to 0 at a rate of $\sqrt(n)$

The reason that DML, specifically, is doubly-robust, is because of the cross-fitting. This step kills the bias of the outcome model.

More formally, here's a paper. On page 10, there's a passage addressing DML itself:

The need for doubly-robust estimators with cross-fitting when using data-adaptive machine learning for nuisance
function estimation arises from two terms in the Von Mises expansion of the estimator [1]. The first term, which is
described by an empirical process term in the expansion, can be controlled by either restricting the complexity of
the nuisance models (e.g., by requiring them to be in the Donsker class) or through cross-fitting. Because it can be
difficult or impossible to verify that a given machine learning method is in the Donsker class, cross-fitting provides
a simple and attractive alternative. The second term is the second-order remainder, and it converges to zero as
the sample size increases. For valid inference, it is desirable for this remainder term to converge as a function of
O(n^{−1/2}), referred to as as root-n convergence. Convergence rates are not a computational issue, but rather a feature
of the estimator itself. Unfortunately, data-adaptive algorithms often have slower convergence rates as a result of
their flexibility. However, because the second-order remainder term of doubly-robust estimators is the product of the
approximation errors of the treatment and outcome nuisance models, doubly-robust estimators only require that the
product of the convergence rates for nuisance models be O(n^{−1/2}). To summarize, cross-fitting permits the use of highly
complex nuisance models, while doubly-robust estimators permit the use of slowly converging nuisance models. Used
together, these approaches allow one to use a wide class of data-adaptive machine learning methods to estimate
causal effects.

[james-belmonte]

thanks so much for the response and the linked paper!

[kbattocchi]

Unfortunately, I don't believe this is correct; DML is robust with respect to bias in the Y model, but not in the T model, though it is true that if they both converge at the relatively slow rate n^(-1/4) then the estimate of the effect will be at the product of those rates, n^(-1/2).

You can see this with a simple example:

from econml.dml import LinearDML
from sklearn.linear_model import LinearRegression
import numpy as np

np.random.seed(1)
n=100000
W = np.random.normal(size=(n,5))
theta = 2
T = W @ np.random.normal(size=(5,)) + np.random.normal(size=(n,))
y = W @ np.random.normal(size=(5,)) + T * theta + np.random.normal(size=(n,))

class BiasedEstimator:
  def __init__(self, estimator, bias):
    self.estimator = estimator
    self.bias = bias

  def fit(self, X, y):
    self.estimator.fit(X, y)
    return self

  def predict(self, X):
    return self.estimator.predict(X) + self.bias

bt, by, ok = [], [], []

bt_dml = LinearDML(model_y=LinearRegression(), model_t=BiasedEstimator(LinearRegression(), 0.1))
by_dml = LinearDML(model_t=LinearRegression(), model_y=BiasedEstimator(LinearRegression(), 0.1))
ok_dml = LinearDML(model_y=LinearRegression(), model_t=LinearRegression())
ns = np.linspace(10000,100000,10,dtype=int)
for n in ns:
  print(n)
  bt.append(bt_dml.fit(y[:n], T[:n], W=W[:n,:]).const_marginal_effect_inference())
  by.append(by_dml.fit(y[:n], T[:n], W=W[:n,:]).const_marginal_effect_inference())
  ok.append(ok_dml.fit(y[:n], T[:n], W=W[:n,:]).const_marginal_effect_inference())

import matplotlib.pyplot as plt

plt.figure()

for lst in [bt,by,ok]:
  intervals = np.array([inf.conf_int() for inf in lst]).squeeze(2)
  plt.plot(ns, [inf.point_estimate for inf in lst])
  plt.fill_between(ns, intervals[:,0], intervals[:,1], alpha=0.2)
plt.legend(['Biased T', 'Biased Y', 'Unbiased'])

[james-belmonte]

@kbattocchi thanks for following up

Your example makes sense, can you give an intuitive explanation why DML would be robust in Y but not T?

[kbattocchi]

Sorry, I should have been more precise in my answer: DML can refer to a number of different estimators; what I actually meant is that the specific DML algorithm we currently implement as part of our package is not doubly robust. If you look at one of the major papers on DML, they propose two different score functions (equations 4.2 and 4.3) that can be used for DML estimators for the treatment effects problem we're trying to solve here. We use the second of these, which is equivalent to regressing the outcome residuals on the treatment residuals in the second stage.

I think one way to see why this is not robust to the treatment model is to look at the moment equation that we're trying to solve:

(T-f(X))ᵀ(Y-h(X)-(T-f(X))θ)=0

where f is our model for T and h is our model for Y. Then, if we perturb the Y model h to h+h' on the left hand side, the score will decrease by (T-f(X))ᵀh', which if the T model is accurate should be roughly zero, and so the estimated value for θ should be stable. On the other hand, if we perturb the T model f to f+f', the score will change by (Y-h(X)-2(T-f(X))θ)ᵀf', which will not necessarily be approximately zero even if the Y model is accurate because of the T-f(X) component, and so the θ we compute will not be stable.

However, the reason that we don't just use the first score function instead is that there the first-stage Y model is no longer a general predictive model, but rather a semi-linear model where T has to enter linearly, which doesn't map well to the paradigm where you can use any sklearn model of your choice.

[james-belmonte]

that is very helpful, thank you!