######  28.2 Observational studies with a multi-valued treatment

Theorem 28.1 is a special case of the following theorem for unconfounded
observational studies with multiple treatment levels (Imai and Van Dyk,, 2004;
Cattaneo,, 2010). Below, I state the general problem and theorem.

Consider an observational study with a multi-valued treatment
$Z\in\\{1,\ldots,K\\}$, covariates $X$, and outcome $Y$. Unit $i$ has $K$
potential outcomes $Y_{i}(1),\ldots,Y_{i}(K)$ corresponding to the $K$
treatment levels. Causal effects can be defined as comparisons of the
potential outcomes. In general, we can define causal effect in terms of
contrasts of the potential outcomes:

$\tau_{c}=\sum_{k=1}^{K}c_{k}E\\{Y(k)\\}$

where $\sum_{k=1}^{K}c_{k}=0$. The canonical choice of the pairwise comparison

$\tau_{k,k^{\prime}}=E\\{Y(k)-Y(k^{\prime})\\}.$

Therefore, the key is to identify and estimate the means of the potential
outcomes $\mu_{k}=E\\{Y(k)\\}$ under the ignorability assumption below based
on IID data of $(Z_{i},X_{i},Y_{i})_{i=1}^{n}$.

###### Assumption 28.2

$Z\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}\\{Y(1),\ldots,Y(K)\\}\mid X$.

Define the generalized propensity score as

$e_{k}(X)=\textup{pr}(Z=k\mid X),$

and define the conditional outcome mean as

$\mu_{k}(X)=E(Y\mid Z=k,X)$

for $k=1,\ldots,K.$ We have the following theorem.

###### Theorem 28.2

Under Assumption 28.2, we have

$\mu_{k}=E\\{\mu_{k}(X)\\}$

or

$\mu_{k}=E\left\\{\frac{I(Z=k)Y}{e_{k}(X)}\right\\}.$

Moreover, based on the working models $e_{k}(X,\alpha)$ and
$\mu_{k}(X,\beta)$, we have the doubly robust formula

$\mu_{k}^{\textup{dr}}=E\\{\mu_{k}(X,\beta)\\}+E\left[\frac{I(Z=k)\\{Y-\mu_{k}(X,\beta)\\}}{e_{k}(X,\alpha)}\right],$

which equals $\mu_{k}$ if either $e_{k}(X,\alpha)=e_{k}(X)$ or
$\mu_{k}(X,\beta)=\mu_{k}(X)$.

Prove Theorem 28.2.

Remark: Theorem 28.1 is a special case of Theorem 28.2 if we view the $(Z,M)$
in Theorem 28.1 as a treatment with four levels. The $\textsc{cde}(m)$ is a
special case of $\tau_{c}$.

######  28.3 cde in the linear outcome model

Show that under Assumption 28.1, if $E(Y\mid
Z,M,X)=\theta_{0}+\theta_{1}Z+\theta_{2}M+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X$,
then

$\textsc{cde}(m)=\theta_{1}$

for all $m$; if $E(Y\mid
Z,M,X)=\theta_{0}+\theta_{1}Z+\theta_{2}M+\theta_{3}ZM+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X$,
then

$\textsc{cde}(m)=\theta_{1}+\theta_{3}m.$

######  28.4 cde in the logit outcome model

Show that for a binary outcome, under Assumption 28.1, if

$\textup{logit}\\{\textup{pr}(Y=1\mid
Z,M,X)\\}=\theta_{0}+\theta_{1}Z+\theta_{2}M+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X,$

then

$\textsc{cde}(m)=E\\{\textup{expit}(\theta_{0}+\theta_{1}+\theta_{2}m+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X)-\textup{expit}(\theta_{0}+\theta_{2}m+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X)\\};$

if

$\textup{logit}\\{\textup{pr}(Y=1\mid
Z,M,X)\\}=\theta_{0}+\theta_{1}Z+\theta_{2}M+\theta_{3}ZM+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X,$

then

$\textsc{cde}(m)=E\\{\textup{expit}(\theta_{0}+\theta_{1}+\theta_{2}m+\theta_{3}m+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X)-\textup{expit}(\theta_{0}+\theta_{2}m+\theta_{4}^{\mkern-1.5mu\mathsf{T}}X)\\}.$

######  28.5 Recommended reading

Nguyen et al., (2021) provided a friendly review of of the topics in Chapters
27 and 29.

###  Chapter 29 Time-Varying Treatment and Confounding

Studies with time-varying treatments are common in biomedical and social
sciences. James Robins championed the research in biostatistics. A classic
example is that HIV patients may take the azidothymidine, an antiretroviral
medication, on and off over time (Robins et al.,, 2000; Hernán et al.,, 2000).
Similar problems also exist in other fields. In education, a classic example
is that students may receive different types of instructions over time (Hong
and Raudenbush,, 2008). In political science, a classic example is that
candidates continuously recalibrate their campaign strategy based on time-
varying polls and opponent actions (Blackwell,, 2013).

Causal inference with a time-varying treatment is not a simple extension of
causal inference with a treatment at a single time point. The main challenge
is time-varying confounding. Even if we assume all time-varying confounders
are observed, we still face statistical challenges in adjusting for those
confounders. One the one hand, we should stratify on these confounders to
adjust for confounding; on the other hand, stratifying on post-treatment
variables will cause bias. Due to these two conflicting goals, causal
inference with time-varying treatments and confounding requires more
sophisticated statistical methods. It is the main topic of this chapter.

To minimize the notational burden, I will use the setting with a treatment at
two time points to convey the most important ideas. Extensions to treatments
at multiple time points can be conceptually straightforward although technical
complexities will arise in finite samples. I will discuss the complexities and
relegate general results to Problems 29.6–29.9.

####  147 Basic setup and sequential ignorability

Start with a treatment at two time points. The temporal order of the variables
with two time points is below:

$X_{0}\rightarrow Z_{1}\rightarrow X_{1}\rightarrow Z_{2}\rightarrow Y$

where

  * •

$X_{0}$ denotes the baseline pre-treatment covariates;

  * •

$Z_{1}$ denotes the treatment at time point 1;

  * •

$X_{1}$ denotes the time-varying covariates between the treatments at time
points 1 and 2;

  * •

$Z_{2}$ denotes the treatment at time point 2;

  * •

$Y$ denotes the outcome.

With binary treatment $(Z_{1},Z_{2})$, each unit has four potential outcomes

$Y(z_{1},z_{2})\text{ for }z_{1},z_{2}=0,1.$

The observed outcome equals

$Y=Y(Z_{1},Z_{2})=\sum_{z_{1}=0,1}\sum_{z_{2}=0,1}1(Z_{1}=z_{1})1(Z_{2}=z_{2})Y(z_{1},z_{2}).$

I will focus on the canonical setting with sequential ignorability, that is,
the treatments are sequentially randomized given the observed history.

######  Assumption 29.1 (sequential ignorability)

(1) $Z_{1}$ is randomized given $X_{0}$:

$Z_{1}\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}Y(z_{1},z_{2})\mid X_{0}\text{ for }z_{1},z_{2}=0,1.$

(2) $Z_{2}$ is randomized given $(Z_{1},X_{1},X_{0})$:

$Z_{2}\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}Y(z_{1},z_{2})\mid(Z_{1},X_{1},X_{0})\text{ for
}z_{1},z_{2}=0,1.$

Figure 42 is a simple causal diagram corresponding to Assumption 29.1, which
does not contains any unmeasured confounding.

