基础：线性回归
考虑一个经典线性高斯模型：
y = a x + u y=ax+u y=ax+u
其中U服从标准高斯分布，a是回归系数，那么回归的目的是找到一个a，使得x与u独立，即
c o v ( y ? a x , x ) = 0 ? c o v ( y , x ) ? a c o v ( x , x ) = 0 ? a = c o v ( y , x ) / c o v ( x , x ) \begin{aligned} & cov( y-ax,x) =0\\ \Longrightarrow & cov( y,x) -acov( x,x) =0\\ \Longrightarrow & a=cov( y,x) /cov( x,x) \end{aligned} ???cov(y?ax,x)=0cov(y,x)?acov(x,x)=0a=cov(y,x)/cov(x,x)?
这些知识在后面会反复用到。
Double/debiased machine learning
我们首先考虑这个因果模型(Partially Linear Regression)：
Y = D θ 0 + g 0 ( X ) + U , E [ U ∣ X , D ] = 0 , D = m 0 ( X ) + V , E [ V ∣ X ] = 0 , \begin{array}{ c l } Y=D\theta _{0} +g_{0} (X)+U, & \mathrm{E} [U\mid X,D]=0,\\ D=m_{0} (X)+V, & \mathrm{E} [V\mid X]=0, \end{array} Y=Dθ0?+g0?(X)+U,D=m0?(X)+V,?E[U∣X,D]=0,E[V∣X]=0,?
其中 Y Y Y 是outcome, D D D 是policy/treatment,且X是个高维的变量
X = ( X 1 , … , X p ) X=( X_{1} ,\dotsc ,X_{p}) X=(X1?,…,Xp?)
在这里，我们最关心的是
θ
0
\displaystyle \theta _{0}
θ0?，因为当我们给定某个X的时候，
θ
0
\displaystyle \theta _{0}
θ0?就表示了在这个X的群体中，D这个treatment的因果效应。那么估计
θ
0
\displaystyle \theta _{0}
θ0?最简单的方法就是，用机器学习去估计，这里我们先将数据随机分成两份，分别是
I
,
I
c
\displaystyle I,I^{c}
I,Ic，不妨假设
g
^
0
\displaystyle \hat{g}_{0}
g^?0?是通过ML估计的函数，于是给定
g
^
0
\displaystyle \hat{g}_{0}
g^?0?，
θ
^
0
\displaystyle \hat{\theta }_{0}
θ^0?可以用线性回归得到
θ
^
0
=
c
o
v
(
D
,
Y
?
g
^
0
(
X
)
)
v
a
r
(
D
)
=
1
n
∑
i
∈
I
D
i
(
Y
i
?
g
^
0
(
X
i
)
)
1
n
∑
i
∈
I
D
i
2
(1)
\hat{\theta }_{0} =\frac{cov( D,Y-\hat{g}_{0} (X))}{var( D)} =\frac{\frac{1}{n}\sum _{i\in I} D_{i}( Y_{i} -\hat{g}_{0} (X_{i} ))}{\frac{1}{n}\sum _{i\in I} D^{2}_{i}} \tag{1}
θ^0?=var(D)cov(D,Y?g^?0?(X))?=n1?∑i∈I?Di2?n1?∑i∈I?Di?(Yi??g^?0?(Xi?))?(1)
这样的估计量其实很容易有bias，主要原因的它非常依赖于
g
^
0
\displaystyle \hat{g}_{0}
g^?0?的准确度，万一它有一点偏差就会产生很大的影响，我们来分析一下：
n ( θ ^ 0 ? θ 0 ) = n 1 n ∑ i ∈ I D i ( Y i ? g ^ 0 ( X i ) ) 1 n ∑ i ∈ I D i 2 ? ( n 1 n ∑ i ∈ I D i ( Y i ? g 0 ( X i ) ) 1 n ∑ i ∈ I D i 2 ? n 1 n ∑ i ∈ I D i U i 1 n ∑ i ∈ I D i 2 ) = ( 1 n ∑ i ∈ I D i 2 ) ? 1 1 n ∑ i ∈ I D i U i ? : = a + ( 1 n ∑ i ∈ I D i 2 ) ? 1 1 n ∑ i ∈ I D i ( g 0 ( X i ) ? g ^ 0 ( X i ) ) ? : = b . \begin{aligned} \sqrt{n}(\hat{\theta }_{0} -\theta _{0}) & =\sqrt{n}\frac{\frac{1}{n}\sum _{i\in I} D_{i}( Y_{i} -\hat{g}_{0} (X_{i} ))}{\frac{1}{n}\sum _{i\in I} D^{2}_{i}} -\left(\sqrt{n}\frac{\frac{1}{n}\sum _{i\in I} D_{i}( Y_{i} -g_{0} (X_{i} ))}{\frac{1}{n}\sum _{i\in I} D^{2}_{i}} -\sqrt{n}\frac{\frac{1}{n}\sum _{i\in I} D_{i} U_{i}}{\frac{1}{n}\sum _{i\in I} D^{2}_{i}}\right)\\ & =\underbrace{\left(\frac{1}{n}\sum _{i\in I} D^{2}_{i}\right)^{-1}\frac{1}{\sqrt{n}}\sum _{i\in I} D_{i} U_{i}}_{:=a} +\underbrace{\left(\frac{1}{n}\sum _{i\in I} D^{2}_{i}\right)^{-1}\frac{1}{\sqrt{n}}\sum _{i\in I} D_{i}( g_{0}( X_{i}) -\hat{g}_{0}( X_{i}))}_{:=b} . \end{aligned} n?(θ^0??θ0?)?=n?n1?∑i∈I?Di2?n1?∑i∈I?Di?(Yi??g^?0?(Xi?))??(n?n1?∑i∈I?Di2?n1?∑i∈I?Di?(Yi??g0?(Xi?))??n?n1?∑i∈I?Di2?n1?∑i∈I?Di?Ui??)=:=a (n1?i∈I∑?Di2?)?1n?1?i∈I∑?Di?Ui???+:=b (n1?i∈I∑?Di2?)?1n?1?i∈I∑?Di?(g0?(Xi?)?g^?0?(Xi?))??.?
误差分为两部分，一部分是U和D的独立性，如果不独立的话会形成误差，即
θ 0 = c o v ( D , Y ? g 0 ( X ) ) v a r ( D ) ? c o v ( D , U ) v a r ( D ) ? ≠ 0 \theta _{0} =\frac{cov( D,Y-g_{0} (X))}{var( D)} - \underbrace{\frac{cov( D,U)}{var( D)}}_{\ne0} θ0?=var(D)cov(D,Y?g0?(X))???=0 var(D)cov(D,U)???
当两者独立的话，第一项的收敛其实很快，服从 a ～ N ( 0 , Σ ￣ ) a\sim N(0,\overline{\Sigma } ) a～N(0,Σ). 问题在 b b b 项，如果X维度很高，那么我们可能会在拟合g时加点正则项（如lasso），就会形成正则化误差，使得g有错误的估计，于是导致b项发散，
b = ( E [ D i 2 ] ) ? 1 1 n ∑ i ∈ I m 0 ( X i ) ( g 0 ( X i ) ? g ^ 0 ( X i ) ) + o P ( 1 ) b=\left(E\left[ D^{2}_{i}\right]\right)^{-1}\frac{1}{\sqrt{n}}\sum _{i\in I} m_{0}( X_{i})( g_{0}( X_{i}) -\hat{g}_{0}( X_{i})) +o_{P} (1) b=(E[Di2?])?1n?1?i∈I∑?m0?(Xi?)(g0?(Xi?)?g^?0?(Xi?))+oP?(1)
b项可以进一步写为上式，可以看到主要问题出在 m 0 ( X i ) ( g 0 ( X i ) ? g ^ 0 ( X i ) ) \displaystyle m_{0}( X_{i})( g_{0}( X_{i}) -\hat{g}_{0}( X_{i})) m0?(Xi?)(g0?(Xi?)?g^?0?(Xi?))这个地方， m 0 ( X i ) \displaystyle m_{0}( X_{i}) m0?(Xi?)的大小将会极大的影响误差的大小，而 m 0 m_0 m0?是数据的性质，无法改变，就会导致我们的估计很不稳定（数据会决定bias的大小）。
有没有可能消掉这一项呢？我们先分析下 m 0 ( X i ) \displaystyle m_{0}( X_{i}) m0?(Xi?)出现的最主要的原因是我们用D来回归，然而D包含了X的信息，于是才会有 m 0 ( X i ) \displaystyle m_{0}( X_{i}) m0?(Xi?)这一项出现。那有没有可能不用D来回归呢？一个重要的观察是，V其实可以看做工具变量，
于是D和Y之前的因果效应可以被以下式子估计：
θ 0 = c o v ( V , Y ? g 0 ( X ) ) c o v ( V , D ) = c o v ( D ? m 0 ( X ) , Y ? g 0 ( X ) ) c o v ( D ? m 0 ( X ) , D ) = c o v ( D ? m 0 ( X ) , D θ 0 + U ) c o v ( D ? m 0 ( X ) , D ) \theta _{0} =\frac{cov( V,Y-g_{0}( X))}{cov( V,D)} =\frac{cov( D-m_{0}( X) ,Y-g_{0}( X))}{cov( D-m_{0}( X) ,D)} =\frac{cov( D-m_{0}( X) ,D\theta _{0} +U)}{cov( D-m_{0}( X) ,D)} θ0?=cov(V,D)cov(V,Y?g0?(X))?=cov(D?m0?(X),D)cov(D?m0?(X),Y?g0?(X))?=cov(D?m0?(X),D)cov(D?m0?(X),Dθ0?+U)?
因此，为了求V，我们可以求：
V ^ = D ? m ^ 0 ( X ) \hat{V} =D-\hat{m}_{0}( X) V^=D?m^0?(X)
这个
m
^
0
(
X
)
\displaystyle \hat{m}_{0}( X)
m^0?(X)可以用X对D回归来得到，于是得到了新的一种估计
θ
?
0
=
1
n
∑
i
∈
I
V
^
i
(
Y
i
?
g
^
0
(
X
i
)
)
1
n
∑
i
∈
I
V
i
^
D
i
(2)
\breve{\theta }_{0} =\frac{\frac{1}{n}\sum _{i\in I}\hat{V}_{i}( Y_{i} -\hat{g}_{0} (X_{i} ))}{\frac{1}{n}\sum _{i\in I}\widehat{V_{i}} D_{i}} \tag{2}
θ?0?=n1?∑i∈I?Vi?
?Di?n1?∑i∈I?V^i?(Yi??g^?0?(Xi?))?(2)
直观来看，因为V不包含X的信息，所以在上面b项中
m
0
(
X
i
)
\displaystyle m_{0}( X_{i})
m0?(Xi?)带来的误差可以被消去，事实上，在这个新的估计量下，这个新的b项将变为
b ? = ( E [ D i 2 ] ) ? 1 1 n ∑ i ∈ I ( m ^ 0 ( X i ) ? m 0 ( X i ) ) ( g 0 ( X i ) ? g ^ 0 ( X i ) ) b^{*} =\left(E\left[ D^{2}_{i}\right]\right)^{-1}\frac{1}{\sqrt{n}}\sum _{i\in I}(\hat{m}_{0}( X_{i}) -m_{0}( X_{i}))( g_{0}( X_{i}) -\hat{g}_{0}( X_{i})) b?=(E[Di2?])?1n?1?i∈I∑?(m^0?(Xi?)?m0?(Xi?))(g0?(Xi?)?g^?0?(Xi?))
显然，这将能得到一个更为robust估计，因为现在bias只取决于回归误差了。
PS：从 θ 0 = c o v ( D ? m 0 ( X ) , D θ 0 + U ) c o v ( D ? m 0 ( X ) , D ) \displaystyle \theta _{0} =\frac{cov( D-m_{0}( X) ,D\theta _{0} +U)}{cov( D-m_{0}( X) ,D)} θ0?=cov(D?m0?(X),D)cov(D?m0?(X),Dθ0?+U)?这个式子来看，其实我们完全可以将下面的D也换成 θ 0 = c o v ( D ? m 0 ( X ) , D θ 0 + U ) c o v ( D ? m 0 ( X ) , D ? m 0 ( X ) ) \displaystyle \theta _{0} =\frac{cov( D-m_{0}( X) ,D\theta _{0} +U)}{cov( D-m_{0}( X) ,D-m_{0}( X))} θ0?=cov(D?m0?(X),D?m0?(X))cov(D?m0?(X),Dθ0?+U)?，这并不影响我们的结果，这样的估计量最早由Robinson (1988) [1] 提出.
Neyman orthogonality and moment conditions.
上文只是介绍了一个Partially Linear Regression构造robust估计量的方法，那这套方法能不能被抽象出来，使其适用于更多的场景？
那我们先总结一些性质，
- 我们有一个目标参数 θ 0 \displaystyle \theta _{0} θ0?，然后有其他不太关心的但可以用各种机器学习方法得到的回归函数， g , m \displaystyle g,m g,m
- 估计 θ 0 \displaystyle \theta _{0} θ0?有多种不同的方法，有的好有的坏，而坏的原因在于某个用机器学习拟合的 g ^ \displaystyle \hat{g} g^?如果出现误差，则误差项的均值将不再为0，从而不收敛。
- 为此，我们需要构造一种“正交”的估计量，使得即使g或者m出现误差也能保证误差足够小。
针对第一个性质，我们可以将公式1的估计量写成：
1 n ∑ i ∈ I φ ( W ; θ ^ 0 , g ^ 0 ) = 0 \frac{1}{n}\sum _{i\in I} \varphi ( W;\hat{\theta }_{0} ,\hat{g}_{0}) =0 n1?i∈I∑?φ(W;θ^0?,g^?0?)=0
在上面的例子，就是线性回归中，找到噪声与x协方差独立，即
c o v ( Y ? θ D ? g ( X ) , D ) = 0 ? θ = c o v ( Y ? g ( X ) , D ) / v a r ( D ) cov( Y-\theta D-g( X) ,D) =0\\ \Longrightarrow \theta =cov( Y-g( X) ,D) /var( D) cov(Y?θD?g(X),D)=0?θ=cov(Y?g(X),D)/var(D)
因此
φ
(
W
;
θ
,
g
)
=
(
Y
?
θ
D
?
g
(
X
)
)
D
\varphi ( W;\theta ,g) =( Y-\theta D-g( X)) D
φ(W;θ,g)=(Y?θD?g(X))D
相当于
c
o
v
(
Y
?
θ
D
?
g
(
X
)
,
D
)
=
0
\displaystyle cov( Y-\theta D-g( X) ,D) =0
cov(Y?θD?g(X),D)=0.
然而在这个估计量的问题是，当g不太准的时候，将会有较大的误差，这个误差理解为g一些微小的扰动而导致
φ
\varphi
φ发生较大的变化。那么怎么衡量这一变化呢？其实就是Gateaux dervative导数：
?
g
E
[
φ
(
W
;
θ
0
,
g
0
)
]
[
g
?
g
0
]
≠
0
\partial _{g} E[ \varphi ( W;\theta _{0} ,g_{0})][ g-g_{0}] \neq 0
?g?E[φ(W;θ0?,g0?)][g?g0?]?=0
不等于0。这个Gateaux dervative的定义是
? g E [ φ ( W ; θ 0 , g 0 ) ] = lim ? r → 0 + E [ φ ( W ; θ 0 , g 0 + r ( g ? g 0 ) ) ] ? E [ φ ( W ; θ 0 , g 0 ) ] r \partial _{g} E[ \varphi ( W;\theta _{0} ,g_{0})] =\lim _{r\rightarrow 0^{+}}\frac{E[ \varphi ( W;\theta _{0} ,g_{0} +r( g-g_{0}))] -E[ \varphi ( W;\theta _{0} ,g_{0})]}{r} ?g?E[φ(W;θ0?,g0?)]=r→0+lim?rE[φ(W;θ0?,g0?+r(g?g0?))]?E[φ(W;θ0?,g0?)]?
直观上，它刻画了 φ \displaystyle \varphi φ对 g \displaystyle g g变化的敏感度，那如果这个导数等于0，则意味着我们找到了一个更为robust的估计量，事实上，我们上文中的第二种方法正是这样的估计量，公式2可以写成：
1 n ∑ i ∈ I ψ ( W ; θ ? 0 , η ^ 0 ) = 0 \frac{1}{n}\sum _{i\in I} \psi ( W;\breve{\theta }_{0} ,\hat{\eta }_{0}) =0 n1?i∈I∑?ψ(W;θ?0?,η^?0?)=0
其中 η ^ 0 = ( m ^ 0 , g ^ 0 ) \displaystyle \hat{\eta }_{0} =(\hat{m}_{0} ,\hat{g}_{0}) η^?0?=(m^0?,g^?0?)，被称为nuisance parameter，并且
ψ ( W ; θ , η ) = ( Y ? θ D ? g ( x ) ) ( D ? m ( x ) ) \psi ( W;\theta ,\eta ) =( Y-\theta D-g( x))( D-m( x)) ψ(W;θ,η)=(Y?θD?g(x))(D?m(x))
相当于
c
o
v
(
Y
?
θ
D
?
g
(
x
)
,
D
?
m
(
x
)
)
\displaystyle cov( Y-\theta D-g( x) ,D-m( x))
cov(Y?θD?g(x),D?m(x))。可以证明这个估计量的Gateaux dervative为0：
?
η
E
[
φ
(
W
;
θ
0
,
η
0
)
]
[
η
?
η
0
]
=
0
(3)
\partial _{\eta } E[ \varphi ( W;\theta _{0} ,\eta _{0})][ \eta -\eta _{0}] =0 \tag{3}
?η?E[φ(W;θ0?,η0?)][η?η0?]=0(3)
我们把这个性质称为Neyman orthogonality.
那有没有一种通用的方法来构造满足Neyman orthogonality性质的这样一个 ψ \displaystyle \psi ψ?
不妨设
max ? θ , β E p [ l ( W ; θ , β ) ] \max_{\theta ,\beta } E_{p}[ l( W;\theta ,\beta )] θ,βmax?Ep?[l(W;θ,β)]
我们的目标可以通过最大化某个函数得到（通常为对数似然度），于是在极值下的导数为0：
E p [ ? θ l ( W ; θ 0 , β 0 ) ] = 0 , E p [ ? β l ( W ; θ 0 , β 0 ) ] = 0 E_{p}[ \partial _{\theta } l( W;\theta _{0} ,\beta _{0})] =0,E_{p}[ \partial _{\beta } l( W;\theta _{0} ,\beta _{0})] =0 Ep?[?θ?l(W;θ0?,β0?)]=0,Ep?[?β?l(W;θ0?,β0?)]=0
因此，自然地我们可以令 ? ( W ; θ , β ) = ? θ l ( W ; θ , β ) \displaystyle \phi ( W;\theta ,\beta ) =\partial _{\theta } l( W;\theta ,\beta ) ?(W;θ,β)=?θ?l(W;θ,β)，然而这样的函数不一定满足Neyman orthogonality，所以我们可以自己构造一个估计量
ψ ( W ; θ , β ) = ? θ l ( W ; θ , β ) ? μ ? β l ( W ; θ , β ) \psi ( W;\theta ,\beta ) =\partial _{\theta } l( W;\theta ,\beta ) -\mu \partial _{\beta } l( W;\theta ,\beta ) ψ(W;θ,β)=?θ?l(W;θ,β)?μ?β?l(W;θ,β)
只要能找到一个合适的 μ \displaystyle \mu μ，使得 ψ \displaystyle \psi ψ，满足Neyman orthogonality就可以了，即
? β ψ ( W ; θ , β ) = ? θ β l ( W ; θ , β ) ? μ ? β β l ( W ; θ , β ) = 0 \partial _{\beta } \psi ( W;\theta ,\beta ) =\partial _{\theta \beta } l( W;\theta ,\beta ) -\mu \partial _{\beta \beta } l( W;\theta ,\beta ) =0 ?β?ψ(W;θ,β)=?θβ?l(W;θ,β)?μ?ββ?l(W;θ,β)=0
我们可以得出解析解， μ = J θ β J β β ? 1 \displaystyle \mu =J_{\theta \beta } J^{-1}_{\beta \beta } μ=Jθβ?Jββ?1?
参考资料
[1] Robinson, Peter M. “Root-N-consistent semiparametric regression.” Econometrica: Journal of the Econometric Society (1988): 931-954.
[2] Double/debiased machine learning for treatment and structural parameters
[3] https://observablehq.com/@herbps10/one-step-estimators-and-pathwise-derivatives
[4] Visually Communicating and Teaching Intuition for Influence Functions
[5] Approximation Theorems of Mathematical Statistics.
[6] https://vanderlaan-lab.org/2019/12/24/cv-tmle-and-double-machine-learning/
[7] https://academic.oup.com/biostatistics/article/21/2/353/5631845