$\textstyle{Z_{1}\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{X_{1}\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{Z_{2}\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{Y}$
Figure 42: Without unmeasured confounding $U$ between $X_{1}$ and $Y$. The
causal diagram conditions on the pre-treatment covariates $X_{0}$.

Figure 43 is a more complex causal diagram corresponding to Assumption 29.1.
Sequential ignorability rules out only the confounding between the treatment
$(Z_{1},Z_{2})$ and the outcome $Y$, but allows for unmeasured confounding
between the time-varying covariate $X_{1}$ and the outcome $Y$. The possible
existence of $U$ causes many subtle issues even under sequential ignorability.

---
$\textstyle{Z_{1}\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{X_{1}\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{Z_{2}\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{Y}$$\textstyle{U\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$
Figure 43: With unmeasured confounding between $X_{1}$ and $Y$. The causal
diagram conditions on the pre-treatment covariates $X_{0}$.

####  148 g-formula and outcome modeling

Recall the outcome-based identification formula with a treatment at a single
time point:

$E\\{Y(z)\\}=E\\{E(Y\mid Z=z,X)\\}.$

With discrete $X$, it reduces to

$E\\{Y(z)\\}=\sum_{x}E(Y\mid Z=z,X=x)\textup{pr}(X=x);$

with continuous $X$, it reduces to

$E\\{Y(z)\\}=\int E(Y\mid Z=z,X=x)f_{X}(x)\text{d}x.$

The following result extends it to the setting with a treatment at two time
points.

###### Theorem 29.1

Under Assumption 29.1,

$E\\{Y(z_{1},z_{2})\\}=E\Big{[}E\\{E(Y\mid z_{2},z_{1},X_{1},X_{0})\mid
z_{1},X_{0}\\}\Big{]}.$ (126)

In Theorem 29.1, I simplify the notation ‘$Z_{2}=z_{2}$” to “$z_{2}$” for
simplicity. To void complex formulas in this Chapter, I will use the lower
case letter to represent the event that the random variable takes the
corresponding value. With discrete $X_{0}$ and $X_{1}$, the identification
formula (126) reduces to

$E\\{Y(z_{1},z_{2})\\}=\sum_{x_{0}}\sum_{x_{1}}E(Y\mid
z_{2},z_{1},x_{1},x_{0})\textup{pr}(x_{1}\mid z_{1},x_{0})\textup{pr}(x_{0});$
(127)

with continuous $X_{0}$ and $X_{1}$, the identification formula (126) reduces
to

$E\\{Y(z_{1},z_{2})\\}=\int\int E(Y\mid z_{2},z_{1},x_{1},x_{0})f(x_{1}\mid
z_{1},x_{0})f(x_{0})\text{d}x_{1}\text{d}x_{0}.$ (128)

Compare (127) with the formula based on the law of total probability to gain
more insights:

$\displaystyle E(Y)$ $\displaystyle=$
$\displaystyle\sum_{x_{0}}\sum_{z_{1}}\sum_{x_{1}}\sum_{z_{2}}E(Y\mid
z_{2},z_{1},x_{1},x_{0})$ (129) $\displaystyle\qquad\textup{pr}(z_{1}\mid
z_{1},x_{1},x_{0})\textup{pr}(x_{1}\mid z_{1},x_{0})\textup{pr}(z_{1}\mid
x_{0})\textup{pr}(x_{0}).$

Erasing the probabilities of $z_{2}$ and $z_{1}$ in (129), we can obtain the
formula (128). This is intuitive because the potential outcome
$Y(z_{1},z_{2})$ has the meaning of fixing $Z_{1}$ and $Z_{2}$ at $z_{1}$ and
$z_{2}$, respectively.

Robins called (127) and (128) the g-formulas. Now I will prove Theorem 29.1.

Proof of Theorem 29.1: By the tower property,

$E\\{Y(z_{1},z_{2})\\}=E\Big{[}E\\{Y(z_{1},z_{2})\mid X_{0}\\}\Big{]},$

so we focus on $E\\{Y(z_{1},z_{2})\mid X_{0}\\}$. By Assumption 29.1(1) and
the tower property,

$\displaystyle E\\{Y(z_{1},z_{2})\mid X_{0}\\}$ $\displaystyle=$
$\displaystyle E\\{Y(z_{1},z_{2})\mid z_{1},X_{0}\\}$ $\displaystyle=$
$\displaystyle E\Big{[}E\\{Y(z_{1},z_{2})\mid z_{1},X_{1},X_{0}\\}\mid
z_{1},X_{0}\Big{]}.$

By Assumption 29.1(2),

$\displaystyle E\\{Y(z_{1},z_{2})\mid X_{0}\\}$ $\displaystyle=$
$\displaystyle E\Big{[}E\\{Y(z_{1},z_{2})\mid z_{2},z_{1},X_{1},X_{0}\\}\mid
z_{1},X_{0}\Big{]}$ $\displaystyle=$ $\displaystyle E\Big{[}E\\{Y\mid
z_{2},z_{1},X_{1},X_{0}\\}\mid z_{1},X_{0}\Big{]}.$

The formula (126) follows. $\square$

#####  148.1 Plug-in estimation based on outcome modeling

The g-formulas (127) and (128) suggest that to estimate the means of the
potential outcomes, we need to model $E(Y\mid z_{2},z_{1},x_{1},x_{0})$,
$\textup{pr}(x_{1}\mid z_{1},x_{0})$ and $\textup{pr}(x_{0})$. With these
fitted models, we can plug them in the g-formulas.

With some special functional forms, this task can be simplified. Example 29.1
below gives the results under a linear model for the outcome.

###### Example 29.1

Assume a linear outcome model

$E(Y\mid
z_{2},z_{1},x_{1},x_{0})=\beta_{0}+\beta_{1}z_{2}+\beta_{2}z_{1}+\beta_{3}x_{1}+\beta_{4}x_{0}.$

We can verify that

$\displaystyle E\\{Y(z_{1},z_{2})\\}$ $\displaystyle=$
$\displaystyle\sum_{x_{0}}\sum_{x_{1}}(\beta_{0}+\beta_{1}z_{2}+\beta_{2}z_{1}+\beta_{3}x_{1}+\beta_{4}x_{0})\textup{pr}(x_{1}\mid
z_{1},x_{0})\textup{pr}(x_{0})$ $\displaystyle=$
$\displaystyle\beta_{0}+\beta_{1}z_{2}+\beta_{2}z_{1}+\beta_{3}\sum_{x_{0}}E(X_{1}\mid
z_{1},x_{0})\textup{pr}(x_{0})+\beta_{4}E(X_{0}).$

Define

$E\\{X_{1}(z_{1})\\}=\sum_{x_{0}}E(X_{1}\mid z_{1},x_{0})\textup{pr}(x_{0})$
(130)

to simplify the formula as

$E\\{Y(z_{1},z_{2})\\}=\beta_{0}+\beta_{1}z_{2}+\beta_{2}z_{1}+\beta_{3}E\\{X_{1}(z_{1})\\}+\beta_{4}E(X_{0}).$

In (130), I introduce the potential outcome of $X_{1}$ under the treatment
$Z_{1}=z_{1}$ at time point 1. It is reasonable because the right-hand side of
(130) is the identification formula of $E\\{X_{1}(z_{1})\\}$ under
ignorability
$X_{1}(z_{1})\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}Z_{1}\mid X_{0}$ for $z_{1}=0,1$. We do not really need the
potential outcome $X_{1}(z_{1})$ and the ignorability, but it is a convenient
notation and matches with our discussion before.

Define $\tau_{Z_{1}\rightarrow X_{1}}=E\\{X_{1}(1)-X_{1}(0)\\}$. We can verify
that

$\displaystyle E\\{Y(1,0)-Y(0,0)\\}$ $\displaystyle=$
$\displaystyle\beta_{2}+\beta_{3}\tau_{Z_{1}\rightarrow X_{1}},$
$\displaystyle E\\{Y(0,1)-Y(0,0)\\}$ $\displaystyle=$
$\displaystyle\beta_{1},$ $\displaystyle E\\{Y(1,1)-Y(0,0)\\}$
$\displaystyle=$
$\displaystyle\beta_{1}+\beta_{2}+\beta_{3}\tau_{Z_{1}\rightarrow X_{1}}.$

Therefore, we can estimate the effect of $(Z_{1},Z_{2})$ on $Y$ based on the
above formulas by first estimating the regression coefficients $\beta$s and
the average causal effect of $Z_{1}$ on $X_{1}$ using standard methods.

However, Robins and Wasserman, (1997) pointed out a surprising drawback of the
plug-in estimation based on outcome modeling. They showed that with model
misspecification in this strategy, data analyzers may falsely reject the null
hypothesis of zero causal effect of $(Z_{1},Z_{2})$ on $Y$ even when the true
effect is zero in the data-generating process. They called it the g-null
paradox. Perhaps surprisingly, they show that the g-null paradox may even
arise in the simple linear outcome model in Example 29.1. McGrath et al.,
(2021) revisited this paradox. See Problem 29.1 for more details.

#####  148.2 Recursive estimation based on outcome modeling

The plug-in estimation in Section 148.1 involves modeling the time-varying
confounder $X_{1}$ and causes the unpleasant g-null paradox. It is not a
desirable method.

Recall the outcome regression estimator with a treatment at a single time
based on $E\\{Y(z)\\}=E\\{E(Y\mid Z=z,X)\\}$. We first fit a model of $Y$ on
$X$ using the subset of the data with $Z=z$, and obtain the fitted values
$\hat{Y}_{i}(z)$ for all units. We then obtain the estimator

$\hat{E}\\{Y(z)\\}=n^{-1}\sum_{i=1}^{n}\hat{Y}_{i}(z).$

Similarly, the recursive expectation formula in (126) motivates a simpler
method for estimation. Start from the inner conditional expectation, denoted
by

$\tilde{Y}_{2}(z_{1},z_{2})=E(Y\mid Z_{2}=z_{2},Z_{1}=z_{1},X_{1},X_{0}).$

We can fit a model of $Y$ on $(X_{1},X_{0})$ using the subset of the data with
$(Z_{2}=z_{2},Z_{1}=z_{1})$, and obtain the fitted values
$\hat{Y}_{2i}(z_{1},z_{2})$ for all units. Move on to outer conditional
expectation, denoted by

$\tilde{Y}_{1}(z_{1},z_{2})=E\\{\tilde{Y}_{2}(z_{1},z_{2})\mid
Z_{1}=z_{1},X_{0}\\}.$

We can fit a model of $\hat{Y}_{2}(z_{1},z_{2})$ on $X_{0}$ using the subset
of data with $Z_{1}=z_{1}$, and obtain the fitted values
$\hat{Y}_{1i}(z_{1},z_{2})$ for all units. The final estimator for
$E\\{Y(z_{1},z_{2})\\}$ is then

$\hat{E}\\{Y(z_{1},z_{2})\\}=n^{-1}\sum_{i=1}^{n}\hat{Y}_{1i}(z_{1},z_{2}).$

The above recursive estimation does not involve fitting a model for $X_{1}$
and avoids the g-null paradox. See Problem 29.2 for a special case.

####  149 Inverse propensity score weighting

Recall the propensity-score-based identification formula with a treatment at a
single time point:

$E\\{Y(z)\\}=E\left\\{\frac{1(Z=z)Y}{\textup{pr}(Z=z\mid X)}\right\\}.$

The following result extends it to the setting with a treatment at two time
points. Define

$e(z_{1},X_{0})=\textup{pr}(Z_{1}=z_{1}\mid X_{0})$

and

$e(z_{2},Z_{1},X_{1},X_{0})=\textup{pr}(Z_{2}=z_{2}\mid Z_{1},X_{1},X_{0})$

as the propensity scores at time points 1 and 2, respectively.

###### Theorem 29.2

Under Assumption 29.1,

$E\\{Y(z_{1},z_{2})\\}=E\left\\{\frac{1(Z_{1}=z_{1})1(Z_{2}=z_{2})Y}{e(z_{1},X_{0})e(z_{2},Z_{1},X_{1},X_{0})}\right\\}.$
(131)

Theorem 29.2 reveals the omitted overlap assumption:

$0<e(z_{1},X_{0})<1,\quad 0<e(z_{2},Z_{1},X_{1},X_{0})<1$

for all $z_{1}$ and $z_{2}$. If some propensity scores are 0 or 1, then the
identification formula (131) blows up to infinity.

Proof of Theorem 29.2: Conditioning on $(Z_{1},X_{1},X_{0})$ and using
Assumption 29.1(2), we can simplify the right-hand side of (131) as

$\displaystyle
E\left\\{\frac{1(Z_{1}=z_{1})1(Z_{2}=z_{2})Y(z_{1},z_{2})}{\textup{pr}(Z_{1}=z_{1}\mid
X_{0})\textup{pr}(Z_{2}=z_{2}\mid Z_{1},X_{1},X_{0})}\right\\}$ (132)
$\displaystyle=$ $\displaystyle
E\left\\{\frac{1(Z_{1}=z_{1})\textup{pr}(Z_{2}=z_{2}\mid
Z_{1},X_{1},X_{0})E(Y(z_{1},z_{2})\mid
Z_{1},X_{1},X_{0})}{\textup{pr}(Z_{1}=z_{1}\mid
X_{0})\textup{pr}(Z_{2}=z_{2}\mid Z_{1},X_{1},X_{0})}\right\\}$
$\displaystyle=$ $\displaystyle
E\left\\{\frac{1(Z_{1}=z_{1})}{\textup{pr}(Z_{1}=z_{1}\mid
X_{0})}E(Y(z_{1},z_{2})\mid Z_{1},X_{1},X_{0})\right\\}$ $\displaystyle=$
$\displaystyle E\left\\{\frac{1(Z_{1}=z_{1})}{\textup{pr}(Z_{1}=z_{1}\mid
X_{0})}Y(z_{1},z_{2})\right\\},$

where (132) follows from the tower property.

Conditioning on $X_{0}$ and using Assumption 29.1(1), we can simplify the
right-hand side of (132) as

$\displaystyle E\left\\{\frac{\textup{pr}(Z_{1}=z_{1}\mid
X_{0})}{\textup{pr}(Z_{1}=z_{1}\mid X_{0})}E(Y(z_{1},z_{2})\mid
X_{0})\right\\}$ $\displaystyle=$ $\displaystyle E\left\\{E(Y(z_{1},z_{2})\mid
X_{0})\right\\}$ $\displaystyle=$ $\displaystyle E\\{Y(z_{1},z_{2})\\},$

where, again, the last line follows from the tower property. $\square$

The estimator based on IPW is much simpler which only involves modeling two
binary indicators. First, we can fit a model of $Z_{1}$ on $X_{0}$ to obtain
the fitted values $\hat{e}_{1}(z_{1},X_{0i})$ and fit a model of $Z_{2}$ on
$(Z_{1},X_{1},X_{0})$ to obtain the fitted values
$\hat{e}_{2}(z_{2},Z_{1i},X_{1i},X_{0i})$ for all units. Then, we obtain the
following IPW estimator:

$\hat{E}^{\textup{ht}}\\{Y(z_{1},z_{2})\\}=n^{-1}\sum_{i=1}^{n}\frac{1(Z_{1i}=z_{1})1(Z_{2i}=z_{2})Y_{i}}{\hat{e}_{1}(z_{1},X_{0i})\hat{e}_{2}(z_{2},Z_{1i},X_{1i},X_{0i})}.$

Similar to the discussion in Chapter 11, the Horvitz–Thompson-type estimator
is not invariant to location shift of the outcome and suffers from instability
in finite samples. A modified Hajek-type estimator is
$\hat{E}^{\textup{haj}}\\{Y(z_{1},z_{2})\\}=\hat{E}^{\textup{ht}}\\{Y(z_{1},z_{2})\\}/\hat{1}^{\textup{ht}}(z_{1},z_{2})$,
where

$\hat{1}^{\textup{ht}}(z_{1},z_{2})=n^{-1}\sum_{i=1}^{n}\frac{1(Z_{1i}=z_{1})1(Z_{2i}=z_{2})}{\hat{e}_{1}(z_{1},X_{0i})\hat{e}_{2}(z_{2},Z_{1i},X_{1i},X_{0i})}.$

####  150 Multiple time points

Extending the estimation strategies in Sections 148 and 149 is not immediate
with multiple time points. Even with a binary treatment and $K$ time points,
the number of treatment combination grows exponentially with $K$ (for example,
$2^{5}=32$ and $2^{10}=1024$). Consequently, the outcome regression and IPW
estimators in Sections 148 and 149 are not feasible in finite samples.

#####  150.1 Marginal structural model

A powerful approach is based on the marginal structural model (MSM) (Robins et
al.,, 2000; Hernán et al.,, 2000). For simplicity of notation, I will only
present the MSM with $K=2$ although its main use is in the general case.

######  Definition 29.1 (MSM)

The marginal mean of $Y(z_{1},z_{2})$ equals

$E\\{Y(z_{1},z_{2})\\}=f(z_{1},z_{2};\beta).$

A leading example of Definition 29.1 is
$E\\{Y(z_{1},z_{2})\\}=\beta_{0}+\beta_{1}z_{1}+\beta_{2}z_{2}$. It is also
straightforward to include the baseline covariates in the model. Definition
29.2 below extends Definition 29.1.

######  Definition 29.2 (MSM with baseline covariates)

The mean of $Y(z_{1},z_{2})$ conditional on $X_{0}$ equals

$E\\{Y(z_{1},z_{2})\mid X_{0}\\}=f(z_{1},z_{2},X_{0};\beta).$

A leading example of Definition 29.2 is

$E\\{Y(z_{1},z_{2})\mid
X_{0}\\}=\beta_{0}+\beta_{1}z_{1}+\beta_{2}z_{2}+\beta_{3}^{\mkern-1.5mu\mathsf{T}}X_{0}.$
(133)

If we observe all potential outcomes, we can solve $\beta$ from the following
minimization problem:

$\beta=\arg\min_{b}\sum_{z_{2}}\sum_{z_{1}}E\\{Y(z_{1},z_{2})-f(z_{1},z_{2},X_{0};b)\\}^{2}.$

For simplicity, I focus on the least squares formulation. We can also extend
the discussion to a general loss function.

Under sequential ignorability, we can solve $\beta$ from the following
minimization problem that only involves observables.

######  Theorem 29.3 (IPW under MSM)

Under Assumption 29.1 and Definition 29.2,

$\beta=\arg\min_{b}\sum_{z_{2}}\sum_{z_{1}}E\left[\frac{1(Z_{1}=z_{1})1(Z_{2}=z_{2})}{e(z_{1},X_{0})e(z_{2},Z_{1},X_{1},X_{0})}\\{Y-f(z_{1},z_{2},X_{0};b)\\}^{2}\right].$

The proof of Theorem 29.3 is similar to that of Theorem 29.2. I relegate it to
Problem 29.3.

Theorem 29.3 implies a simple estimation strategy based on weighted
regressions. For instance, under (133), we can fit WLS of $Y_{i}$ on
$(1,Z_{1i},Z_{2i},X_{0i})$ with weights
$\hat{e}_{1}^{-1}(Z_{1i},X_{0i})\hat{e}_{2i}^{-1}(Z_{2i},Z_{i1},X_{1i},X_{0i})$.

#####  150.2 Structural nested model

A key problem of IPW is that it is not applicable if the overlap assumption is
violated. To address this challenge, Robins proposed the structural nested
model. Again, to simplify the presentation, I only review the version with two
time points.

######  Definition 29.3 (structural nested model)

The conditional effect at time point 1 is

$E\\{Y(z_{1},0)-Y(0,0)\mid Z_{1}=z_{1},X_{0}\\}=g_{1}(z_{1},X_{0};\beta)\text{
for all }z_{1}$

and the conditional effect at time point 2 is

$E\\{Y(z_{1},z_{2})-Y(z_{1},0)\mid
Z_{2}=z_{2},Z_{2}=z_{1},X_{1},X_{0}\\}=g_{2}(z_{2},z_{1},X_{1},X_{0};\beta)\text{
for all }z_{1},z_{2}.$

In Definition 29.3, two logical restrictions are

$g_{1}(0,X_{0};\beta)=0$

and

$g_{2}(0,z_{1},X_{1},X_{0};\beta)=0\text{ for all }z_{1}.$

Two leading choices of Definition 29.3 are below.

###### Example 29.2

Assume

$\begin{cases}g_{1}(z_{1},X_{0};\beta)=\beta_{1}z_{1},\\\
g_{2}(z_{2},z_{1},X_{1},X_{0};\beta)=(\beta_{2}+\beta_{3}z_{1})z_{2}.\end{cases}$

###### Example 29.3

Assume

$\begin{cases}g_{1}(z_{1},X_{0};\beta)=(\beta_{1}+\beta^{\mkern-1.5mu\mathsf{T}}X_{0})z_{1},\\\
g_{2}(z_{2},z_{1},X_{1},X_{0};\beta)=(\beta_{2}+\beta_{3}z_{1}+\beta_{4}^{\mkern-1.5mu\mathsf{T}}X_{1})z_{2}.\end{cases}$

Compare Definitions 29.2 and 29.3. The structural nested model allows for
adjusting for the time-varying covariates whereas the marginal structural
model only allows for adjusting for baseline covariates. The estimation under
Definition 29.3 is more involved. A strategy is to estimate the parameter
based on estimating equations.

I first introduce two important building blocks for the discussing the
estimation. Define

$U_{2}(\beta)=Y-g_{2}(Z_{2},Z_{1},X_{1},X_{0};\beta)$

and

$U_{1}(\beta)=Y-g_{2}(Z_{2},Z_{1},X_{1},X_{0};\beta)-g_{1}(Z_{1},X_{0};\beta).$

They are not directly computable from the data because they depend on the true
value of the parameter $\beta$. At the true value, they have the following
properties.

###### Lemma 29.1

Under Assumption 29.1 and Definition 29.3, we have

$\displaystyle E\\{U_{2}(\beta)\mid Z_{2},Z_{1},X_{1},X_{0}\\}$
$\displaystyle=$ $\displaystyle E\\{U_{2}(\beta)\mid Z_{1},X_{1},X_{0}\\}$
$\displaystyle=$ $\displaystyle E\\{Y(Z_{1},0)\mid Z_{1},X_{1},X_{0}\\}$

and

$\displaystyle E\\{U_{1}(\beta)\mid Z_{1},X_{0}\\}$ $\displaystyle=$
$\displaystyle E\\{U_{1}(\beta)\mid X_{0}\\}$ $\displaystyle=$ $\displaystyle
E\\{Y(0,0)\mid X_{0}\\}.$

Lemma 29.1 involves a subtle notation $Y(Z_{1},0)$ because $Z_{1}$ is random.
It should be read as $Y(Z_{1},0)=Z_{1}Y(1,0)+(1-Z_{1})Y(0,0)$. Based on the
definitions and Lemma 29.1, $U_{1}(\beta)$ acts as the control potential
outcome before receiving any treatment and $U_{2}(\beta)$ acts as the control
potential outcome after receiving the treatment at time point 1.

Proof of Lemma 29.1: First, we have

$\displaystyle E\\{U_{2}(\beta)\mid Z_{2}=1,Z_{1},X_{1},X_{0}\\}$
$\displaystyle=$ $\displaystyle
E\\{Y(Z_{1},1)-g_{2}(1,Z_{1},X_{1},X_{0};\beta)\mid
Z_{2}=1,Z_{1},X_{1},X_{0}\\}$ $\displaystyle=$ $\displaystyle
E\\{Y(Z_{1},0)\mid Z_{2}=1,Z_{1},X_{1},X_{0}\\}$ $\displaystyle
E\\{U_{2}(\beta)\mid Z_{2}=0,Z_{1},X_{1},X_{0}\\}$ $\displaystyle=$
$\displaystyle E\\{Y(Z_{1},0)-g_{2}(0,Z_{1},X_{1},X_{0};\beta)\mid
Z_{2}=0,Z_{1},X_{1},X_{0}\\}$ $\displaystyle=$ $\displaystyle
E\\{Y(Z_{1},0)\mid Z_{2}=0,Z_{1},X_{1},X_{0}\\}$

so

$E\\{U_{2}(\beta)\mid Z_{2},Z_{1},X_{1},X_{0}\\}=E\\{Y(Z_{1},0)\mid
Z_{2},Z_{1},X_{1},X_{0}\\}=E\\{Y(Z_{1},0)\mid Z_{1},X_{1},X_{0}\\}$

where the last identity follows from sequential ignorability. Since the last
term does not depend on $Z_{2}$, we also have

$E\\{U_{2}(\beta)\mid Z_{2},Z_{1},X_{1},X_{0}\\}=E\\{U_{2}(\beta)\mid
Z_{1},X_{1},X_{0}\\}.$

Using the above results, we have

$\displaystyle E\\{U_{1}(\beta)\mid Z_{1},X_{0}\\}$ $\displaystyle=$
$\displaystyle E\\{U_{2}(\beta)-g_{1}(Z_{1},X_{0};\beta)\mid Z_{1},X_{0}\\}$
$\displaystyle=$ $\displaystyle
E\left[E\\{U_{2}(\beta)-g_{1}(Z_{1},X_{0};\beta)\mid X_{1},Z_{1},X_{0}\\}\mid
Z_{1},X_{0}\right]$ $\displaystyle=$ $\displaystyle
E\left[E\\{Y(Z_{1},0)-g_{1}(Z_{1},X_{0};\beta)\mid X_{1},Z_{1},X_{0}\\}\mid
Z_{1},X_{0}\right]$ $\displaystyle=$ $\displaystyle
E\\{Y(Z_{1},0)-g_{1}(Z_{1},X_{0};\beta)\mid Z_{1},X_{0}\\}$ $\displaystyle=$
$\displaystyle E\\{Y(0,0)\mid Z_{1},X_{0}\\}$ $\displaystyle=$ $\displaystyle
E\\{Y(0,0)\mid X_{0}\\}$

where the last identity follows from sequential ignorability. Since the last
term does not depend on $Z_{1}$, we also have

$E\\{U_{1}(\beta)\mid Z_{1},X_{0}\\}=E\\{U_{1}(\beta)\mid X_{0}\\}.$

$\square$

With Lemma 29.1, we can prove Theorem 29.4 below.

###### Theorem 29.4

Under Assumption 29.1 and Definition 29.3,

$E\Big{[}h_{2}(Z_{1},X_{1},X_{0})\\{Z_{2}-e(1,Z_{1},X_{1},X_{0})\\}U_{2}(\beta)\Big{]}=0$

and

$E\Big{[}h_{1}(X_{0})\\{Z_{1}-e(1,X_{0})\\}U_{1}(\beta)\Big{]}=0.$

for any functions $h_{1}$ and $h_{2}$, provided that the moments exist.

Proof of Theorem 29.2: Use the tower property by conditioning on
$(Z_{2},Z_{1},X_{1},X_{0})$ and Lemma 29.1 to obtain

$\displaystyle
E\left[h_{2}(Z_{1},X_{1},X_{0})\\{Z_{2}-e(1,Z_{1},X_{1},X_{0})\\}E\\{U_{2}(\beta)\mid
Z_{2},Z_{1},X_{1},X_{0}\\}\right]$ $\displaystyle=$ $\displaystyle
E\left[h_{2}(Z_{1},X_{1},X_{0})\\{Z_{2}-e(1,Z_{1},X_{1},X_{0})\\}E\\{U_{2}(\beta)\mid
Z_{1},X_{1},X_{0}\\}\right].$

Use the tower property by conditioning on $(Z_{1},X_{1},X_{0})$ to show that
the last identity equals 0.

Similarly, use the tower property by conditioning on $(Z_{1},X_{0})$ and Lemma
29.1 to obtain

$\displaystyle E\left[h_{1}(X_{0})\\{Z_{1}-e(1,X_{0})\\}E\\{U_{1}(\beta)\mid
Z_{1},X_{0}\\}\right]$ $\displaystyle=$ $\displaystyle
E\left[h_{1}(X_{0})\\{Z_{1}-e(1,X_{0})\\}E\\{U_{1}(\beta)\mid
X_{0}\\}\right].$

Use the tower property by conditioning on $X_{0}$ to show that the last
identity equals 0. $\square$

To use Theorem 29.4, we must specify $h_{1}$ and $h_{2}$ to ensure that there
are enough equations for solving $\beta$. Example 29.4 below revisits Example
29.2.

###### Example 29.4

Under Example 29.2, we can choose $h_{1}=1$ and $h_{2}=(1,Z_{1})$ to obtain

$\displaystyle
E\left[\\{Z_{2}-e(1,Z_{1},X_{1},X_{0})\\}\\{Y-(\beta_{2}+\beta_{3}Z_{1})Z_{2}\\}\right]$
$\displaystyle=$ $\displaystyle 0,$ $\displaystyle
E\left[Z_{1}\\{Z_{2}-e(1,Z_{1},X_{1},X_{0})\\}\\{Y-(\beta_{2}+\beta_{3}Z_{1})Z_{2}\\}\right]$
$\displaystyle=$ $\displaystyle 0,$ $\displaystyle
E\left[\\{Z_{1}-e(1,X_{0})\\}\\{Y-(\beta_{2}+\beta_{3}Z_{1})Z_{2}-\beta_{1}Z_{1}\\}\right]$
$\displaystyle=$ $\displaystyle 0.$

We can then solve for the $\beta$’s from the above linear equations; see
Problem 29.5. A natural question is that whether alternative choices of
$(h_{1},h_{2})$ can lead to more efficient estimators. The answer is yes. For
example, we can choose many $(h_{1},h_{2})$ and use the generalized method of
moment (Hansen,, 1982). The technical details are beyond this book.

Naimi et al., (2017) and Vansteelandt and Joffe, (2014) provided tutorials on
the structural nested models.

####  151 Homework problems

######  29.1 g-null paradox

Consider the simple causal diagram in Figure 44 without pre-treatment
covariates $X_{0}$ and without the arrows from $(Z_{1},Z_{2})$ to $Y$. So the
effect of $(Z_{1},Z_{2})$ on $Y$ is zero.

---
$\textstyle{Z_{1}\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{X_{1}\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$$\textstyle{Z_{2}}$$\textstyle{Y}$$\textstyle{U\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces\ignorespaces}$
Figure 44: With unmeasured confounding between $X_{1}$ and $Y$. The causal
diagram ignores the pre-treatment covariates $X_{0}$.

Revisit Example 29.1. Show that the expectation $E\\{Y(z_{1},z_{2})\\}$ does
not depend on $(z_{1},z_{2})$ if

$\beta_{1}=\beta_{3}=0\text{ and }\beta_{2}=0$

or

$\beta_{1}=\beta_{3}=0\text{ and }E\\{X_{1}(z_{1})\\}\text{ does not depend on
}z_{1}.$

holds.

Remark: However, $\beta_{2}=0$ in the first condition rules out the dependence
of $Y$ on $X_{1}$, contradicting with the existence of unmeasured confounder
$U$ between $X_{1}$ and $Y$; the independence of $E\\{X_{1}(z_{1})\\}$ on
$z_{1}$ rules out the dependence of $X_{1}$ on $Z_{1}$, contradicting with the
existence of the arrow from $Z_{1}$ on $X_{1}$. That is, if there is an
unmeasured confounder $U$ between $X_{1}$ and $Y$ and there is an arrow from
$Z_{1}$ on $X_{1}$, then the formula of $E\\{Y(z_{1},z_{2})\\}$ in Example
29.1 must depend on $(z_{1},z_{2})$, which leads to a contradiction with the
absence of arrows from $(Z_{1},Z_{2})$ to $Y$.

######  29.2 Recursive estimation under the null model

Consider the recursive estimation method in 148.2 under the causal diagram in
Problem 29.1. Show that based linear models, the estimator converges to 0.

######  29.3 IPW under MSM

Prove Theorem 29.3.

######  29.4 Structural nested model with a single time point

Recall the standard setting of observational studies with IID data draw from
$\\{X,Z,Y(1),Y(0)\\}$. Define the propensity score as
$e(X)=\textup{pr}(Z=1\mid X)$. Assume

$Z\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}Y(0)\mid X$

and the following structural nested model.

######  Definition 29.4 (structural nested model with a single time point)

The conditional mean of the individual effect is

$E\\{Y(z)-Y(0)\mid Z=z,X\\}=g(z,X;\beta).$

In Definition 29.4, a logical restriction is $g(0,X;\beta)=0$. Prove the
following results.

  1. 1.

We have

$E\\{Y-g(Z,X;\beta)\mid X,Z\\}=E\\{Y-g(Z,X;\beta)\mid X\\}=E\\{Y(0)\mid X\\}.$

  2. 2.

We have

$\displaystyle E\Big{[}h(X)\\{Z-e(X)\\}\\{Y-g(Z,X;\beta)\\}\Big{]}=0$ (134)

for any function $h$, provided that the moment exists.

Remark: (134) is the basis for parameter estimation. Consider a special case
of Definition 29.4 with $g(z,X;\beta)=\beta z$. Choose $h(X)=1$ to obtain

$E\\{(Z-e(X))(Y-\beta Z)\\}=0.$

Solve for $\beta$ to obtain

$\beta=\frac{E\\{(Z-e(X))Y\\}}{E\\{(Z-e(X))Z\\}}.$

That is, $\beta$ equals the coefficient of $Z$ in the two-stage least squares
of $Y$ on $Z$ with $Z-e(X)$ being the instrument variable for $Z$

Consider a special case of Definition 29.4 with
$g(z,X;\beta)=(\beta_{0}+\beta_{1}^{\mkern-1.5mu\mathsf{T}}X)z$. Choose
$h(X)=(1,X)$ to obtain

$E\left\\{\begin{pmatrix}Z-e(X)\\\
(Z-e(X))X\end{pmatrix}(Y-\beta_{0}Z-\beta_{1}^{\mkern-1.5mu\mathsf{T}}XZ)\right\\}=0.$

That is, $(\beta_{0},\beta_{1})$ equal the coefficients in the two-stage least
squares of $Y$ on $(Z,XZ)$ with $(Z-e(X),(Z-e(X))X)$ being the instrument
variable for $(Z,XZ)$.

######  29.5 Estimation under Example 29.4

We can estimate the $\beta$’s by solving the empirical version of the
estimating equations in Example 29.4. We first estimate the two propensity
scores and obtain the centered treatment

$\check{Z}_{1i}=Z_{1i}-\hat{e}(1,X_{0i})$

at time point 1 and

$\check{Z}_{2i}=Z_{2i}-\hat{e}(1,Z_{1i},X_{1i},X_{0i})$

at time point 2.

Show that we can estimate $\beta_{2}$ and $\beta_{3}$ by running two-stage
least squares of $Y_{i}$ on $(Z_{2i},Z_{1i}Z_{2i})$ with
$(\check{Z}_{2i},Z_{1i}\check{Z}_{2i})$ as the instrumental variable for
$(Z_{2i},Z_{1i}Z_{2i})$, and then we can estimate $\beta_{1}$ by running two-
stage least squares of $Y_{i}-(\hat{\beta}_{2}+\hat{\beta}_{3}Z_{1i})Z_{2i}$
on $Z_{1i}$ with $\check{Z}_{1i}$ as the instrumental variable for $Z_{1i}$.

######  29.6 g-formula with a treatment at multiple time points

Extend the discussion to the setting with $K$ time points. The temporal
ordering of the variables is

$X_{0}\rightarrow Z_{1}\rightarrow X_{1}\rightarrow Z_{2}\rightarrow\cdots
X_{K-1}\rightarrow Z_{K}.$

Introduce the notation $\overline{Z}_{k}=(Z_{1},\ldots,Z_{k})$ and
$\overline{X}_{k}=(X_{0},X_{1},\ldots,X_{k})$ with lower case
$\overline{z}_{k}$ and $\overline{x}_{k}$ denoting the corresponding realized
values. With $k=0$, we have $\overline{X}_{0}=X_{0}$ and $\overline{Z}_{0}$ is
empty. Each unit has $2^{K}$ potential outcomes:

$Y(\overline{z}_{K})\text{ for all }z_{1},\ldots,z_{K}=0,1.$

Assume sequential ignorability below.

######  Assumption 29.2 (sequential ignorability at multiple time points)

We have

$Z_{k}\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}Y(\overline{z}_{K})\mid(\overline{Z}_{k-1},\overline{X}_{k-1})$

for all $k=1,\ldots,K$ and all $z_{1},\ldots,z_{K}=0,1.$

Prove Theorem 29.5 below.

######  Theorem 29.5 (g-formula with multiple time points)

Under Assumption 29.2,

$E\\{Y(\overline{z}_{K})\\}=E\left[\cdots
E\\{E(Y\mid\overline{z}_{K},\overline{X}_{K-1})\mid\overline{z}_{K-1},\overline{X}_{K-2}\\}\cdots\mid
z_{1},X_{0}\right].$

Remark: In Theorem 29.5, I use the simplified notation “$\overline{z}_{k}$”
for “$\overline{Z}_{k}=\overline{z}_{k}$.” With discrete $X$, Theorem 29.5
reduces to

$\displaystyle E\\{Y(\overline{z}_{K})\\}$ $\displaystyle=$
$\displaystyle\sum_{x_{0}}\sum_{x_{1}}\cdots\sum_{x_{K-1}}E(Y\mid\overline{z}_{K},\overline{x}_{K-1})$
$\displaystyle\qquad\qquad\cdot\textup{pr}(x_{K-1}\mid\overline{z}_{K-1},\overline{x}_{K-2})\cdots\textup{pr}(x_{1}\mid
z_{1},x_{0})\textup{pr}(x_{0});$

with continuous $X$, Theorem 29.5 reduces to

$\displaystyle E\\{Y(\overline{z}_{K})\\}$ $\displaystyle=$ $\displaystyle\int
E(Y\mid\overline{z}_{K},\overline{x}_{K-1})$ $\displaystyle\cdot
f(x_{K-1}\mid\overline{z}_{K-1},\overline{x}_{K-2})\cdots f(x_{1}\mid
z_{1},x_{0})f(x_{0})\textnormal{d}\overline{x}_{K-1}.$

######  29.7 IPW with a treatment at multiple time points

Inherit the setting of Problem 29.6. Define the propensity score at $K$ time
points as

$\displaystyle e(z_{1},X_{0})$ $\displaystyle=$
$\displaystyle\textup{pr}(Z_{1}=z_{1}\mid X_{0}),$ $\displaystyle\vdots$
$\displaystyle e(z_{k},\overline{Z}_{k-1},\overline{X}_{k-1})$
$\displaystyle=$
$\displaystyle\textup{pr}(Z_{k}=z_{k}\mid\overline{Z}_{k-1},\overline{X}_{k-1}),$
$\displaystyle\vdots$ $\displaystyle
e(z_{K},\overline{Z}_{K-1},\overline{X}_{K-1})$ $\displaystyle=$
$\displaystyle\textup{pr}(Z_{K}=z_{K}\mid\overline{Z}_{K-1},\overline{X}_{K-1}).$

Prove Theorem 29.7 below assuming overlap implicitly.

######  Theorem 29.6 (IPW with multiple time points)

Under Assumption 29.2,

$E\\{Y(\overline{z}_{K})\\}=E\left\\{\frac{1(Z_{1}=z_{1})\cdots
1(Z_{K}=z_{K})Y}{e(z_{1},X_{0})\cdots
e(z_{K},\overline{Z}_{K-1},\overline{X}_{K-1})}\right\\}.$

Based on Theorem 29.7, construct the Horvitz–Thompson and Hajek estimators.

######  29.8 MSM with a treatment at multiple time points

The number of potential outcomes grows exponentially with $K$. The formulas in
Problems 29.6 and 29.7 are not directly applicable in finite samples. We can
impose the following structural assumptions on the potential outcomes.

######  Definition 29.5 (MSM with multiple time points)

Assume

$E\\{Y(\overline{z}_{K})\mid X_{0}\\}=f(\overline{z}_{K},X_{0};\beta).$

Two leading examples of Definition 29.5 are $E\\{Y(\overline{z}_{K})\mid
X_{0}\\}=\beta_{0}+\beta_{1}\sum_{k=1}^{K}z_{k}+\beta_{2}^{\mkern-1.5mu\mathsf{T}}X_{0}$
and $E\\{Y(\overline{z}_{K})\mid
X_{0}\\}=\beta_{0}+\sum_{k=1}^{K}\beta_{k}z_{k}+\beta_{K+1}^{\mkern-1.5mu\mathsf{T}}X_{0}.$

If we know all the potential outcomes, we can solve $\beta$ from the following
minimization problem:

$\beta=\arg\min_{b}\sum_{\overline{z}_{K}}E\\{Y(\overline{z}_{K})-f(\overline{z}_{K},X_{0};\beta)\\}^{2}.$

Theorem 29.7 below shows that under Assumption 29.2, we can solve $\beta$ from
a minimization problem that only involves observables.

######  Theorem 29.7 (IPW for MSM with multiple time points)

Under Assumption 29.2,

$\beta=\arg\min_{b}\sum_{\overline{z}_{K}}E\left[\frac{1(Z_{1}=z_{1})\cdots
1(Z_{K}=z_{K})}{e(z_{1},X_{0})\cdots
e(z_{K},\overline{Z}_{K-1},\overline{X}_{K-1})}\\{Y-f(\overline{z}_{K},X_{0};\beta)\\}^{2}\right].$

######  29.9 Structural nested model with a treatment at multiple time points

Inherit the setting from Problem 29.6 and the notation from Problem 29.7. This
problem presents a general structural nested model.

######  Definition 29.6 (structural nested model with multiple time points)

The conditional effect at time $k$ is

$E\\{Y(\overline{z}_{k},0)-Y(\overline{z}_{k-1},0)\mid\overline{z}_{k},\overline{X}_{k-1}\\}=g_{k}(\overline{z}_{k},\overline{X}_{k-1};\beta)$

for all $\overline{z}_{k}$ and all $k=1,\ldots,K$.

In Definition 29.6, a logical restriction is

$g_{k}(0,\overline{z}_{k-1},\overline{X}_{k-1};\beta)=0$

for all $\overline{z}_{k-1}$ and all $k=1,\ldots,K$.

Define

$U_{k}(\beta)=Y-\sum_{s=1}^{k}g_{s}(\overline{Z}_{s},\overline{X}_{s-1};\beta)$

for all $k=1,\ldots,K$. Theorem 29.8 below extends Theorem 29.4.

###### Theorem 29.8

Under Assumption 29.2 and Definition 29.6,

$E\left[h_{k}(\overline{Z}_{k-1},\overline{X}_{k-1})\\{Z_{k}-e(1,\overline{Z}_{k-1},\overline{X}_{k-1})\\}U_{k}(\beta)\right]=0$

for all $k=1,\ldots,K$.

Remark: Choosing appropriate $h_{k}$’s, we can estimate $\beta$ by solving the
empirical version of Theorem 29.8.

######  29.10 Recommended reading

Robins et al., (2000) reviewed the MSM. Naimi et al., (2017) reviewed the
g-methods.

##  Part VII Appendices

###  Chapter A1 Probability and Statistics

####  152 Probability

#####  152.1 Tower property and variance decomposition

Given random variables or vectors $A,B,C$, we have

$E(A)=E\\{E(A\mid B)\\}$

and

$E(A\mid C)=E\\{E(A\mid B,C)\mid C\\}.$

Given a random variable $A$ and random variables or vectors $B,C$, we have

$\textup{var}(A)=E\\{\textup{var}(A\mid B)\\}+\textup{var}\\{E(A\mid B)\\}$

and

$\textup{var}(A\mid C)=E\\{\textup{var}(A\mid B,C)\mid
C\\}+\textup{var}\\{E(A\mid B,C)\mid C\\}.$

Similarly, we can decompose the covariance as

$\textup{cov}(A_{1},A_{2})=E\left\\{\textup{cov}(A_{1},A_{2}\mid
B)\right\\}+\textup{cov}\\{E(A_{1}\mid B),E(A_{2}\mid B)\\}$

and

$\textup{cov}(A_{1},A_{2}\mid C)=E\left\\{\textup{cov}(A_{1},A_{2}\mid
B,C)\mid C\right\\}+\textup{cov}\\{E(A_{1}\mid B,C),E(A_{2}\mid B,C)\mid
C\\}.$

#####  152.2 Limiting theorems

######  Definition A1.1 (convergence in probability)

A sequence of random variables $(X_{n})_{n\geq 1}$ converges to $X$ in
probability, if for every $\varepsilon>0$, we have

$\textup{pr}(|X_{n}-X|>\varepsilon)\rightarrow 0$

as $n\rightarrow\infty$.

######  Definition A1.2 (convergence in distribution)

A sequence of random variables $(X_{n})_{n\geq 1}$ converges to $X$ in
distribution, if

$\textup{pr}(X_{n}\leq x)\rightarrow\textup{pr}(X\leq x)$

for all continuity point $x$ of $\textup{pr}(X\leq x)$, as
$n\rightarrow\infty$.

Convergence in probability is stronger than convergence in distribution.
Definitions A1.1 and A1.2 are useful for stating the following two fundamental
theorems in probability theory.

######  Theorem A1.1 (law of large numbers)

If $X_{1},\ldots,X_{n}\stackrel{{\scriptstyle\text{IID}}}{{\sim}}X$ with
$E|X|<\infty$, then $\bar{X}=n^{-1}\sum_{i=1}^{n}X_{i}\rightarrow E(X)$ in
probability.

The law of large numbers in Theorem A1.1 states that the sample average is
close to the population mean in the limit.

######  Theorem A1.2 (central limit theorem)

If $X_{1},\ldots,X_{n}\stackrel{{\scriptstyle\text{IID}}}{{\sim}}X$ with
$\textup{var}(X)<\infty$, then

$\frac{\bar{X}-E(X)}{\sqrt{\textup{var}(X)/n}}\rightarrow\textsc{N}(0,1)$

in distribution.

The central limit theorem in Theorem A1.2 states that the standardized sample
average is close to a standard Normal random variable in the limit.

Theorems A1.1 and A1.2 assume IID random variables for convenience. There are
also many law of large numbers and central limit theorems for the sample mean
of independent random variable (e.g., Durrett,, 2019).

#####  152.3 Delta method

Delta method is a power tool to derive asymptotic Normality of nonlinear
functions of an asymptotically Normal random vector. I review a special case
of delta method below.

######  Theorem A1.3 (delta method)

Assume $\sqrt{n}(X_{n}-\mu)\rightarrow\textsc{N}(0,\Sigma)$ in distribution
and the function $g(x)$ has non-zero derivative $\nabla g(\mu)$ at $\mu$. Then

$\sqrt{n}\\{g(X_{n})-g(\mu)\\}\rightarrow\textsc{N}\left(0,(\nabla
g(\mu)^{\mkern-1.5mu\mathsf{T}}\Sigma\nabla g(\mu)\right)$

in distribution.

I will omit the proof of Theorem A1.3. It is intuitive based on the first-
order Taylor expansion:

$g(X_{n})-g(\mu)\cong(\nabla g(\mu)^{\mkern-1.5mu\mathsf{T}}(X_{n}-\mu).$

A leading example of delta method is to obtain the asymptotic Normality of
ratio.

######  Example A1.1 (asymptotic normality for ratio)

Assume

$\sqrt{n}\begin{pmatrix}Y_{n}-\mu_{Y}\\\
X_{n}-\mu_{X}\end{pmatrix}\rightarrow\textsc{N}\left(\begin{pmatrix}0\\\
0\end{pmatrix},\begin{pmatrix}\sigma_{Y}^{2}&\sigma_{YX}\\\
\sigma_{YX}&\sigma_{X}^{2}\end{pmatrix}\right)$ (135)

in distribution with $\mu_{X}\neq 0.$ Apply Theorem A1.3 to obtain that

$\sqrt{n}\left(\frac{Y_{n}}{X_{n}}-\frac{\mu_{Y}}{\mu_{X}}\right)\rightarrow\textsc{N}\left(0,\frac{\sigma_{Y}^{2}}{\mu_{X}^{2}}+\frac{\mu_{Y}^{2}\sigma_{X}^{2}}{\mu_{X}^{4}}-\frac{2\mu_{Y}\sigma_{YX}}{\mu_{X}^{3}}\right)$
(136)

in distribution. In the special case that $X_{n}$ and $Y_{n}$ are
asymptotically independent, the asymptotic variance of $Y_{n}/X_{n}$
simplifies to
$\sigma_{Y}^{2}/\mu_{X}^{2}+\mu_{Y}^{2}\sigma_{X}^{2}/\mu_{X}^{4}$. I leave
the details to Problem A1.2.

The asymptotic variance in Example A1.1 is a little cumbersome. An easier way
to memorize it is based on the following approximation:

$\frac{Y_{n}}{X_{n}}-\frac{\mu_{Y}}{\mu_{X}}=\frac{Y_{n}-\mu_{Y}/\mu_{X}\cdot
X_{n}}{X_{n}}\cong\frac{Y_{n}-\mu_{Y}/\mu_{X}\cdot X_{n}}{\mu_{X}},$ (137)

so the asymptotic variance of the ratio equals the asymptotic variance of

$\frac{Y_{n}-\mu_{Y}/\mu_{X}\cdot X_{n}}{\mu_{X}},$

which is a linear combination of $Y_{n}$ and $X_{n}$. Slutsky’s theorem can
make the approximation in (137) rigorous; it is beyond this book.

######  Example A1.2 (asymptotic normality for product)

Assume (135). Apply Theorem A1.3 to obtain that

$\sqrt{n}\left(X_{n}Y_{n}-\mu_{X}\mu_{Y}\right)\rightarrow\textsc{N}\left(0,\mu_{Y}^{2}\sigma_{X}^{2}+\mu_{X}^{2}\sigma_{Y}^{2}+2\mu_{X}\mu_{Y}\sigma_{XY}\right)$
(138)

in distribution. In the special case that $X_{n}$ and $Y_{n}$ are
asymptotically independent, the asymptotic variance of $X_{n}Y_{n}$ simplifies
to $\mu_{Y}^{2}\sigma_{X}^{2}+\mu_{X}^{2}\sigma_{Y}^{2}$. I leave the details
to Problem A1.3.

####  153 Statistical inference

#####  153.1 Point estimation

Assume that $\theta$ is the parameter of interest. Oftentimes, the problem
also contain other parameters not of interest, denoted by $\eta$.
Statisticians call $\eta$ the nuisance parameter. Based on data, we can
compute an estimator $\hat{\theta}$. Throughout this book, we take the
frequentist’s perspective by assuming that $\theta$ is a fixed number and
$\hat{\theta}$ is random due to the randomness of data. Two basic requirements
for an estimator are below.

######  Definition A1.3 (unbiasedness)

The estimator $\hat{\theta}$ is unbiased for $\theta$ if

$E(\hat{\theta})=\theta$

for all possible values of $\theta$ and $\eta$.

######  Definition A1.4 (consistency)

The estimator $\hat{\theta}$ is consistent for $\theta$ if

$\hat{\theta}\rightarrow\theta$

in probability as the sample size approaches to infinity, for all possible
values of $\theta$ and $\eta$.

Unbiasedness requires that the mean of the estimator is identical to the
parameter of interest. Consistency requires that the estimator is close to the
true parameter in the limit. Unbiased does not imply consistency, and
consistency does not imply unbiasedness either. Unbiasedness can be
restrictive because it is impossible even in some simple statistics problems.
Consistency is often the basic requirement in most statistics problems.

#####  153.2 Confidence interval

A point estimator $\hat{\theta}$ is a random variable which differs from the
true parameter. Statisticians are often interested in finding an interval that
covers the true parameter with certain given probability. This interval is
computed based on the data, and it is random.

######  Definition A1.5 (confidence interval)

A data-dependent interval
$[\hat{\theta}_{\textsc{l}},\hat{\theta}_{\textsc{u}}]$ is a confidence
interval for $\theta$ with coverage probability $1-\alpha$ if

$\textup{pr}(\hat{\theta}_{\textsc{l}}\leq\theta\leq\hat{\theta}_{\textsc{u}})\geq
1-\alpha.$

######  Definition A1.6 (asymptotic confidence interval)

A data-dependent interval
$[\hat{\theta}_{\textsc{l}}],\hat{\theta}_{\textsc{u}}]$ is an asymptotic
confidence interval for $\theta$ with coverage probability $1-\alpha$ if

$\textup{pr}(\hat{\theta}_{\textsc{l}}\leq\theta\leq\hat{\theta}_{\textsc{u}})\rightarrow
1-\alpha^{\prime}$

with $\alpha^{\prime}\geq\alpha$, as $n\rightarrow\infty$.

A standard choice is $\alpha=0.05$. In Definitions A1.5 and A1.6, the coverage
probabilities can be larger than the nominal level $1-\alpha$. That is, the
definitions allow for over overage but do not allow for under coverage. With
over coverage, we say that the confidence interval is conservative. Of course,
we hope the confidence interval to be as narrow as possible. Otherwise, the
definition of confidence interval can be arbitrary.

#####  153.3 Hypothesis testing

Many applied problems can be formulated as testing a hypothesis:

$H_{0}:\theta=0.$

The decision rule $\phi$ is a binary function of the data: $\phi=1$ if we
reject $H_{0}$; $\phi=0$ if we fail to reject $H_{0}$. The type one error rate
of the test is the probability of rejection if the null hypothesis holds. I
review the definition below.

###### Definition A1.7

When $H_{0}$ holds, define the type one error rate of the test $\phi$ as the
maximum possible value of the probability

$\textup{pr}(\phi=1).$

A standard choice is to make sure that the type one error rate is below
$\alpha=0.05$. The type two error rate of the test is the probability of no
rejection if the null hypothesis does not hold. I review the definition below.

###### Definition A1.8

When $H_{0}$ does not hold, define the type two error rate of the test $\phi$
as the probability

$\textup{pr}(\phi=0).$

Given the control of the type one error rate under $H_{0}$, we hope the type
two error rate is as low as possible when $H_{0}$ does not hold.

#####  153.4 Wald-type confidence interval and test

Many statistics problems have the following structure. The parameter of
interest is $\theta$. We first find a consistent estimator $\hat{\theta}$ that
converges in probability to $\theta$, and show that it is asymptotically
Normal with mean $\theta$ and variance $v$ which may depends on $\theta$ as
well as other parameters. We then find a consistent estimator $\hat{v}$ for
$v$, based on analytic formulas or the bootstrap reviewed in Chapter 156. We
finally construct the Wald-type confidence interval for $\theta$ as

$\hat{\theta}\pm z_{1-\alpha/2}\sqrt{\hat{v}}$

which covers $\theta$ with probability approximately $1-\alpha$. When this
interval excludes a particular $c$, for example, $c=0$, we reject the null
hypothesis $H_{0}(c):\theta=c$, which is called the Wald test.

#####  153.5 Duality between constructing confidence sets and testing null
hypotheses

Consider the statistical inference problem for a scalar parameter $\theta$. A
fundamental result in statistics is that constructing confidence sets for
$\theta$ is equivalent to testing null hypotheses about $\theta$. This is
often called the duality between constructing confidence sets and testing null
hypotheses.

Section 153.4 has reviewed the duality based on the Wald-type confidence
interval and test. The duality also holds in general. Assume that
$\hat{\Theta}$ is a $(1-\alpha)$-level confidence set for $\theta$:

$\textup{pr}(\theta\in\hat{\Theta})=1-\alpha.$

Then we can reject the null hypothesis $H_{0}(c):\theta=c$ if $c$ is not in
the set $\hat{\Theta}$. This is a valid test because when $\theta$ indeed
equals $c$, we have correct type one error rate
$\textup{pr}(\theta\not\in\hat{\Theta})=\alpha$. Conversely, if we test a
sequence of null hypotheses $H_{0}(c):\theta=c$, we can obtain the
corresponding $p$-values, $p(c)$, as a function of $c$. Then the values of $c$
that we fail to reject at level $\alpha$ form a confidence set for $\theta$:

$\hat{\Theta}=\\{c:p(c)\geq\alpha\\}=\\{c:\text{ fail to reject
}H_{0}(c)\text{ at level }\alpha\\}.$

It is a valid confidence set because

$\textup{pr}(\theta\in\hat{\Theta})=\textup{pr}\\{\text{fail to reject
}H_{0}(\theta)\text{ at level }\alpha\\}=1-\alpha.$

Here I use “confidence set” instead of “confidence interval” because
$\hat{\Theta}$ based on inverting tests may not be an interval. See the use of
the duality in Sections 155.2 and 15.1.

####  154 Inference with $2\times 2$ tables

#####  154.1 Fisher’s exact test

Fisher proposed an exact test for $H_{0}:p_{1}=p_{0}$ under the statistical
model:

$n_{11}\sim\text{Binomial}(n_{1},p_{1}),\quad
n_{01}\sim\text{Binomial}(n_{0},p_{0}),\quad
n_{11}\begin{picture}(9.0,8.0)\put(0.0,0.0){\line(1,0){9.0}}
\put(3.0,0.0){\line(0,1){8.0}} \put(6.0,0.0){\line(0,1){8.0}}
\end{picture}n_{01}.$

The table below summarizes the data.

| 1 | 0 | row sum
---|---|---|---
sample 1 | $n_{11}$ | $n_{10}$ | $n_{1}$
sample 0 | $n_{01}$ | $n_{00}$ | $n_{0}$
column sum | $n_{\cdot 1}$ | $n_{\cdot 0}$ | $n$

He argued that the sum $n_{11}+n_{01}\equiv n_{\cdot 1}$ contains little
information for the difference between $p_{1}$ and $p_{0}$, and $n_{11}$
conditioning on the sum has Hypergeometric distribution that does not depend
on the unknown parameter $p_{1}=p_{0}$ under $H_{0}$:

$\textup{pr}(n_{11}=k)=\frac{\binom{n_{\cdot 1}}{k}\binom{n-n_{\cdot
1}}{n_{1}-k}}{\binom{n}{n_{1}}}.$

In R, the function fisher.test implement this test.

#####  154.2 Estimation with $2\times 2$ tables

Based on the model in Section 154.1, we can estimate the parameters $p_{1}$
and $p_{0}$ by sample frequencies:

$\hat{p}_{1}=\frac{n_{11}}{n_{1}},\quad\hat{p}_{0}=\frac{n_{01}}{n_{0}}.$

Therefore, we can estimate the risk difference, log risk ratio, and log odds
ratio

rd $\displaystyle=$ $\displaystyle p_{1}-p_{0},$
$\displaystyle\log\textsc{rr}$ $\displaystyle=$
$\displaystyle\log\frac{p_{1}}{p_{0}},$ $\displaystyle\log\textsc{or}$
$\displaystyle=$ $\displaystyle\log\frac{p_{1}/(1-p_{1})}{p_{0}/(1-p_{0})}$

by the sample analogues

$\displaystyle\hat{\textsc{rd}}$ $\displaystyle=$
$\displaystyle\hat{p}_{1}-\hat{p}_{0},$ $\displaystyle\log\hat{\textsc{rr}}$
$\displaystyle=$ $\displaystyle\log\frac{\hat{p}_{1}}{\hat{p}_{0}},$
$\displaystyle\log\hat{\textsc{or}}$ $\displaystyle=$
$\displaystyle\log\frac{\hat{p}_{1}/(1-\hat{p}_{1})}{\hat{p}_{0}/(1-\hat{p}_{0})}=\log\frac{n_{11}n_{00}}{n_{10}n_{01}}.$

Based on the asymptotic approximation (see Problem A1.4), the estimated
variance for the above parameters are

$\displaystyle\frac{\hat{p}_{1}(1-\hat{p}_{1})}{n_{1}}+\frac{\hat{p}_{0}(1-\hat{p}_{0})}{n_{0}},$
$\displaystyle\frac{1-\hat{p}_{1}}{n_{1}\hat{p}_{1}}+\frac{1-\hat{p}_{0}}{n_{0}\hat{p}_{0}},$
$\displaystyle\frac{1}{n_{1}\hat{p}_{1}(1-\hat{p}_{1})}+\frac{1}{n_{0}\hat{p}_{0}(1-\hat{p}_{0})},$

respectively. The log transformation above yields better Normal approximations
because the risk ratio and odds ratio are always positive.

####  155 Two famous problems in statistics

#####  155.1 Behrens–Fisher problem

Consider the two-sample problem with $n_{1}$ units under the treatment and
$n_{0}$ units under the control, respectively. Assume the outcomes under the
treatment $\\{Y_{i}:Z_{i}=1\\}$ are IID from
$\textsc{N}(\mu_{1},\sigma_{1}^{2})$ and the outcomes under the control
$\\{Y_{i}:Z_{i}=0\\}$ are IID from $\textsc{N}(\mu_{0},\sigma_{0}^{2})$,
respectively. The goal is to test $H_{0}:\mu_{1}=\mu_{0}$.

Start with the easier case with $\sigma_{1}^{2}=\sigma_{0}^{2}.$ Coherent with
Chapter 3, let $\hat{\bar{Y}}(1)$ and $\hat{\bar{Y}}(0)$ denote the sample
means of the outcomes under the treatment and control, respectively. A
standard result is that

$t_{\textup{equal}}=\frac{\hat{\bar{Y}}(1)-\hat{\bar{Y}}(0)}{\sqrt{\frac{n}{n_{1}n_{0}(n-2)}\left[\sum_{Z_{i}=1}\\{Y_{i}-\hat{\bar{Y}}(1)\\}^{2}+\sum_{Z_{i}=0}\\{Y_{i}-\hat{\bar{Y}}(0)\\}^{2}\right]}}\sim
t_{n-2}.$

Based on $t_{\textup{equal}}$, we can construct a test for $H_{0}.$

Now consider the more difficult case with possibly different $\sigma_{1}^{2}$
and $\sigma_{0}^{2}.$ The distribution of $t_{\textup{equal}}$ is no longer
$t_{n-2}.$ Estimating the variances separately, we can also define

$t_{\textup{unequal}}=\frac{\hat{\bar{Y}}(1)-\hat{\bar{Y}}(0)}{\sqrt{\frac{\hat{S}^{2}(1)}{n_{1}}+\frac{\hat{S}^{2}(0)}{n_{0}}}},$

where

$\hat{S}^{2}(1)=(n_{1}-1)^{-1}\sum_{Z_{i}=1}\\{Y_{i}-\hat{\bar{Y}}(1)\\}^{2},\quad\hat{S}^{2}(0)=(n_{0}-1)^{-1}\sum_{Z_{i}=0}\\{Y_{i}-\hat{\bar{Y}}(0)\\}^{2}$

are the sample variances of the outcomes under the treatment and control,
respectively. Unfortunately, the exact distribution of $t_{\textup{unequal}}$
depends on the known variances. Testing $H_{0}$ without assuming equal
variances is the famous Behrens–Fisher problem. With large sample sizes
$n_{1}$ and $n_{0}$, the central limit theorem ensures that
$t_{\textup{unequal}}$ is approximately $\textsc{N}(0,1)$. So we can construct
approximate test for $H_{0}.$

#####  155.2 Fieller–Creasy problem

Consider the two-sample problem with $n_{1}$ units under the treatment and
$n_{0}$ units under the control, respectively. Assume the outcomes under the
treatment $\\{Y_{i}:Z_{i}=1\\}$ are IID from $\textsc{N}(\mu_{1},1)$ and the
outcomes under the control $\\{Y_{i}:Z_{i}=0\\}$ are IID from
$\textsc{N}(\mu_{0},1)$, respectively. The goal is to estimate
$\gamma=\mu_{1}/\mu_{0}$. We can use
$\hat{\gamma}=\hat{\bar{Y}}(1)/\hat{\bar{Y}}(0)$ to estimate $\gamma$. But the
point estimator has a complicated distribution, which does not yield a simple
procedure to construct the confidence interval for $\gamma$.

Fieller’s confidence interval can be formulated as inverting tests for a
sequence of null hypotheses: $H_{0}(c):\gamma=c$. Under $H_{0}(c)$, we have

$\frac{\hat{\bar{Y}}(1)-c\hat{\bar{Y}}(0)}{\sqrt{1/n_{1}+c^{2}/n_{0}}}\sim\textsc{N}(0,1)$

which motivates the confidence interval

$\left\\{c:\Big{|}\frac{\hat{\bar{Y}}(1)-c\hat{\bar{Y}}(0)}{\sqrt{1/n_{1}+c^{2}/n_{0}}}\Big{|}\leq
z_{\alpha}\right\\}$

where $z_{\alpha}$ is the upper $1-\alpha/2$ quantile of a standard Normal
random variable.

####  156 Bootstrap

It is often very tedious to derive the variance formulas for complex
estimators. Efron, (1979) proposed the bootstrap as a general tool for
variance estimation. There are many versions of the bootstrap (Davison and
Hinkley,, 1997). In this book, we only need the most basic one: the
nonparametric bootstrap, which will be simply called the bootstrap.

Consider the generic setting with

$Y_{1},\ldots,Y_{n}\stackrel{{\scriptstyle\text{IID}}}{{\sim}}Y,$

where $Y_{i}$ can be a general random element denoting the observed data for
unit $i$. An estimator $\hat{\theta}$ is a function of the observed data:
$\hat{\theta}=T(Y_{1},\ldots,Y_{n})$. When $T$ is a complex function, it may
not be easy to obtain the variance or asymptotic variance of $\hat{\theta}$.

The uncertainty of $\hat{\theta}$ is driven by the IID sampling of
$Y_{1},\ldots,Y_{n}$ from the true distribution. Although the true
distribution is unknown, it can be well approximated by its empirical version

$\hat{F}_{n}(y)=n^{-1}\sum_{i=1}^{n}I(Y_{i}\leq y),$

when the sample size $n$ is large. If we believe this approximation, we can
simulate $\hat{\theta}$ by sampling

$(Y_{1}^{*},\ldots,Y_{n}^{*})\stackrel{{\scriptstyle\text{IID}}}{{\sim}}\hat{F}_{n}(y).$

Because $\hat{F}_{n}(y)$ is a discrete distribution with mass $1/n$ on each
observed data point, the simulation of $\hat{\theta}$ reduces to the following
procedure:

  1. 1.

sample $(Y_{1}^{*},\ldots,Y_{n}^{*})$ from $\\{Y_{1},\ldots,Y_{n}\\}$ with
replacement;

  2. 2.

compute $\hat{\theta}^{*}=T(Y_{1}^{*},\ldots,Y_{n}^{*})$;

  3. 3.

repeat the above two steps $B$ times to obtain the bootstrap replicates
$\\{\hat{\theta}^{*}_{1},\ldots,\hat{\theta}^{*}_{B}\\}$.

We can then approximate the (asymptotic) variance of $\hat{\theta}$ by the
sample variance of the bootstrap replicates:

$\hat{V}_{\text{boot}}=(B-1)^{-1}\sum_{b=1}^{B}(\hat{\theta}^{*}_{b}-\bar{\theta}^{*})^{2},$

where $\bar{\theta}^{*}=B^{-1}\sum_{b=1}^{B}\hat{\theta}^{*}_{b}$. The
bootstrap confidence interval based on the Normal approximation is then

$\hat{\theta}\pm z_{1-\alpha/2}\sqrt{\hat{V}_{\text{boot}}},$

where $z_{1-\alpha/2}$ is the $1-\alpha/2$ upper quantile of
$\textsc{N}(0,1)$.

####  157 Homework problems

######  A1.1 Independent but not IID data

Assume that the $X_{i}$’s are independent with mean $\mu_{i}$ and variances
$\sigma_{i}^{2}$ for $i=1,\ldots,n$. The parameter of interest is
$\mu=n^{-1}\sum_{i=1}^{n}\mu_{i}.$ Show that
$\hat{\mu}=n^{-1}\sum_{i=1}^{n}X_{i}$ is unbiased for $\mu$ and find its
variance. Show that the usual variance estimator for IID data

$\hat{v}=\\{n(n-1)\\}^{-1}\sum_{i=1}^{n}(X_{i}-\hat{\mu})^{2}$

is a conservative estimator for the variance of $\hat{\mu}$ in the sense that

$E(\hat{v})-\textup{var}(\hat{\mu})=\\{n(n-1)\\}^{-1}\sum_{i=1}^{n}(\mu_{i}-\mu)^{2}\geq
0.$

Remark: Consider a simpler case with $\mu_{i}=\mu$ and
$\sigma_{i}^{2}=\sigma^{2}$ for all $i=1,\ldots,n$. The sample mean is
unbiased for $\mu$ with variance $\sigma^{2}/n$. Moreover, an unbiased
estimator for the variance $\sigma^{2}/n$ is $\hat{\sigma}^{2}/n=\hat{v}$,
where $\hat{\sigma}^{2}=(n-1)^{-1}\sum_{i=1}^{n}(X_{i}-\hat{\mu})^{2}$.

######  A1.2 Asymptotic Normality of ratio

Prove (136).

######  A1.3 Asymptotic Normality of product

Prove (138).

######  A1.4 Variance estimators in two-by-two tables

Use delta method to show the variance estimators in Section 154.2.

###  Chapter A2 Linear and Logistic Regressions

####  158 Population ordinary least squares

Assume that
$(x_{i},y_{i})_{i=1}^{n}\stackrel{{\scriptstyle\text{IID}}}{{\sim}}(x,y)$,
where $x$ is a $p$-dimensional random scalar or vector and $y$ is a random
scalar. Below I will use $(x,y)$ to denote a general observation, dropping the
subscript $i$ for simplicity. Define the population ordinary least squares
(OLS) coefficient as

$\beta=\arg\min_{b}E\left\\{(y-x^{\mkern-1.5mu\mathsf{T}}b)^{2}\right\\}.$

The objective function is quadratic in $b$, so we can show that the minimizer
is

$\beta=\left\\{E\left(xx^{\mkern-1.5mu\mathsf{T}}\right)\right\\}^{-1}E\left(xy\right)$

if the moments exist and $E\left(xx^{\mkern-1.5mu\mathsf{T}}\right)$ is
invertible.

With $\beta$, we can define

$\varepsilon=y-x^{\mkern-1.5mu\mathsf{T}}\beta$ (139)

as the population residual. By the definition of $\beta$, we can verify that

$E(x\varepsilon)=E\left\\{x(y-x^{\mkern-1.5mu\mathsf{T}}\beta)\right\\}=E(xy)-E(xx^{\mkern-1.5mu\mathsf{T}})\beta=0.$

######  Example A2.1 (population OLS with an intercept)

If we include $1$ as a component of $x$, then

$E(\varepsilon)=E(y-x^{\mkern-1.5mu\mathsf{T}}\beta)=0$

which further implies that $\textup{cov}(x,\varepsilon)=0$. So with an
intercept in $\beta$, the mean of the population residual must be zero, and it
is uncorrelated with other covariates by construction.

######  Example A2.2 (univariate population OLS with an intercept)

An important special case is that for scalars $x$ and $y$, we can define

$(\alpha,\beta)=\arg\min_{a,b}E\\{(y-a-bx)^{2}\\}$

which have explicit formulas

$\beta=\frac{\textup{cov}(x,y)}{\textup{var}(x)},\quad\alpha=E(y)-\beta E(x).$

######  Example A2.3 (univariate population OLS without an intercept)

Without intercept, we can define

$\gamma=\arg\min_{c}E\\{(y-cx)^{2}\\}$

which equals

$\gamma=\frac{E(xy)}{E(x^{2})}.$

When $x$ has mean zero, $\beta=\gamma$ in the above two population OLS.

We can also rewrite (139) as

$y=x^{\mkern-1.5mu\mathsf{T}}\beta+\varepsilon,$ (140)

which holds by the definition of the population OLS coefficient and residual
without any modeling assumption. We call (140) the population OLS
decomposition.

####  159 Sample OLS

Based on data
$(x_{i},y_{i})_{i=1}^{n}\stackrel{{\scriptstyle\text{IID}}}{{\sim}}(x,y)$, we
can easily obtain the moment estimator for the population OLS coefficient

$\hat{\beta}=\left(n^{-1}\sum_{i=1}^{n}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\right)^{-1}\left(n^{-1}\sum_{i=1}^{n}x_{i}y_{i}\right),$

and the residuals
$\hat{\varepsilon}_{i}=y_{i}-x_{i}^{\mkern-1.5mu\mathsf{T}}\hat{\beta}.$ This
is called the sample OLS or simply the OLS. The OLS coefficient $\hat{\beta}$
minimizes the residual sum of squares

$\hat{\beta}=\arg\min_{b}n^{-1}\sum_{i=1}^{n}(y_{i}-x_{i}^{\mkern-1.5mu\mathsf{T}}b)^{2},$

which satisfies the following Normal equation:

$\sum_{i=1}^{n}x_{i}(y_{i}-x_{i}^{\mkern-1.5mu\mathsf{T}}\hat{\beta})=0.$

The fitted values equal

$\hat{y}_{i}=x_{i}^{\mkern-1.5mu\mathsf{T}}\hat{\beta}\quad(i=1,\ldots,n).$

Using the matrix notation

$X=\begin{pmatrix}x_{1}^{\mkern-1.5mu\mathsf{T}}\\\ \vdots\\\
x_{n}^{\mkern-1.5mu\mathsf{T}}\end{pmatrix},\quad Y=\begin{pmatrix}y_{1}\\\
\vdots\\\ y_{n}\end{pmatrix},$

we can write the OLS coefficient as

$\hat{\beta}=(X^{\mkern-1.5mu\mathsf{T}}X)^{-1}X^{\mkern-1.5mu\mathsf{T}}Y$

and the fitted vector as

$\hat{Y}=X\hat{\beta}=X(X^{\mkern-1.5mu\mathsf{T}}X)^{-1}X^{\mkern-1.5mu\mathsf{T}}Y.$

Define the hat matrix as

$H=X(X^{\mkern-1.5mu\mathsf{T}}X)^{-1}X^{\mkern-1.5mu\mathsf{T}}.$

Then we also have $\hat{Y}=HY$, justifying the name “hat matrix.”

Assuming finite fourth moments of $(x,y)$, we can use the law of large numbers
and the central limit theorem to show that

$\sqrt{n}(\hat{\beta}-\beta)\rightarrow\text{N}(0,V=B^{-1}MB^{-1})$

in distribution, where $B=E(xx^{\mkern-1.5mu\mathsf{T}})$ and
$M=E(\varepsilon^{2}xx^{\mkern-1.5mu\mathsf{T}})$. So a moment estimator for
the asymptotic variance of $\hat{\beta}$ is

$\hat{V}_{\textsc{ehw}}=n^{-1}\left(n^{-1}\sum_{i=1}^{n}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\right)^{-1}\left(n^{-1}\sum_{i=1}^{n}\hat{\varepsilon}_{i}^{2}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\right)\left(n^{-1}\sum_{i=1}^{n}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\right)^{-1},$

which is called the Eicker–Huber–White (EHW) robust covariance estimator
(Eicker,, 1967; Huber,, 1967; White,, 1980). We can show that
$n\hat{V}_{\textsc{ehw}}\rightarrow V$ in probability. Based on $\hat{\beta}$
and $\hat{V}_{\textsc{ehw}}$, we can make inference about the population OLS
coefficient $\beta$. In R, the lm function can compute $\hat{\beta}$, and the
hccm function in the package car can compute $\hat{V}_{\textsc{ehw}}$.

There are many variants of the EHW robust covariance estimator (Long and
Ervin,, 2000). In particular, the HC1 variant modifies
$\hat{\varepsilon}_{i}^{2}$ to $\hat{\varepsilon}_{i}^{2}/(n-p)$, the HC2
variant modifies $\hat{\varepsilon}_{i}^{2}$ to
$\hat{\varepsilon}_{i}^{2}/(1-h_{ii})$, the HC3 variant modifies
$\hat{\varepsilon}_{i}^{2}$ to $\hat{\varepsilon}_{i}^{2}/(1-h_{ii})^{2}$, in
the definition of $\hat{V}_{\textsc{ehw}}$, where $h_{ii}$ is the $(i,i)$th
diagonal element of $H$, also called the leverage scores.

####  160 Frisch–Waugh–Lovell Theorem

The Frisch–Waugh–Lovell (FWL) theorem has two versions: one at the population
level and the other at the sample level. It reduces multivariate OLS to
univariate OLS and therefore facilitate the understanding and calculation of
the OLS coefficients. Below I will present special cases of the FWL Theorem
which are enough for this book.

######  Theorem A2.1 (population FWL)

The coefficient of $x_{1}$ in the OLS fit of $y$ on
$(x_{1},x_{2},\ldots,x_{p})$ equals the coefficient of $\tilde{x}_{1}$ in the
OLS fit of $y$ or $\tilde{y}$ on $\tilde{x}_{1}$, where $\tilde{y}$ is the
residual from the OLS fit of $y$ on $(x_{2},\ldots,x_{p})$ and $\tilde{x}_{1}$
is the residual from the OLS fit of $x_{1}$ on $(x_{2},\ldots,x_{p})$.

In Theorem A2.1, residualizing $x_{1}$ is crucial but residualizing $y$ is
not.

######  Theorem A2.2 (sample FWL)

With data $(Y,X_{1},X_{2},\ldots,X_{p})$ containing column vectors, the
coefficient of $X_{1}$ equals the coefficient of $\tilde{X}_{1}$ in the OLS
fit of $Y$ or $\tilde{Y}$ on $\tilde{X}_{1}$, where $\tilde{Y}$ is the
residual vector from the OLS fit of $Y$ on $(X_{2},\ldots,X_{p})$ and
$\tilde{X}_{1}$ is the residual from the OLS fit of $X_{1}$ on
$(X_{2},\ldots,X_{p})$.

Again, in Theorem A2.2, residualizing $X_{1}$ is crucial but residualizing $Y$
is not.

####  161 Linear model

Sometimes, we impose a stronger model assumption which requires the
conditional mean of $y$ given $x$ is linear:

$E(y\mid x)=x^{\mkern-1.5mu\mathsf{T}}\beta$

or, equivalently,

$y=x^{\mkern-1.5mu\mathsf{T}}\beta+\varepsilon,\qquad E(\varepsilon\mid x)=0,$

which is called the restricted mean model. Under this model, the population
OLS coefficient is the true parameter of interest:

$\displaystyle\left\\{E(xx^{\mkern-1.5mu\mathsf{T}})\right\\}^{-1}E(xy)$
$\displaystyle=\left\\{E(xx^{\mkern-1.5mu\mathsf{T}})\right\\}^{-1}E\left\\{xE(y\mid
x)\right\\}$
$\displaystyle=\left\\{E(xx^{\mkern-1.5mu\mathsf{T}})\right\\}^{-1}E(xx^{\mkern-1.5mu\mathsf{T}}\beta)$
$\displaystyle=\beta.$

Moreover, the population OLS coefficient does not depend on the distribution
of $x$. The asymptotic inference in Section 158 applies to this model too.

In the special case with $\textup{var}(\varepsilon\mid x)=\sigma^{2}$, the
asymptotic variance of the OLS coefficient reduces to

$V=\sigma^{2}\\{E(xx^{\mkern-1.5mu\mathsf{T}})\\}^{-1}$

so a simpler moment estimator for the asymptotic variance of $\hat{\beta}$ is

$\hat{V}_{\textsc{ols}}=\hat{\sigma}^{2}\left(\sum_{i=1}^{n}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\right)^{-1}$

where $\hat{\sigma}^{2}=(n-p)^{-1}\sum_{i=1}^{n}\hat{\varepsilon}_{i}^{2}$.
This is the standard covariance estimator from the lm function.

####  162 Weighted least squares

Assuming that
$(w_{i},x_{i},y_{i})\stackrel{{\scriptstyle\text{IID}}}{{\sim}}(w,x,y)$ with
$w\neq 0$. At the population level, we can define weighted least squares (WLS)
coefficient as

$\beta_{w}=\arg\min_{b}E\\{w(y-x^{\mkern-1.5mu\mathsf{T}}b)^{2}\\},$

which satisfies

$E\\{wx(y-x^{\mkern-1.5mu\mathsf{T}}\beta_{w})\\}=0$

and thus equals

$\beta_{w}=\\{E(wxx^{\mkern-1.5mu\mathsf{T}})\\}^{-1}E(wxy)$

if $E(wxx^{\mkern-1.5mu\mathsf{T}})$ is invertible.

At the sample level, we can define the WLS coefficient as

$\hat{\beta}_{w}=\arg\min_{b}\sum_{i=1}^{n}w_{i}(y_{i}-x_{i}^{\mkern-1.5mu\mathsf{T}}b)^{2},$

which satisfies

$\sum_{i=1}^{n}w_{i}x_{i}(y_{i}-x_{i}^{\mkern-1.5mu\mathsf{T}}\hat{\beta}_{w})=0$

and thus equals

$\hat{\beta}_{w}=\left(n^{-1}\sum_{i=1}^{n}w_{i}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\right)^{-1}\left(n^{-1}\sum_{i=1}^{n}w_{i}x_{i}y_{i}\right)$

if $\sum_{i=1}^{n}w_{i}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}$ is invertible.

####  163 Logistic regression

#####  163.1 Model

Technically, we can use apply the OLS procedure even if the outcome $y$ is
binary. However, it is a little awkward to have predicted probabilities
outside the range of $[0,1]$. This motivates us to consider the following
model:

$\textup{pr}(y_{i}=1\mid x_{i})=g(x_{i}^{\mkern-1.5mu\mathsf{T}}\beta),$

where $g(\cdot):\mathbb{R}\rightarrow[0,1]$ is a monotone function, and its
inverse is often called the link function. The $g(\cdot)$ function can be any
distribution function of a random variable, but we will focus on the logistic
form:

$g(z)=\frac{e^{z}}{1+e^{z}}=(1+e^{-z})^{-1}.$

We can also write the logistic model as

$\textup{pr}(y_{i}=1\mid
x_{i})\equiv\pi(x_{i},\beta)=\frac{e^{x_{i}^{\mkern-1.5mu\mathsf{T}}\beta}}{1+e^{x_{i}^{\mkern-1.5mu\mathsf{T}}\beta}},$

or, equivalently,

$\text{logit}\left\\{\textup{pr}(y_{i}=1\mid
x_{i})\right\\}\equiv\log\frac{\textup{pr}(y_{i}=1\mid
x_{i})}{1-\textup{pr}(y_{i}=1\mid
x_{i})}=x_{i}^{\mkern-1.5mu\mathsf{T}}\beta.$

Assume that $x_{i1}$ is binary. Under the logistic model, we have

$\displaystyle\beta_{1}$ $\displaystyle=$
$\displaystyle\text{logit}\left\\{\textup{pr}(y_{i}=1\mid
x_{i1}=1,\ldots)\right\\}-\text{logit}\left\\{\textup{pr}(y_{i}=1\mid
x_{i1}=0,\ldots)\right\\}$ $\displaystyle=$
$\displaystyle\log\frac{\textup{pr}(y_{i}=1\mid
x_{i1}=1,\ldots)/\textup{pr}(y_{i}=0\mid
x_{i1}=1,\ldots)}{\textup{pr}(y_{i}=1\mid
x_{i1}=0,\ldots)/\textup{pr}(y_{i}=0\mid x_{i1}=0,\ldots)},$

where $\ldots$ contains all other regressor $x_{i2},\ldots,x_{ip}$. Therefore,
the coefficient $\beta_{1}$ equals the log odds ratio of $x_{i1}$ on $y_{i}$
conditional on other regressors.

#####  163.2 Maximum likelihood estimate

To estimate the parameter $\beta$, we can maximize the following likelihood
function:

$\displaystyle L(\beta)$
$\displaystyle=\prod_{i=1}^{n}\left\\{\pi(x_{i},\beta)\right\\}^{y_{i}}\left\\{1-\pi(x_{i},\beta)\right\\}^{1-y_{i}}$
$\displaystyle=\prod_{i=1}^{n}\left\\{\frac{\pi(x_{i},\beta)}{1-\pi(x_{i},\beta)}\right\\}^{y_{i}}\left\\{1-\pi(x_{i},\beta)\right\\}$
$\displaystyle=\prod_{i=1}^{n}\left(e^{x_{i}^{\mkern-1.5mu\mathsf{T}}\beta}\right)^{y_{i}}\frac{1}{1+e^{x_{i}^{\mkern-1.5mu\mathsf{T}}\beta}}$
$\displaystyle=\prod_{i=1}^{n}\frac{e^{y_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}\beta}}{1+e^{x_{i}^{\mkern-1.5mu\mathsf{T}}\beta}}.$

Let $\hat{\beta}$ denote the maximizer, which is called the maximum likelihood
estimate (MLE). Taking the log of $L(\beta)$ and differentiating it with
respect to $\beta$, we can show that the MLE must satisfy the first order
condition:

$\sum_{i=1}^{n}x_{i}\\{y_{i}-\pi(x_{i},\hat{\beta})\\}=0.$

So if $x_{i}$ contains an intercept, the MLE must satisfy

$\sum_{i=1}^{n}\\{y_{i}-\pi(x_{i},\hat{\beta})\\}=0,$

that is, the average of the observed $y_{i}$’s must be identical to the
average of the fitted probabilities $\pi(x_{i},\hat{\beta})$’s.

Using the general theory for the MLE, we can show that it is consistent for
the true parameter $\beta$ and is asymptotically normal:

$\sqrt{n}(\hat{\beta}-\beta)\rightarrow\text{N}(0,V)$

in distribution, where
$V=E\left[\pi(x_{i},\beta)\\{1-\pi(x_{i},\beta)\\}xx^{\mkern-1.5mu\mathsf{T}}\right].$
So we can approximate the covariance matrix of $\hat{\beta}$ by

$n^{-1}\sum_{i=1}^{n}\pi(x_{i},\hat{\beta})\\{1-\pi(x_{i},\hat{\beta})\\}x_{i}x_{i}^{\mkern-1.5mu\mathsf{T}}.$

In R, the glm function can find the MLE and report the estimated covariance
matrix.

#####  163.3 Extension to the case-control study

In case-control studies, sampling is conditional on the binary outcome, that
is, units with outcomes $y_{i}=1$ and $y_{i}=0$ are sampled with different
probabilities. Let $s_{i}$ be the sampling indicator. In case control studies,
we have

$\textup{pr}(s_{i}=1\mid x_{i},y_{i})=\textup{pr}(s_{i}=1\mid y_{i})$

as a function of $y_{i}$, and we only observe units with $s_{i}=1$.

Prentice and Pyke, (1979) showed that logistic regression is applicable in
case-control studies although the above discussion assume IID sampling.

#####  163.4 Logistic regression with weights

Sometimes, unit $i$ has weight $w_{i}$, then we can fit a weighted logistic
regression by solving

$\sum_{i=1}^{n}w_{i}x_{i}\\{y_{i}-\pi(x_{i},\hat{\beta})\\}=0.$

####  164 Homework problems

######  A2.1 Sample OLS with intercept

Assume the regressor $x_{i}$ contains an intercept. Show that

$\bar{y}=\bar{x}^{\mkern-1.5mu\mathsf{T}}\hat{\beta}.$ (141)

######  A2.2 Univariate weighed least squares

As a special case of WLS, define

$(\hat{\alpha}_{w},\hat{\beta}_{w})=\arg\min_{(a,b)}\sum_{i=1}^{n}w_{i}(y_{i}-a-bx_{i})^{2}$

where $w_{i}\geq 0$. Show that

$\displaystyle\hat{\beta}_{w}=\frac{\sum_{i=1}^{n}w_{i}(x_{i}-\bar{x}_{w})(y_{i}-\bar{y}_{w})}{\sum_{i=1}^{n}w_{i}(x_{i}-\bar{x}_{w})^{2}}$
(142)

and

$\displaystyle\hat{\alpha}_{w}=\bar{y}_{w}-\hat{\beta}_{w}\bar{x}_{w},$ (143)

where $\bar{x}_{w}=\sum_{i=1}^{n}w_{i}x_{i}/\sum_{i=1}^{n}w_{i}$ and
$\bar{y}_{w}=\sum_{i=1}^{n}w_{i}y_{i}/\sum_{i=1}^{n}w_{i}$ are the weighted
averages of the $x_{i}$’s and $y_{i}$’s.

Further assume that the $x_{i}$’s are binary. Show that

$\hat{\beta}_{w}=\frac{\sum_{i=1}^{n}w_{i}x_{i}y_{i}}{\sum_{i=1}^{n}w_{i}x_{i}}-\frac{\sum_{i=1}^{n}w_{i}(1-x_{i})y_{i}}{\sum_{i=1}^{n}w_{i}(1-x_{i})}.$

That is, if the regressor is binary in the univariate WLS, the coefficient of
the regressor equals the difference in the weighted means.

Hint: Think about an appropriate reparametrization of the WLS problem.
Otherwise, the derivation can be tedious.

######  A2.3 OLS with orthogonal regressors

Consider sample OLS fit of an $n$-vector $Y$ on an $n\times p$ matrix $X$,
with coefficient $\hat{\beta}$. Partition $X$ into $X=(X_{1},X_{2})$, where
$X_{1}$ is an $n\times k$ matrix and $X_{2}$ is an $n\times l$ matrix, with
$p=k+l$. Correspondingly, partition $\hat{\beta}$ into

$\hat{\beta}=\begin{pmatrix}\hat{\beta}_{1}\\\ \hat{\beta}_{2}\end{pmatrix}.$

Assume $X_{1}$ and $X_{2}$ are orthogonal, that is,
$X_{1}^{\mkern-1.5mu\mathsf{T}}X_{2}=0$. Show that $\hat{\beta}_{1}$ equals
the coefficient from OLS of $Y$ on $X_{1}$ and $\hat{\beta}_{2}$ equals the
coefficient from OLS of $Y$ on $X_{2}$, respectively.

######  A2.4 OLS with a non-degenerate transformation of the regressors

Define $\hat{\beta}$ as the coefficient from the sample OLS fit of an
$n$-vector $Y$ on an $n\times p$ matrix $X$. Let $\Gamma$ be a $p\times p$
non-degenerate matrix, and define $X^{\prime}=X\Gamma$. Define
$\hat{\beta}^{\prime}$ as the coefficient from the sample OLS fit of $Y$ on
$X^{\prime}$.

Show that

$\hat{\beta}=\Gamma\hat{\beta}^{\prime}.$

###  Chapter A3 Some Useful Lemmas for Simple Random Sampling

Simple Random Sampling

####  165 Lemmas

Simple random sampling is a basic topic in standard survey sampling textbooks
(e.g., Cochran,, 1953). Below I review some results for simple random sampling
that are useful for design-based inference in the CRE in Chapters 3 and 4.

A simple random sample of size $n_{1}$ consists of a subset from a finite
population of $n$ units indexed by $i=1,\ldots,n$. Let
$\bm{Z}=(Z_{1},\ldots,Z_{n})$ be the inclusion indicators of the $n$ units
with $Z_{i}=1$ if unit $i$ is sampled and $Z_{i}=0$ otherwise. The vector
$\bm{Z}$ can take $\binom{n}{n_{1}}$ possible permutations of a vector of
$n_{1}$ 1’s and $n_{0}$ 0’s, and each has equal probability. The following
lemma summarizes the first two moments of the inclusion indicators.

###### Lemma A3.1

Under simple random sampling, we have

$E(Z_{i})=\frac{n_{1}}{n},\quad\textup{var}(Z_{i})=\frac{n_{1}n_{0}}{n^{2}},\quad\textup{cov}(Z_{i},Z_{j})=-\frac{n_{1}n_{0}}{n^{2}(n-1)}.$

In more compact forms, we have

$E(\bm{Z})=\frac{n_{1}}{n}\bm{1}_{n},\quad\textup{cov}(\bm{Z})=\frac{n_{1}n_{0}}{n(n-1)}\bm{P}_{n},$

where $\bm{1}_{n}$ is a $n$-dimensional vector of $1$’s, and
$\bm{P}_{n}=\bm{I}_{n}-n^{-1}\bm{1}_{n}\bm{1}_{n}^{\top}$ is the $n\times n$
projection matrix orthogonal to $\bm{1}_{n}.$

Let $\\{c_{1},\ldots,c_{n}\\}$ be a finite population with mean
$\bar{c}=\sum_{i=1}^{n}c_{i}/n$ and variance

$S_{c}^{2}=(n-1)^{-1}\sum_{i=1}^{n}(c_{i}-\bar{c})^{2};$

let $\\{d_{1},\ldots,d_{n}\\}$ be another finite population with mean
$\bar{d}=\sum_{i=1}^{n}d_{i}/n$ and variance

$S_{d}^{2}=(n-1)^{-1}\sum_{i=1}^{n}(d_{i}-\bar{d})^{2};$

their covariance is

$S_{cd}=(n-1)^{-1}\sum_{i=1}^{n}(c_{i}-\bar{c})(d_{i}-\bar{d}).$

Based on the simple random sample, the sample means are

$\hat{\bar{c}}=n_{1}^{-1}\sum_{i=1}^{n}Z_{i}c_{i},\quad\hat{\bar{d}}=n_{1}^{-1}\sum_{i=1}^{n}Z_{i}d_{i};$

sample variances are

$\hat{S}_{c}^{2}=(n_{1}-1)^{-1}\sum_{i=1}^{n}Z_{i}(c_{i}-\hat{\bar{c}})^{2},\quad\hat{S}_{d}^{2}=(n_{1}-1)^{-1}\sum_{i=1}^{n}Z_{i}(d_{i}-\hat{\bar{d}})^{2};$

the sample covariance is

$\hat{S}_{cd}=(n_{1}-1)^{-1}\sum_{i=1}^{n}Z_{i}(c_{i}-\hat{\bar{c}})(d_{i}-\hat{\bar{d}}).$

Lemma A3.2 below gives the moments of the sample means $\hat{\bar{c}}$ and
$\hat{\bar{d}}$.

###### Lemma A3.2

The sample means are unbiased for the population means:

$E(\hat{\bar{c}})=\bar{c},\quad E(\hat{\bar{d}})=\bar{d}.$

Their variances and covariance are

$\displaystyle\textup{var}\left(\hat{\bar{c}}\right)=\frac{n_{0}}{nn_{1}}S_{c}^{2},\quad\textup{var}\left(\hat{\bar{d}}\right)=\frac{n_{0}}{nn_{1}}S_{d}^{2},\quad\textup{cov}\left(\hat{\bar{c}},\hat{\bar{d}}\right)=\frac{n_{0}}{nn_{1}}S_{cd}.$

In the variance formula in Lemma A3.2, the coefficient
$n_{0}/(nn_{1})=1/n_{1}\times(1-n_{1}/n)$ in Lemma A3.2 is different from
$1/n_{1}$ under IID sampling. The additional factor $1-n_{1}/n=n_{0}/n$ is
called the finite population correction.

Lemma A3.3 below gives the unbiasedness of the sample variances and covariance
for estimating the population analogs.

###### Lemma A3.3

The sample variances and covariance are unbiased for their population
versions:

$E(\hat{S}_{c}^{2})=S_{c}^{2},\quad E(\hat{S}_{d}^{2})=S_{d}^{2},\quad
E(\hat{S}_{cd})=S_{cd}.$

An important practical question is to make inference about $\bar{c}$ based on
the simple random sample. This requires a more precise characterization of the
distribution of its unbiased estimator $\hat{\bar{c}}$. The finite-sample
exact distribution of $\hat{\bar{c}}$ depends on the whole finite population
$\\{c_{1},\ldots,c_{n}\\}$, which is intractable in general. The following
finite population central limit theorem characterizes the asymptotic
distribution of $\hat{\bar{c}}$ based on its first two moments.

######  Lemma A3.4 (finite population central limit theorem)

As $n\rightarrow\infty$, if

$\frac{\max_{1\leq i\leq
n}(c_{i}-\bar{c})^{2}}{\min(n_{1},n_{0})S_{c}^{2}}\rightarrow 0,$

then

$\frac{\hat{\bar{c}}-\bar{c}}{\sqrt{\frac{n_{0}}{nn_{1}}S_{c}^{2}}}\rightarrow\textsc{N}(0,1)$

in distribution, and $\hat{S}_{c}^{2}/S_{c}^{2}\rightarrow 1$ in probability.

Lemma A3.4 justifies the Wald-type $1-\alpha$ confidence interval for
$\bar{c}$:

$\hat{\bar{c}}\pm z_{1-\alpha/2}\sqrt{\frac{n_{0}}{nn_{1}}\hat{S}_{c}^{2}}$

where $z_{1-\alpha/2}$ is the $1-\alpha/2$ upper quantile of the standard
Normal random variable.

####  166 Proofs

Proof of Lemma A3.1: By symmetry, the $Z_{i}$’s have the same mean, so

$n_{1}=\sum_{i=1}^{n}Z_{i}=E\left(\sum_{i=1}^{n}Z_{i}\right)=nE(Z_{i})\Longrightarrow
E(Z_{i})=n_{1}/n.$

Because $Z_{i}$ is a Bernoulli random variable, its variance is

$\textup{var}(Z_{i})=\frac{n_{1}}{n}\left(1-\frac{n_{1}}{n}\right)=\frac{n_{1}n_{0}}{n^{2}}.$

By symmetry again, the $Z_{i}$’s have the same variance and the pairs
$(Z_{i},Z_{j})$’s have the same covariance, so

$0=\textup{var}\left(\sum_{i=1}^{n}Z_{i}\right)=n\textup{var}(Z_{i})+n(n-1)\textup{cov}(Z_{i},Z_{j})$

which implies that

$\textup{cov}(Z_{i},Z_{j})=-\frac{n_{1}n_{0}}{n^{2}(n-1)}\quad(i\neq j).$

$\square$

Proof of Lemma A3.1: The unbiasedness of the sample mean follows from
linearity. For example,

$E(\hat{\bar{c}})=E\left(\frac{1}{n_{1}}\sum_{i=1}^{n}Z_{i}c_{i}\right)=\frac{1}{n_{1}}\sum_{i=1}^{n}E(Z_{i})c_{i}=\bar{c}.$

The covariance of the sample means is

$\displaystyle\textup{cov}(\hat{\bar{c}},\hat{\bar{d}})$ $\displaystyle=$
$\displaystyle\textup{cov}\left\\{\frac{1}{n_{1}}\sum_{i=1}^{n}Z_{i}(c_{i}-\bar{c}),\frac{1}{n_{1}}\sum_{i=1}^{n}Z_{i}(d_{i}-\bar{d})\right\\}$
$\displaystyle=$
$\displaystyle\frac{1}{n_{1}^{2}}\left[\sum_{i=1}^{n}\textup{var}(Z_{i})(c_{i}-\bar{c})(d_{i}-\bar{d})+\sum_{i\neq
j}\textup{cov}(Z_{i},Z_{j})(c_{i}-\bar{c})(d_{j}-\bar{d})\right]$
$\displaystyle=$
$\displaystyle\frac{1}{n_{1}^{2}}\left[\frac{n_{1}n_{0}}{n^{2}}\sum_{i=1}^{n}(c_{i}-\bar{c})(d_{i}-\bar{d})-\frac{n_{1}n_{0}}{n^{2}(n-1)}\sum_{i\neq
j}(c_{i}-\bar{c})(d_{j}-\bar{d})\right].$

Because

$0=\sum_{i=1}^{n}(c_{i}-\bar{c})\sum_{i=1}^{n}(d_{i}-\bar{d})=\sum_{i=1}^{n}(c_{i}-\bar{c})(d_{i}-\bar{d})+\sum_{i\neq
j}(c_{i}-\bar{c})(d_{j}-\bar{d}),$

the covariance of the sample means reduces to

$\displaystyle\textup{cov}(\hat{\bar{c}},\hat{\bar{d}})$ $\displaystyle=$
$\displaystyle\frac{1}{n_{1}^{2}}\left[\frac{n_{1}n_{0}}{n^{2}}\sum_{i=1}^{n}(c_{i}-\bar{c})(d_{i}-\bar{d})+\frac{n_{1}n_{0}}{n^{2}(n-1)}\sum_{i=1}^{n}(c_{i}-\bar{c})(d_{i}-\bar{c})\right]$
$\displaystyle=$ $\displaystyle\frac{n_{0}}{nn_{1}}S_{cd}.$

The variance formulas are just special cases with
$\hat{\bar{c}}=\hat{\bar{d}}$. $\square$

Proof of Lemma A3.3: We prove only the sample covariance term, because the
formulas for sample variances are special cases. We have the following
decomposition:

$\displaystyle(n_{1}-1)\hat{S}_{cd}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{n}Z_{i}(c_{i}-\hat{\bar{c}})(d_{i}-\hat{\bar{d}})$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{n}Z_{i}\\{(c_{i}-\bar{c})-(\hat{\bar{c}}-\bar{c})\\}\\{(d_{i}-\bar{d})-(\hat{\bar{d}}-\bar{d})\\}$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{n}Z_{i}(c_{i}-\bar{c})(d_{i}-\bar{d})-n_{1}(\hat{\bar{c}}-\bar{c})(\hat{\bar{d}}-\bar{d}).$

Taking expectation on both sides, we have

$\displaystyle E\\{(n_{1}-1)\hat{S}_{cd}\\}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{n}E(Z_{i})(c_{i}-\bar{c})(d_{i}-\bar{d})-n_{1}E\\{(\hat{\bar{c}}-\bar{c})(\hat{\bar{d}}-\bar{d})\\}$
$\displaystyle=$
$\displaystyle\frac{n_{1}}{n}\sum_{i=1}^{n}(c_{i}-\bar{c})(d_{i}-\bar{d})-n_{1}\frac{n_{0}}{nn_{1}}S_{cd}$
$\displaystyle=$ $\displaystyle
S_{cd}\left\\{\frac{n_{1}(n-1)}{n}-\frac{n_{0}}{n}\right\\}$ $\displaystyle=$
$\displaystyle(n_{1}-1)S_{cd},$

and the conclusion follows by dividing both sides by $n_{1}-1$. $\square$

Proof of Lemma A3.4: Hájek, (1960) gave a proof of the central limit theorem
for simple random sampling, and Lehmann, (1975) gave a more accessible version
of the proof. Li and Ding, (2017) modified the central limit theorem as
presented in Lemma A3.4, and proved the consistency of the sample variance.
Due to the technical complexities, I omit the proof. $\square$

####  167 Comments on the literature

Survey sampling and experimental design are deeply connected ever since
Neyman, (1934, 1935)’s seminal work. Li and Ding, (2017) and Mukerjee et al.,
(2018) made many theoretical ties between these two areas.

####  168 Homework Problems

######  A3.1 Vector form of the results

Assume the $c_{i}$’s are vectors and modify

$S_{c}^{2}=(n-1)^{-1}\sum_{i=1}^{n}(c_{i}-\bar{c})(c_{i}-\bar{c})^{\mkern-1.5mu\mathsf{T}},\quad\hat{S}_{c}^{2}=(n_{1}-1)^{-1}\sum_{i=1}^{n}Z_{i}(c_{i}-\hat{\bar{c}})(c_{i}-\hat{\bar{c}})^{\mkern-1.5mu\mathsf{T}}.$

Show that

$E(\hat{\bar{c}})=\bar{c},\quad\textup{cov}(\hat{\bar{c}})=\frac{n_{0}}{nn_{1}}S_{c}^{2},\quad
E(\hat{S}_{c}^{2})=S_{c}^{2}.$

### References

  * Abadie and Imbens, (2006) Abadie, A. and Imbens, G. W. (2006).  Large sample properties of matching estimators for average treatment effects.  Econometrica, 74:235–267.
  * Abadie and Imbens, (2008) Abadie, A. and Imbens, G. W. (2008).  On the failure of the bootstrap for matching estimators.  Econometrica, 76:1537–1557.
  * Abadie and Imbens, (2011) Abadie, A. and Imbens, G. W. (2011).  Bias-corrected matching estimators for average treatment effects.  Journal of Business and Economic Statistics, 29:1–11.
  * Abadie and Imbens, (2016) Abadie, A. and Imbens, G. W. (2016).  Matching on the estimated propensity score.  Econometrica, 84:781–807.
  * Alwin and Hauser, (1975) Alwin, D. F. and Hauser, R. M. (1975).  The decomposition of effects in path analysis.  American Sociological Review, 40:37–47.
  * Amarante et al., (2016) Amarante, V., Manacorda, M., Miguel, E., and Vigorito, A. (2016).  Do cash transfers improve birth outcomes? evidence from matched vital statistics, program, and social security data.  American Economic Journal: Economic Policy, 8:1–43.
  * Anderson and Rubin, (1950) Anderson, T. W. and Rubin, H. (1950).  The asymptotic properties of estimates of the parameters of a single equation in a complete system of stochastic equations.  Annals of Mathematical Statistics, 21:570–582.
  * Angrist et al., (2009) Angrist, J., Lang, D., and Oreopoulos, P. (2009).  Incentives and services for college achievement: Evidence from a randomized trial.  American Economic Journal: Applied Economics, 1:136–163.
  * Angrist and Lavy, (2009) Angrist, J. and Lavy, V. (2009).  The effects of high stakes high school achievement awards: Evidence from a randomized trial.  American Economic Review, 99:1384–1414.
  * Angrist, (1990) Angrist, J. D. (1990).  Lifetime earnings and the Vietnam era draft lottery: evidence from social security administrative records.  American Economic Review, 80:313–336.
  * Angrist, (1998) Angrist, J. D. (1998).  Estimating the labor market impact of voluntary military service using social security data on military applicants.  Econometrica, 66:249–288.
  * Angrist and Evans, (1998) Angrist, J. D. and Evans, W. N. (1998).  Children and their parents’ labor supply: Evidence from exogenous variation in family size.  American Economic Review, 88:450–477.
  * Angrist and Imbens, (1995) Angrist, J. D. and Imbens, G. W. (1995).  Two-stage least squares estimation of average causal effects in models with variable treatment intensity.  Journal of the American Statistical Association, 90:431–442.
  * Angrist et al., (1996) Angrist, J. D., Imbens, G. W., and Rubin, D. B. (1996).  Identification of causal effects using instrumental variables (with discussion).  Journal of the American Statistical Association, 91:444–455.
  * Angrist and Krueger, (1991) Angrist, J. D. and Krueger, A. B. (1991).  Does compulsory school attendance affect schooling and earnings?  Quarterly Journal of Economics, 106:979–1014.
  * Angrist and Pischke, (2008) Angrist, J. D. and Pischke, J.-S. (2008).  Mostly Harmless Econometrics: An Empiricist’s Companion.  Princeton: Princeton University Press.
  * Angrist and Pischke, (2014) Angrist, J. D. and Pischke, J.-S. (2014).  Mastering’Metrics: The Path from Cause to Effect.  Princeton: Princeton University Press.
  * Aronow et al., (2014) Aronow, P. M., Green, D. P., and Lee, D. K. K. (2014).  Sharp bounds on the variance in randomized experiments.  Annals of Statistics, 42:850–871.
  * Asher and Novosad, (2020) Asher, S. and Novosad, P. (2020).  Rural roads and local economic development.  American Economic Review, 110:797–823.
  * Baker and Lindeman, (1994) Baker, S. G. and Lindeman, K. S. (1994).  The paired availability design: a proposal for evaluating epidural analgesia during labor.  Statistics in Medicine, 13:2269–2278.
  * Balke and Pearl, (1997) Balke, A. and Pearl, J. (1997).  Bounds on treatment effects from studies with imperfect compliance.  Journal of the American Statistical Association, 92:1171–1176.
  * Ball et al., (1973) Ball, S., Bogatz, G., Rubin, D., and Beaton, A. (1973).  Reading with television: An evaluation of the electric company. a report to the children’s television workshop. volumes 1 and 2.
  * Bang and Robins, (2005) Bang, H. and Robins, J. M. (2005).  Doubly robust estimation in missing data and causal inference models.  Biometrics, 61:962–973.
  * Barnard, (1947) Barnard, G. A. (1947).  Significance tests for $2\times 2$ tables.  Biometrika, 34:123–138.
  * Baron and Kenny, (1986) Baron, R. M. and Kenny, D. A. (1986).  The moderator-mediator variable distinction in social psychological research: Conceptual, strategic, and statistical considerations.  Journal of Personality and Social Psychology, 51:1173–1182.
  * Basmann, (1957) Basmann, R. L. (1957).  A generalized classical method of linear estimation of coefficients in a structural equation.  Econometrica, 25:77–83.
  * Bazzano et al., (2003) Bazzano, L. A., He, J., Muntner, P., Vupputuri, S., and Whelton, P. K. (2003).  Relationship between cigarette smoking and novel risk factors for cardiovascular disease in the United States.  Annals of Internal Medicine, 138:891–897.
  * Berk et al., (2013) Berk, R., Pitkin, E., Brown, L., Buja, A., George, E., and Zhao, L. (2013).  Covariance adjustments for the analysis of randomized field experiments.  Evaluation Review, 37:170–196.
  * Bertrand and Mullainathan, (2004) Bertrand, M. and Mullainathan, S. (2004).  Are Emily and Greg more employable than Lakisha and Jamal? A field experiment on labor market discrimination.  American Economic Review, 94:991–1013.
  * Bickel et al., (1975) Bickel, P. J., Hammel, E. A., and O’Connell, J. W. (1975).  Sex bias in graduate admissions: Data from Berkeley.  Science, 187:398–404.
  * Bickel et al., (1993) Bickel, P. J., Klaassen, C. A. J., Ritov, Y., and Wellner, J. A. (1993).  Efficient and Adaptive Estimation for Semiparametric Models.  Baltimore: Johns Hopkins University Press.
  * Bind and Rubin, (2020) Bind, M.-A. C. and Rubin, D. B. (2020).  When possible, report a fisher-exact p value and display its underlying null randomization distribution.  Proceedings of the National Academy of Sciences of the United States of America, 117:19151–19158.
  * Blackwell, (2013) Blackwell, M. (2013).  A framework for dynamic causal inference in political science.  American Journal of Political Science, 57:504–520.
  * Bloniarz et al., (2016) Bloniarz, A., Liu, H., Zhang, C. H., Sekhon, J., and Yu, B. (2016).  Lasso adjustments of treatment effect estimates in randomized experiments.  Proceedings of the National Academy of Sciences of the United States of America, 113:7383–7390.
  * Bloom, (1984) Bloom, H. S. (1984).  Accounting for no-shows in experimental evaluation designs.  Evaluation Review, 8:225–246.
  * Bor et al., (2014) Bor, J., Moscoe, E., Mutevedzi, P., Newell, M.-L., and Bärnighausen, T. (2014).  Regression discontinuity designs in epidemiology: causal inference without randomized trials.  Epidemiology, 25:729.
  * Bowden et al., (2015) Bowden, J., Davey Smith, G., and Burgess, S. (2015).  Mendelian randomization with invalid instruments: effect estimation and bias detection through Egger regression.  International Journal of Epidemiology, 44:512–525.
  * Bowden et al., (2018) Bowden, J., Spiller, W., Del Greco M, F., Sheehan, N., Thompson, J., Minelli, C., and Davey Smith, G. (2018).  Improving the visualization, interpretation and analysis of two-sample summary data mendelian randomization via the radial plot and radial regression.  International Journal of Epidemiology, 47:1264–1278.
  * Bradford Hill, (1965) Bradford Hill, A. (1965).  The environment and disease: association or causation?  Proceedings of the Royal Society of Medicine, 58:295–300.
  * Bradford Hill, (2020) Bradford Hill, A. (2020).  The environment and disease: association or causation? (with discussion).  Observational Studies, 6:1–65.
  * Bruhn and McKenzie, (2009) Bruhn, M. and McKenzie, D. (2009).  In pursuit of balance: Randomization in practice in development field experiments.  American Economic Journal: Applied Economics, 1:200–232.
  * Butler, (1969) Butler, C. C. (1969).  A test for symmetry using the sample distribution function.  Annals of Mathematical Statistics, 40:2209–2210.
  * Cao et al., (2009) Cao, W., Tsiatis, A. A., and Davidian, M. (2009).  Improving efficiency and robustness of the doubly robust estimator for a population mean with incomplete data.  Biometrika, 96:723–734.
  * Card, (1993) Card, D. (1993).  Using geographic variation in college proximity to estimate the return to schooling.  Technical report, National Bureau of Economic Research.
  * Carpenter and Dobkin, (2009) Carpenter, C. and Dobkin, C. (2009).  The effect of alcohol consumption on mortality: regression discontinuity evidence from the minimum drinking age.  American Economic Journal: Applied Economics, 1:164–182.
  * Cattaneo, (2010) Cattaneo, M. D. (2010).  Efficient semiparametric estimation of multi-valued treatment effects under ignorability.  Journal of Econometrics, 155:138–154.
  * Cattaneo et al., (2015) Cattaneo, M. D., Frandsen, B. R., and Titiunik, R. (2015).  Randomization inference in the regression discontinuity design: An application to party advantages in the US Senate.  Journal of Causal Inference, 3:1–24.
  * Chan et al., (2016) Chan, K. C. G., Yam, S. C. P., and Zhang, Z. (2016).  Globally efficient non-parametric inference of average treatment effects by empirical balancing calibration weighting.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 78:673–700.
  * Charig et al., (1986) Charig, C. R., Webb, D. R., Payne, S. R., and Wickham, J. E. (1986).  Comparison of treatment of renal calculi by open surgery, percutaneous nephrolithotomy, and extracorporeal shockwave lithotripsy.  British Medical Journal, 292:879–882.
  * Chen et al., (2007) Chen, H., Geng, Z., and Jia, J. (2007).  Criteria for surrogate end points.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 69:919–932.
  * Cheng and Small, (2006) Cheng, J. and Small, D. S. (2006).  Bounds on causal effects in three-arm trials with non-compliance.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 68:815–836.
  * Chernozhukov et al., (2018) Chernozhukov, V., Chetverikov, D., Demirer, M., Duflo, E., Hansen, C., Newey, W., and Robins, J. (2018).  Double/debiased machine learning for treatment and structural parameters.  Econometrics Journal, 21:C1–C68.
  * Chong et al., (2016) Chong, A., Cohen, I., Field, E., Nakasone, E., and Torero, M. (2016).  Iron deficiency and schooling attainment in peru.  American Economic Journal: Applied Economics, 8:222–55.
  * Cochran, (1938) Cochran, W. G. (1938).  The omission or addition of an independent variate in multiple linear regression.  Supplement to the Journal of the Royal Statistical Society, 5:171–176.
  * Cochran, (1953) Cochran, W. G. (1953).  Sampling Techniques.  New York: Wiley.
  * Cochran, (1957) Cochran, W. G. (1957).  Analysis of covariance: its nature and uses.  Biometrics, 13:261–281.
  * Cochran, (1965) Cochran, W. G. (1965).  The planning of observational studies of human populations (with discussion).  Journal of the Royal Statistical Society: Series A (General), 128:234–266.
  * Cochran, (1968) Cochran, W. G. (1968).  The effectiveness of adjustment by subclassification in removing bias in observational studies.  Biometrics, 24:295–313.
  * Cochran and Rubin, (1973) Cochran, W. G. and Rubin, D. B. (1973).  Controlling bias in observational studies: A review.  Sankhyā, 35:417–446.
  * Cornfield et al., (1959) Cornfield, J., Haenszel, W., Hammond, E. C., Lilienfeld, A. M., Shimkin, M. B., and Wynder, E. L. (1959).  Smoking and lung cancer: recent evidence and a discussion of some questions.  Journal of the National Cancer Institute, 22:173–203.
  * Cox, (1982) Cox, D. R. (1982).  Randomization and concomitant variables in the design of experiments.  In G. Kallianpur, P. R. K. and Ghosh, J. K., editors, Statistics and Probability: Essays in Honor of C. R. Rao, pages 197–202. North-Holland, Amsterdam.
  * Cox, (2007) Cox, D. R. (2007).  On a generalization of a result of W. G. Cochran.  Biometrika, 94:755–759.
  * Crump et al., (2009) Crump, R. K., Hotz, V. J., Imbens, G. W., and Mitnik, O. A. (2009).  Dealing with limited overlap in estimation of average treatment effects.  Biometrika, 96:187–199.
  * Cuzick et al., (1997) Cuzick, J., Edwards, R., and Segnan, N. (1997).  Adjusting for non-compliance and contamination in randomized clinical trials.  Statistics in Medicine, 16:1017–1029.
  * D’Amour et al., (2021) D’Amour, A., Ding, P., Feller, A., Lei, L., and Sekhon, J. (2021).  Overlap in observational studies with high-dimensional covariates.  Journal of Econometrics, 221:644–654.
  * Davey Smith and Ebrahim, (2003) Davey Smith, G. and Ebrahim, S. (2003).  “Mendelian randomization”: can genetic epidemiology contribute to understanding environmental determinants of disease?  International Journal of Epidemiology, 32:1–22.
  * Davison and Hinkley, (1997) Davison, A. C. and Hinkley, D. V. (1997).  Bootstrap Methods and Their Application.  Cambridge: Cambridge University Press.
  * Dawid, (1979) Dawid, A. P. (1979).  Conditional independence in statistical theory.  Journal of the Royal Statistical Society: Series B (Methodological), 41:1–15.
  * Dawid, (2000) Dawid, A. P. (2000).  Causal inference without counterfactuals (with discussion).  Journal of the American Statistical Association, 95:407–424.
  * Dehejia and Wahba, (1999) Dehejia, R. H. and Wahba, S. (1999).  Causal effects in nonexperimental studies: Reevaluating the evaluation of training programs.  Journal of the American statistical Association, 94:1053–1062.
  * Ding, (2016) Ding, P. (2016).  A paradox from randomization-based causal inference (with discussion).  Statistical Science, 32:331–345.
  * Ding, (2021) Ding, P. (2021).  The Frisch–Waugh–Lovell theorem for standard errors.  Statistics and Probability Letters, 168:108945.
  * Ding and Dasgupta, (2016) Ding, P. and Dasgupta, T. (2016).  A potential tale of two by two tables from completely randomized experiments.  Journal of American Statistical Association, 111:157–168.
  * Ding and Dasgupta, (2017) Ding, P. and Dasgupta, T. (2017).  A randomization-based perspective on analysis of variance: a test statistic robust to treatment effect heterogeneity.  Biometrika, 105:45–56.
  * Ding et al., (2019) Ding, P., Feller, A., and Miratrix, L. (2019).  Decomposing treatment effect variation.  Journal of the American Statistical Association, 114:304–317.
  * Ding et al., (2011) Ding, P., Geng, Z., Yan, W., and Zhou, X.-H. (2011).  Identifiability and estimation of causal effects by principal stratification with outcomes truncated by death.  Journal of the American Statistical Association, 106:1578–1591.
  * Ding and Li, (2018) Ding, P. and Li, F. (2018).  Causal inference: A missing data perspective.  Statistical Science, 33:214–237.
  * (78) Ding, P., Li, X., and Miratrix, L. W. (2017a).  Bridging finite and super population causal inference.  Journal of Causal Inference, 5:20160027.
  * Ding and Lu, (2017) Ding, P. and Lu, J. (2017).  Principal stratification analysis using principal scores.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 79:757–777.
  * Ding and Miratrix, (2015) Ding, P. and Miratrix, L. W. (2015).  To adjust or not to adjust? Sensitivity analysis of M-bias and butterfly-bias.  Journal of Causal Inference, 3:41–57.
  * Ding and VanderWeele, (2014) Ding, P. and VanderWeele, T. J. (2014).  Generalized Cornfield conditions for the risk difference.  Biometrika, 101:971–977.
  * Ding and VanderWeele, (2016) Ding, P. and VanderWeele, T. J. (2016).  Sensitivity analysis without assumptions.  Epidemiology, 27:368–377.
  * Ding and Vanderweele, (2016) Ding, P. and Vanderweele, T. J. (2016).  Sharp sensitivity bounds for mediation under unmeasured mediator-outcome confounding.  Biometrika, 103:483–490.
  * (84) Ding, P., VanderWeele, T. J., and Robins, J. M. (2017b).  Instrumental variables as bias amplifiers with general outcome and confounding.  Biometrika, 104:291–302.
  * Doll and Hill, (1950) Doll, R. and Hill, A. B. (1950).  Smoking and carcinoma of the lung.  British Medical Journal, 2:739.
  * Dorn, (1953) Dorn, H. F. (1953).  Philosophy of inferences from retrospective studies.  American Journal of Public Health and the Nations Health, 43:677–683.
  * Durrett, (2019) Durrett, R. (2019).  Probability: Theory and Examples.  Cambridge: Cambridge University Press.
  * Efron, (1979) Efron, B. (1979).  Bootstrap methods: Another look at the jackknife.  The Annals of Statistics, 7:1–26.
  * Efron and Feldman, (1991) Efron, B. and Feldman, D. (1991).  Compliance as an explanatory variable in clinical trials (with discussion).  Journal of the American Statistical Association, 86:9–17.
  * Eicker, (1967) Eicker, F. (1967).  Limit theorems for regressions with unequal and dependent errors.  In Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability, volume 1, pages 59–82. Berkeley, CA: University of California Press.
  * Fan and Gijbels, (1996) Fan, J. and Gijbels, I. (1996).  Local Polynomial Modelling and Its Applications.  New York: Chapman and Hall/CRC.
  * Fieller, (1954) Fieller, E. C. (1954).  Some problems in interval estimation.  Journal of the Royal Statistical Society: Series B (Methodological), 16:175–185.
  * Firth and Bennett, (1998) Firth, D. and Bennett, K. E. (1998).  Robust models in probability sampling (with discussion).  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 60:3–21.
  * Fisher, (1925) Fisher, R. A. (1925).  Statistical Methods for Research Workers.  Edinburgh by Oliver and Boyd, 1st edition.
  * Fisher, (1935) Fisher, R. A. (1935).  The Design of Experiments.  Edinburgh, London: Oliver and Boyd, 1st edition.
  * Fisher, (1957) Fisher, R. A. (1957).  Dangers of cigarette smoking [letter].  British Medical Journal, 2:297–298.
  * (97) Fogarty, C. B. (2018a).  On mitigating the analytical limitations of finely stratified experiments.  Journal of the Royal Statistical Society. Series B (Statistical Methodology), 80:1035–1056.
  * (98) Fogarty, C. B. (2018b).  Regression assisted inference for the average treatment effect in paired experiments.  Biometrika, 105:994–1000.
  * Follmann, (2000) Follmann, D. A. (2000).  On the effect of treatment among would-be treatment compliers: An analysis of the multiple risk factor intervention trial.  Journal of the American Statistical Association, 95:1101–1109.
  * Forastiere et al., (2018) Forastiere, L., Mattei, A., and Ding, P. (2018).  Principal ignorability in mediation analysis: through and beyond sequential ignorability.  Biometrika, 105:979–986.
  * Frangakis and Rubin, (2002) Frangakis, C. E. and Rubin, D. B. (2002).  Principal stratification in causal inference.  Biometrics, 58:21–29.
  * (102) Freedman, D. A. (2008a).  On regression adjustments in experiments with several treatments.  Annals of Applied Statistics, 2:176–196.
  * (103) Freedman, D. A. (2008b).  On regression adjustments to experimental data.  Advances in Applied Mathematics, 40:180–193.
  * (104) Freedman, D. A. (2008c).  Randomization does not justify logistic regression.  Statistical Science, 23:237–249.
  * Freedman and Berk, (2008) Freedman, D. A. and Berk, R. A. (2008).  Weighting regressions by propensity scores.  Evaluation Review, 32:392–409.
  * Funk et al., (2011) Funk, M. J., Westreich, D., Wiesen, C., Stürmer, T., Brookhart, M. A., and Davidian, M. (2011).  Doubly robust estimation of causal effects.  American Journal of Epidemiology, 173:761–767.
  * Gastwirth et al., (1998) Gastwirth, J. L., KRIEGER, A. M., and ROSENBAUM, P. R. (1998).  Cornfield’s inequality.  In Armitage, P. and Colton, T., editors, Encyclopedia of Biostatistics. New York: Wiley.
  * Gerber and Green, (2012) Gerber, A. S. and Green, D. P. (2012).  Field Experiments: Design, Analysis, and Interpretation.  WW Norton.
  * Gerber et al., (2008) Gerber, A. S., Green, D. P., and Larimer, C. W. (2008).  Social pressure and voter turnout: Evidence from a large-scale field experiment.  American Political Science Review, 102:33–48.
  * Gilbert and Hudgens, (2008) Gilbert, P. B. and Hudgens, M. G. (2008).  Evaluating candidate principal surrogate endpoints.  Biometrics, 64:1146–1154.
  * Gould, (1998) Gould, A. L. (1998).  Multi-centre trial analysis revisited.  Statistics in Medicine, 17:1779–1797.
  * Greevy et al., (2004) Greevy, R., Lu, B., Silber, J. H., and Rosenbaum, P. (2004).  Optimal multivariate matching before randomization.  Biostatistics, 5:263–275.
  * Guo and Basse, (2023) Guo, K. and Basse, G. (2023).  The generalized Oaxaca–Blinder estimator.  Journal of American Statistical Association, 118:524–536.
  * Hahn, (1998) Hahn, J. (1998).  On the role of the propensity score in efficient semiparametric estimation of average treatment effects.  Econometrica, 66:315–331.
  * Hahn et al., (2001) Hahn, J., Todd, P., and Van der Klaauw, W. (2001).  Identification and estimation of treatment effects with a regression-discontinuity design.  Econometrica, 69:201–209.
  * Hahn et al., (2020) Hahn, P. R., Murray, J. S., and Carvalho, C. M. (2020).  Bayesian regression tree models for causal inference: regularization, confounding, and heterogeneous effects.  Bayesian Analysis, 15:965–1056.
  * Hainmueller, (2012) Hainmueller, J. (2012).  Entropy balancing for causal effects: A multivariate reweighting method to produce balanced samples in observational studies.  Political Analysis, 20:25–46.
  * Hájek, (1960) Hájek, J. (1960).  Limiting distributions in simple random sampling from a finite population.  Publications of the Mathematics Institute of the Hungarian Academy of Science, 5:361–74.
  * Hájek, (1971) Hájek, J. (1971).  Comment on “an essay on the logical foundations of survey sampling, part one”.  The foundations of survey sampling, 236.
  * Hammond and Horn, (1958) Hammond, E. C. and Horn, D. (1958).  Smoking and death rates: report on forty four months of follow-up of $187,783$ men.  Journal of the American Medicial Association, 166:1159–1172, 1294–1308.
  * Hansen, (1982) Hansen, L. P. (1982).  Large sample properties of generalized method of moments estimators.  Econometrica, 50:1029–1054.
  * Hartley et al., (1969) Hartley, H. O., Rao, J. N. K., and Kiefer, G. (1969).  Variance estimation with one unit per stratum.  Journal of the American Statistical Association, 64:841–851.
  * Hausman, (1978) Hausman, J. A. (1978).  Specification tests in econometrics.  Econometrica, 46:1251–1271.
  * Hearst et al., (1986) Hearst, N., Newman, T. B., and Hulley, S. B. (1986).  Delayed effects of the military draft on mortality.  New England Journal of Medicine, 314:620–624.
  * Heckman and Navarro-Lozano, (2004) Heckman, J. and Navarro-Lozano, S. (2004).  Using matching, instrumental variables, and control functions to estimate economic choice models.  Review of Economics and Statistics, 86:30–57.
  * Heckman, (1979) Heckman, J. J. (1979).  Sample selection bias as a specification error.  Econometrica, 47:153–161.
  * Hennessy et al., (2016) Hennessy, J., Dasgupta, T., Miratrix, L., Pattanayak, C., and Sarkar, P. (2016).  A conditional randomization test to account for covariate imbalance in randomized experiments.  Journal of Causal Inference, 4:61–80.
  * Hernán et al., (2000) Hernán, M. Á., Brumback, B., and Robins, J. M. (2000).  Marginal structural models to estimate the causal effect of zidovudine on the survival of hiv-positive men.  Epidemiology, 11:561–570.
  * Hernán and Robins, (2020) Hernán, M. A. and Robins, J. M. (2020).  Causal Inference: What If.  Boca Raton: Chapman & Hall/CRC.
  * Hill et al., (2002) Hill, J., Waldfogel, J., and Brooks-Gunn, J. (2002).  Differential effects of high-quality child care.  Journal of Policy Analysis and Management, 21:601–627.
  * Hill, (2011) Hill, J. L. (2011).  Bayesian nonparametric modeling for causal inference.  Journal of Computational and Graphical Statistics, 20:217–240.
  * Hirano and Imbens, (2001) Hirano, K. and Imbens, G. W. (2001).  Estimation of causal effects using propensity score weighting: An application to data on right heart catheterization.  Health Services and Outcomes Research Methodology, 2:259–278.
  * Hirano et al., (2000) Hirano, K., Imbens, G. W., Rubin, D. B., and Zhou, X. H. (2000).  Assessing the effect of an influenza vaccine in an encouragement design.  Biostatistics, 1:69–88.
  * Ho et al., (2007) Ho, D. E., Imai, K., King, G., and Stuart, E. A. (2007).  Matching as nonparametric preprocessing for reducing model dependence in parametric causal inference.  Political Analysis, 15:199–236.
  * Ho et al., (2011) Ho, D. E., Imai, K., King, G., and Stuart, E. A. (2011).  Matchit: nonparametric preprocessing for parametric causal inference.  Journal of Statistical Software, 42:1–28.
  * Hodges and Lehmann, (1962) Hodges, J. L. and Lehmann, E. L. (1962).  Rank methods for combination of independent experiments in analysis of variance.  Annals of Mathematical Statistics, 33:482–497.
  * Holland, (1986) Holland, P. W. (1986).  Statistics and causal inference (with discussion).  Journal of the American statistical Association, 81:945–960.
  * Hong and Raudenbush, (2008) Hong, G. and Raudenbush, S. W. (2008).  Causal inference for time-varying instructional treatments.  Journal of Educational and Behavioral Statistics, 33:333–362.
  * Horvitz and Thompson, (1952) Horvitz, D. G. and Thompson, D. J. (1952).  A generalization of sampling without replacement from a finite universe.  Journal of the American statistical Association, 47:663–685.
  * Huber, (1967) Huber, P. J. (1967).  The behavior of maximum likelihood estimates under nonstandard conditions.  In Cam, L. M. L. and Neyman, J., editors, Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability, volume 1, pages 221–233. Berkeley, California: University of California Press.
  * Hyman, (1955) Hyman, H. H. (1955).  Survey Design and Analysis: Principles, Cases, and Procedures.  Glencoe, IL: Free Press.
  * (142) Imai, K. (2008a).  Sharp bounds on the causal effects in randomized experiments with “truncation-by-death”.  Statistics and Probability Letters, 78:144–149.
  * (143) Imai, K. (2008b).  Variance identification and efficiency analysis in randomized experiments under the matched-pair design.  Statistics in Medicine, 27:4857–4873.
  * Imai et al., (2010) Imai, K., Keele, L., and Yamamoto, T. (2010).  Identification, inference and sensitivity analysis for causal mediation effects.  Statistical Science, 25:51–71.
  * Imai and Van Dyk, (2004) Imai, K. and Van Dyk, D. A. (2004).  Causal inference with general treatment regimes: Generalizing the propensity score.  Journal of the American Statistical Association, 99:854–866.
  * Imbens, (2020) Imbens, G. (2020).  Potential outcome and directed acyclic graph approaches to causality: Relevance for empirical practice in economics.  Journal of Economic Literature, 58:1129–1179.
  * Imbens, (2003) Imbens, G. W. (2003).  Sensitivity to exogeneity assumptions in program evaluation.  American Economic Review, 93:126–132.
  * Imbens, (2004) Imbens, G. W. (2004).  Nonparametric estimation of average treatment effects under exogeneity: A review.  Review of Economics and Statistics, 86:4–29.
  * Imbens, (2014) Imbens, G. W. (2014).  Instrumental variables: An econometrician’s perspective.  Statistical Science, 29:323–358.
  * Imbens, (2015) Imbens, G. W. (2015).  Matching methods in practice: Three examples.  Journal of Human Resources, 50:373–419.
  * Imbens and Angrist, (1994) Imbens, G. W. and Angrist, J. D. (1994).  Identification and estimation of local average treatment effects.  Econometrica, 62:467–475.
  * Imbens and Lemieux, (2008) Imbens, G. W. and Lemieux, T. (2008).  Regression discontinuity designs: A guide to practice.  Journal of Econometrics, 142:615–635.
  * Imbens and Manski, (2004) Imbens, G. W. and Manski, C. F. (2004).  Confidence intervals for partially identified parameters.  Econometrica, 72:1845–1857.
  * Imbens and Rubin, (1997) Imbens, G. W. and Rubin, D. B. (1997).  Estimating outcome distributions for compliers in instrumental variables models.  Review of Economic Studies, 64:555–574.
  * Imbens and Rubin, (2015) Imbens, G. W. and Rubin, D. B. (2015).  Causal Inference for Statistics, Social, and Biomedical Sciences: An Introduction.  Cambridge: Cambridge University Press.
  * Investigators et al., (2014) Investigators, I. T. et al. (2014).  Endovascular or open repair strategy for ruptured abdominal aortic aneurysm: 30 day outcomes from improve randomised trial.  British Medical Journal, 348:f7661.
  * Ioannidis et al., (2019) Ioannidis, J. P. A., Tan, Y. J., and Blum, M. R. (2019).  Limitations and misinterpretations of E-values for sensitivity analyses of observational studies.  Annals of Internal Medicine, 170:108–111.
  * Jackson et al., (2006) Jackson, L. A., Jackson, M. L., Nelson, J. C., Neuzil, K. M., and Weiss, N. S. (2006).  Evidence of bias in estimates of influenza vaccine effectiveness in seniors.  International Journal of Epidemiology, 35:337–344.
  * Jiang and Ding, (2020) Jiang, Z. and Ding, P. (2020).  Measurement errors in the binary instrumental variable model.  Biometrika, 107:238–245.
  * Jiang and Ding, (2021) Jiang, Z. and Ding, P. (2021).  Identification of causal effects within principal strata using auxiliary variables.  Statistical Science, 36:493–508.
  * Jiang et al., (2016) Jiang, Z., Ding, P., and Geng, Z. (2016).  Principal causal effect identification and surrogate end point evaluation by multiple trials.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 78:829–848.
  * Jiang et al., (2022) Jiang, Z., Yang, S., and Ding, P. (2022).  Multiply robust estimation of causal effects under principal ignorability.  Journal of the Royal Statistical Society - Series B (Statistical Methodology), 84:1423–1445.
  * Jo and Stuart, (2009) Jo, B. and Stuart, E. A. (2009).  On the use of propensity scores in principal causal effect estimation.  Statistics in Medicine, 28:2857–2875.
  * Jo et al., (2011) Jo, B., Stuart, E. A., MacKinnon, D. P., and Vinokur, A. D. (2011).  The use of propensity scores in mediation analysis.  Multivariate Behavioral Research, 46:425–452.
  * Judd and Kenny, (1981) Judd, C. M. and Kenny, D. A. (1981).  Process analysis estimating mediation in treatment evaluations.  Evaluation Review, 5:602–619.
  * Kang and Schafer, (2007) Kang, J. D. Y. and Schafer, J. L. (2007).  Demystifying double robustness: A comparison of alternative strategies for estimating a population mean from incomplete data.  Statistical Science, 22:523–539.
  * Katan, (1986) Katan, M. B. (1986).  Apoupoprotein E isoforms, serum cholesterol, and cancer.  Lancet, 327:507–508.
  * King and Zeng, (2006) King, G. and Zeng, L. (2006).  The dangers of extreme counterfactuals.  Political Analysis, 14:131–159.
  * Kitagawa, (2015) Kitagawa, T. (2015).  A test for instrument validity.  Econometrica, 83:2043–2063.
  * Koenker and Xiao, (2002) Koenker, R. and Xiao, Z. (2002).  Inference on the quantile regression process.  Econometrica, 70:1583–1612.
  * Künzel et al., (2019) Künzel, S. R., Sekhon, J. S., Bickel, P. J., and Yu, B. (2019).  Metalearners for estimating heterogeneous treatment effects using machine learning.  Proceedings of the National Academy of Sciences of the United States of America, 116:4156–4165.
  * Kurth et al., (2005) Kurth, T., Walker, A. M., Glynn, R. J., Chan, K. A., Gaziano, J. M., Berger, K., and Robins, J. M. (2005).  Results of multivariable logistic regression, propensity matching, propensity adjustment, and propensity-based weighting under conditions of nonuniform effect.  American Journal of Epidemiology, 163:262–270.
  * LaLonde, (1986) LaLonde, R. J. (1986).  Evaluating the econometric evaluations of training programs with experimental data.  American Economic Review, 76:604–620.
  * Lee, (2008) Lee, D. S. (2008).  Randomized experiments from non-random selection in US House elections.  Journal of Econometrics, 142:675–697.
  * Lee, (2009) Lee, D. S. (2009).  Training, wages, and sample selection: Estimating sharp bounds on treatment effects.  Review of Economic Studies, 76:1071–1102.
  * Lee and Lemieux, (2010) Lee, D. S. and Lemieux, T. (2010).  Regression discontinuity designs in economics.  Journal of Economic Literature, 48:281–355.
  * Lee, (2018) Lee, M.-J. (2018).  Simple least squares estimator for treatment effects using propensity score residuals.  Biometrika, 105:149–164.
  * Lee, (2011) Lee, W.-C. (2011).  Bounding the bias of unmeasured factors with confounding and effect-modifying potentials.  Statistics in Medicine, 30:1007–1017.
  * Lehmann, (1975) Lehmann, E. L. (1975).  Nonparametrics: Statistical Methods Based on Ranks.  California: Holden-Day, Inc.
  * Lei and Ding, (2021) Lei, L. and Ding, P. (2021).  Regression adjustment in completely randomized experiments with a diverging number of covariates.  Biometrika, 108:815–828.
  * Li et al., (2015) Li, F., Mattei, A., and Mealli, F. (2015).  Evaluating the causal effect of university grants on student dropout: evidence from a regression discontinuity design using principal stratification.  Annals of Applied Statistics, 9:1906–1931.
  * (182) Li, F., Morgan, K. L., and Zaslavsky, A. M. (2018a).  Balancing covariates via propensity score weighting.  Journal of the American Statistical Association, 113:390–400.
  * Li et al., (2019) Li, F., Thomas, L. E., and Li, F. (2019).  Addressing extreme propensity scores via the overlap weights.  American Journal of Epidemiology, 188:250–257.
  * Li and Ding, (2016) Li, X. and Ding, P. (2016).  Exact confidence intervals for the average causal effect on a binary outcome.  Statistics in Medicine, 35:957–960.
  * Li and Ding, (2017) Li, X. and Ding, P. (2017).  General forms of finite population central limit theorems with applications to causal inference.  Journal of the American Statistical Association, 112:1759–1769.
  * Li and Ding, (2020) Li, X. and Ding, P. (2020).  Rerandomization and regression adjustment.  Journal of the Royal Statistical Society, Series B (Statistical Methodology), 82:241–268.
  * (187) Li, X., Ding, P., and Rubin, D. B. (2018b).  Asymptotic theory of rerandomization in treatment-control experiments.  Proceedings of the National Academy of Sciences of the United States of America, 115:9157–9162.
  * Lin, (2013) Lin, W. (2013).  Agnostic notes on regression adjustments to experimental data: Reexamining Freedman’s critique.  Annals of Applied Statistics, 7:295–318.
  * Lin et al., (2023) Lin, Z., Ding, P., and Han, F. (2023).  Estimation based on nearest neighbor matching: from density ratio to average treatment effect.  Econometrica.
  * Lind, (1753) Lind, J. (1753).  A treatise of the scurvy.  Three Parts. Containing an Inquiry into the Nature, Causes and Cure, of that Disease. Together with a Critical and Chronological View of what has been Published on the Subject.
  * Lipsitch et al., (2010) Lipsitch, M., Tchetgen Tchetgen, E., and Cohen, T. (2010).  Negative controls: a tool for detecting confounding and bias in observational studies.  Epidemiology, 21:383–388.
  * Little and An, (2004) Little, R. and An, H. (2004).  Robust likelihood-based analysis of multivariate data with missing values.  Statistica Sinica, 14:949–968.
  * Liu and Yang, (2020) Liu, H. and Yang, Y. (2020).  Regression-adjusted average treatment effect estimates in stratified randomized experiments.  Biometrika, 107:935–948.
  * Long and Ervin, (2000) Long, J. S. and Ervin, L. H. (2000).  Using heteroscedasticity consistent standard errors in the linear regression model.  American Statistician, 54:217–224.
  * Lu and Ding, (2023) Lu, S. and Ding, P. (2023).  Flexible sensitivity analysis for causal inference in observational studies subject to unmeasured confounding.  https://arxiv.org/abs/2305.17643.
  * Lumley et al., (2011) Lumley, T., Shaw, P. A., and Dai, J. Y. (2011).  Connections between survey calibration estimators and semiparametric models for incomplete data.  International Statistical Review, 79:200–220.
  * Lunceford and Davidian, (2004) Lunceford, J. K. and Davidian, M. (2004).  Stratification and weighting via the propensity score in estimation of causal treatment effects: a comparative study.  Statistics in Medicine, 23:2937–2960.
  * Luo et al., (2021) Luo, X., Dasgupta, T., Xie, M., and Liu, R. Y. (2021).  Leveraging the fisher randomization test using confidence distributions: Inference, combination and fusion learning.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 83:777–797.
  * Manski, (1990) Manski, C. F. (1990).  Nonparametric bounds on treatment effects.  American Economic Review, 2:319–323.
  * Manski, (2003) Manski, C. F. (2003).  Partial Identification of Probability Distributions.  New York: Springer.
  * Mattei et al., (2013) Mattei, A., Li, F., and Mealli, F. (2013).  Exploiting multiple outcomes in bayesian principal stratification analysis with application to the evaluation of a job training program.  Annals of Applied Statistics, 7:2336–2360.
  * McCrary, (2008) McCrary, J. (2008).  Manipulation of the running variable in the regression discontinuity design: A density test.  Journal of Econometrics, 142:698–714.
  * McDonald et al., (1992) McDonald, C. J., Hui, S. L., and Tierney, W. M. (1992).  Effects of computer reminders for influenza vaccination on morbidity during influenza epidemics.  MD Computing: Computers in Medical Practice, 9:304–312.
  * McGrath et al., (2021) McGrath, S., Young, J. G., and Hernán, M. A. (2021).  Revisiting the g-null paradox.  Epidemiology, 33:114–120.
  * Mealli and Pacini, (2013) Mealli, F. and Pacini, B. (2013).  Using secondary outcomes to sharpen inference in randomized experiments with noncompliance.  Journal of the American Statistical Association, 108:1120–1131.
  * Meinert et al., (1970) Meinert, C. L., Knatterud, G. L., Prout, T. E., and Klimt, C. R. (1970).  A study of the effects of hypoglycemic agents on vascular complications in patients with adult-onset diabetes. ii. mortality results.  Diabetes, 19:Suppl–789.
  * Mercatanti and Li, (2014) Mercatanti, A. and Li, F. (2014).  Do debit cards increase household spending? evidence from a semiparametric causal analysis of a survey.  Annals of Applied Statistics, 8:2485–2508.
  * Ming and Rosenbaum, (2000) Ming, K. and Rosenbaum, P. R. (2000).  Substantial gains in bias reduction from matching with a variable number of controls.  Biometrics, 56:118–124.
  * Ming and Rosenbaum, (2001) Ming, K. and Rosenbaum, P. R. (2001).  A note on optimal matching with variable controls using the assignment algorithm.  Journal of Computational and Graphical Statistics, 10:455–463.
  * Miratrix et al., (2013) Miratrix, L. W., Sekhon, J. S., and Yu, B. (2013).  Adjusting treatment effect estimates by post-stratification in randomized experiments.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 75:369–396.
  * Morgan and Rubin, (2012) Morgan, K. L. and Rubin, D. B. (2012).  Rerandomization to improve covariate balance in experiments.  Annals of Statistics, 40:1263–1282.
  * Mukerjee et al., (2018) Mukerjee, R., Dasgupta, T., and Rubin, D. B. (2018).  Using standard tools from finite population sampling to improve causal inference for complex experiments.  Journal of the American Statistical Association, 113:868–881.
  * Naimi et al., (2017) Naimi, A. I., Cole, S. R., and Kennedy, E. H. (2017).  An introduction to g methods.  International Journal of Epidemiology, 46:756–762.
  * Negi and Wooldridge, (2021) Negi, A. and Wooldridge, J. M. (2021).  Revisiting regression adjustment in experiments with heterogeneous treatment effects.  Econometric Reviews, 40:504–534.
  * Neyman, (1923) Neyman, J. (1923).  On the application of probability theory to agricultural experiments. essay on principles (with discussion). section 9 (translated). reprinted ed.  Statistical Science, 5:465–472.
  * Neyman, (1934) Neyman, J. (1934).  On the two different aspects of the representative method: the method of stratified sampling and the method of purposive selection (with discussion).  Journal of the Royal Statistical Society, 97:558–625.
  * Neyman, (1935) Neyman, J. (1935).  Statistical problems in agricultural experimentation (with discussion).  Supplement to the Journal of the Royal Statistical Society, 2:107–180.
  * Nguyen et al., (2021) Nguyen, T. Q., Schmid, I., Ogburn, E. L., and Stuart, E. A. (2021).  Clarifying causal mediation analysis for the applied researcher: Effect identification via three assumptions and five potential outcomes.  Psychological Methods, 26:255–271.
  * Otsu and Rai, (2017) Otsu, T. and Rai, Y. (2017).  Bootstrap inference of matching estimators for average treatment effects.  Journal of the American Statistical Association, 112:1720–1732.
  * Pearl, (1995) Pearl, J. (1995).  Causal diagrams for empirical research (with discussion).  Biometrika, 82:669–688.
  * Pearl, (2000) Pearl, J. (2000).  Causality: Models, Reasoning and Inference.  Cambridge: Cambridge University Press.
  * Pearl, (2001) Pearl, J. (2001).  Direct and indirect effects.  In Breese, J. S. and Koller, D., editors, Proceedings of the 17th Conference on Uncertainty in Artificial Intelligence, pages 411–420. pp. 411–420. San Francisco: Morgan Kaufmann Publishers Inc.
  * Pearl, (2010) Pearl, J. (2010).  On a class of bias-amplifying variables that endanger effect estimates.  In Grunwald, P. and Spirtes, P., editors, Proceedings of the Twenty-Sixth Conference on Uncertainty in Artificial Intelligence (UAI 2010), Corvallis, OR: 425–432. Association for Uncetainty in Artificial Intelligence.
  * Pearl, (2011) Pearl, J. (2011).  Invited commentary: Understanding bias amplification.  American Journal of Epidemiology, 174:1223–1227.
  * Pearl, (2018) Pearl, J. (2018).  Does obesity shorten life? Or is it the soda? On non-manipulable causes.  Journal of Causal Inference, 6:20182001.
  * Pearl and Bareinboim, (2014) Pearl, J. and Bareinboim, E. (2014).  External validity: From do-calculus to transportability across populations.  Statistical Science, 29:579–595.
  * Permutt and Hebel, (1989) Permutt, T. and Hebel, J. R. (1989).  Simultaneous-equation estimation in a clinical trial of the effect of smoking on birth weight.  Biometrics, 45:619–622.
  * Phipson and Smyth, (2010) Phipson, B. and Smyth, G. K. (2010).  Permutation p-values should never be zero: calculating exact p-values when permutations are randomly drawn.  Statistical Applications in Genetics and Molecular Biology, 9:Article39.
  * Pimentel et al., (2015) Pimentel, S. D., Yoon, F., and Keele, L. (2015).  Variable-ratio matching with fine balance in a study of the Peer Health Exchange.  Statistics in Medicine, 34:4070–4082.
  * Poole, (2010) Poole, C. (2010).  On the origin of risk relativism.  Epidemiology, 21:3–9.
  * Popper, (1963) Popper, K. (1963).  Conjectures and Refutations: The Growth of Scientific Knowledge.  Routledge.
  * Powers and Swinton, (1984) Powers, D. E. and Swinton, S. S. (1984).  Effects of self-study for coachable test item types.  Journal of Educational Psychology, 76:266–278.
  * Prentice and Pyke, (1979) Prentice, R. L. and Pyke, R. (1979).  Logistic disease incidence models and case-control studies.  Biometrika, 66:403–411.