Asymptotic Theory of Rerandomization in

arXiv:1604.00698v4 [math.ST] 12 Aug 2017

Treatment-Control Experiments
Xinran Li, Peng Ding and Donald B. Rubin

∗

Abstract
Although complete randomization ensures covariate balance on average, the chance for observing significant differences between treatment and control covariate distributions increases
with many covariates. Rerandomization discards randomizations that do not satisfy a predetermined covariate balance criterion, generally resulting in better covariate balance and more
precise estimates of causal effects. Previous theory has derived finite sample theory for rerandomization under the assumptions of equal treatment group sizes, Gaussian covariate and outcome distributions, or additive causal effects, but not for the general sampling distribution of
the difference-in-means estimator for the average causal effect. To supplement existing results,
we develop asymptotic theory for rerandomization without these assumptions, which reveals a
non-Gaussian asymptotic distribution for this estimator, specifically a linear combination of a
Gaussian random variable and a truncated Gaussian random variable. This distribution follows
because rerandomization affects only the projection of potential outcomes onto the covariate
space but does not affect the corresponding orthogonal residuals. We also demonstrate that,
compared to complete randomization, rerandomization reduces the asymptotic sampling variances and quantile ranges of the difference-in-means estimator. Moreover, our work allows the
construction of accurate large-sample confidence intervals for the average causal effect, thereby
revealing further advantages of rerandomization over complete randomization.

Keywords: Causal inference; Covariate balance; Geometry of rerandomization; Mahalanobis distance; Quantile range; Tiers of covariates.
∗

Xinran Li is Doctoral Candidate, Department of Statistics, Harvard University, Cambridge, MA 02138 (E-mail:
xinranli@fas.harvard.edu). Peng Ding is Assistant Professor, Department of Statistics, University of California,
Berkeley, CA 94720 (E-mail: pengdingpku@berkeley.edu). Donald B. Rubin is Professor, Department of Statistics,
Harvard University, Cambridge, MA 02138 (E-mail: rubin@stat.harvard.edu).

1

1.

Introduction

Ever since Fisher (1925, 1926, 1935)’s seminal work, randomized experiments have become the
“gold standard” for drawing causal inferences. Complete randomization balances the covariate
distributions between treatment groups in expectation, thereby ensuring the existence of unbiased
estimators of the average causal effect. Covariate imbalance, however, often occurs in specific randomized experiments, as recognized by Fisher (1926) and later researchers (e.g., Student 1938;
Greevy et al. 2004; Hansen and Bowers 2008; Keele et al. 2009; Bruhn and McKenzie 2009; Krieger
et al. 2016; Athey and Imbens 2017). The standard approach advocated by Fisher (1935), stratification or blocking, ensures balance with a few discrete covariates; see Cochran and Cox (1992)
and Imbens and Rubin (2015) for detailed discussions.
When a randomized allocation is unbalanced, it is reasonable to discard that allocation and
re-draw another one until a certain pre-determined covariate balance criterion is satisfied. This is
rerandomization, an experimental design hinted by R. A. Fisher (cf. Savage 1962, page 88) and
Cox (1982, 2009), and formally proposed by Rubin (2008) and Morgan and Rubin (2012). Note
that rerandomization is conceptually the same as restricted or constrained randomization (e.g.,
Yates 1948; Grundy and Healy 1950; Youden 1972; Bailey 1983). For more historical discussion,
see Fienberg and Hinkley (1980, page 45), Speed (1992), Lehmann (2011, page 57), and Morgan
and Rubin (2012).
Morgan and Rubin (2012) showed that the difference-in-means estimator is generally unbiased
for the average causal effect under rerandomization with equal-sized treatment groups, and obtained
the sampling variance of this estimator under additional assumptions of Gaussian covariate and
outcome distributions and additive causal effects. When rerandomization is applied when these
assumptions do not hold, statistical inference becomes more challenging, because the Gaussian
distributional theory that is justified by the central limit theorem under complete randomization
(cf. Hájek 1960; Lin 2013) no longer generally holds. Some applied researchers believe that “the
only analysis that we can be completely confident in is a permutation test or rerandomization test”
(Bruhn and McKenzie 2009). However, randomization-based tests require sharp null hypotheses
that all individual causal effects are known from observed values.
Analogous to the repeated sampling properties for complete randomization (Neyman 1923;
2

Imbens and Rubin 2015), we evaluate the sampling properties of the difference-in-means estimator
when rerandomization is used, where all potential outcomes and covariates are regarded as fixed
quantities and all randomness arises solely from the random treatment assignments. The geometry
of rerandomization reveals non-Gaussian asymptotic distributions, which serve as the foundation for
constructing large-sample confidence intervals for average causal effects. Furthermore, we compare
the lengths of quantile ranges of the asymptotic distributions of the difference-in-means estimator
under rerandomization and complete randomization, extending Morgan and Rubin (2012, 2015)’s
comparison of their sampling variances.

2.

Framework, Notation, and Basic Results

2.1.

Covariate imbalance and rerandomization

Inferring the causal effect of some binary treatment on an outcome Y is of central interest in many
studies. We consider an experiment with n units, with n1 assigned to treatment and n0 assigned
to control, n = n1 + n0 , indexed by i = 1, . . . , n. Before conducting the experiment, we collect
K covariates Xi = (X1i , X2i , . . . , XKi ) for each unit, which can possibly include transformations
of basic covariates and their interactions. Let Zi be the indicator variable for unit i assigned to
treatment (Zi = 1 if active treatment level; Zi = 0 if the control level), and Z = (Z1 , Z2 , . . . , Zn )0
be the treatment assignment column vector. In a completely randomized experiment (CRE), the
distribution of Z is such that each value, z = (z1 , . . . , zn )0 , of Z has probability n1 !n0 !/n!, where
Pn
Pn
i=1 zi = n1 and
i=1 (1 − zi ) = n0 , which does not depend on the values of any observed or
unobserved covariates. The difference-in-means vector of the covariates between treatment and
control groups is
τ̂X

n
n
1 X
1 X
=
Zi Xi −
(1 − Zi )Xi .
n1
n0
i=1

i=1

Although on average τ̂X has mean zero over all


n
n1

randomizations, for any realized value of

Z, imbalancedness in covariate distributions between treatment groups often occurs. As pointed
out by Morgan and Rubin (2012), with 10 independent covariates and significance level 5%, the
probability of a significant difference for at least one covariate is 40%.
When significant covariate imbalance arises in a drawn allocation, it is reasonable to discard
3

the unlucky allocation and draw another treatment assignment vector until some a priori covariate balance criterion is satisfied. This is rerandomization, an intuitive experimental design tool
apparently personally advocated by R. A. Fisher (see the discussion by Rubin 2008) and formally
discussed by Morgan and Rubin (2012).
In general, rerandomization entails the following steps:
(1) collect covariate data;
(2) specify a balance criterion to determine whether a randomization is acceptable or not;
(3) randomize the units to treatment and control groups;
(4) if the balance criterion is satisfied, proceed to Step (5); otherwise, return to Step (3);
(5) conduct the experiment using the final randomization obtained in Step (4);
(6) analyze the data taking into account the rerandomization used in Steps (2)–(4).
Although the balance criterion in Step (2) can be general, Morgan and Rubin (2012) suggested
using the Mahalanobis distance between covariate means in treatment and control groups, and
Morgan and Rubin (2015) suggested considering tiers of covariates according to their presumed
importance in predicting the outcomes in this experiment. We will discuss these two types of
rerandomization in detail, and apposite statistical inference after these rerandomizations as implied
by Step (6). We then extend the theory to rerandomizations under more general covariate balance
criteria in Section 5.

2.2.

Potential outcomes and definitions of finite population quantities

We use the potential outcomes framework (sometimes called the Rubin Causal Model; see Holland
1986; Imbens and Rubin 2015) to define causal effects, and let Yi (1) and Yi (0) denote the potential
outcomes of unit i under active treatment and control, respectively. On the difference scale, the
individual causal effect for unit i is τi = Yi (1) − Yi (0), and the average causal effect in the finite
P
P
population of n units is τY = ni=1 τi /n. Let Ȳ (z) = ni=1 Yi (z)/n be the finite population average
of potential outcomes under treatment arm z, X̄ the finite population average of covariates, SY2 (z)
the finite population variance (with divisor n − 1) of the potential outcomes under treatment arm
4

0
z, SY (z),X = SX,Y
(z) the finite population covariance between potential outcomes and covariates,
2 the finite population covariance matrix of covariates. For simplicity, we avoid notation for
and SX

these quantities’ dependence on n. Notice that these quantities are fixed, and are not dependent
on the randomization or rerandomization scheme.

2.3.

Repeated sampling inference in a CRE

The observed outcome for unit i is Yi = Zi Yi (1) + (1 − Zi )Yi (0), a function of treatment assignment
and potential outcomes. In a CRE, Neyman (1923) showed that, for estimating τY , the differencein-means estimator
n
n
1 X
1 X
τ̂Y =
Zi Yi −
(1 − Zi )Yi
n1
n0
i=1

i=1

is unbiased (the expectation of τ̂Y over all randomizations is τY ), and obtained its sampling variance over all randomizations for constructing a large-sample confidence interval for τY . However,
Neyman (1923)’s interval is not accurate if rerandomization is used, except in an asymptotic conservative sense.
Let r1 = n1 /n and r0 = n0 /n be the proportions of units receiving treatment and control.
According to the finite population central limit theorem (Hájek 1960), under some regularity con√
0 ) is Gaussian
ditions, the large n sampling distribution, over all randomizations, of n(τ̂Y − τY , τ̂X
with mean zero and covariance matrix V , where

 Vτ τ
V =
Vxτ





−1 2
r1 SY (1)

r0−1 SY2 (0)

Sτ2

−
Vτ x 
+
=
r1−1 SX,Y (1) + r0−1 SX,Y (0)
Vxx

r1−1 SY (1),X



+

r0−1 SY (0),X 

2
(r1 r0 )−1 SX

.

Note again that we are conducting randomization-based inference, where all the covariates and
potential outcomes are fixed numbers, and randomness comes solely from the treatment assignment.
We embed n units into an infinite sequence of finite populations with increasing sizes, and a sufficient
√
0 ) is as follows (Li and Ding 2016).
condition for the asymptotic Gaussianity of n(τ̂Y − τY , τ̂X
Condition 1. As n → ∞, for z = 0, 1,
(i) rz , the proportion of units under treatment arm z, has positive limits,
2 and S
(ii) the finite population variances and covariances SY2 (z) , Sτ2 , SX
X,Y (z) have limiting values,

5

(iii) max1≤i≤n |Yi (z) − Ȳ (z)|2 /n → 0 and max1≤i≤n kXi − X̄k22 /n → 0.
.

We introduce the notation ∼ for two sequences of random vectors converging weakly to the
√
.
0 ) ∼
(A, B 0 ), where
same distribution. Therefore, under CRE and Condition 1, n(τ̂Y − τY , τ̂X
(A, B 0 ) is a random vector from N (0, V ).

3.

Rerandomization using the Mahalanobis distance

3.1.

Mahalanobis distance

The Mahalanobis distance between the covariate means in treatment and control groups is
0
τ̂X
{Var(τ̂X )}−1 τ̂X =

√

0 −1 √

nτ̂X Vxx
nτ̂X ,

2 is a fixed and known K × K matrix in our finite population
recalling that Vxx = (r1 r0 )−1 SX

setting. A rerandomization scheme proposed by Morgan and Rubin (2012) accepts only those
randomizations with the Mahalanobis distance less than or equal to a, a pre-specified threshold.
Let
−1
M = {µ : µ0 Vxx
µ ≤ a}

√
denote the acceptance region for nτ̂X ; that is, a treatment assignment vector Z is accepted if
√
and only if the corresponding nτ̂X ∈ M. Below we use ReM to denote rerandomization using
this criterion.
Several practical issues with ReM are worth mentioning. First, if we include transformations
and interactions of X, then ReM can incorporate a wide class of rerandomization schemes. Second, for small sample sizes, it can be that there does not exist any randomization satisfying the
balance criterion. However, according to the finite population central limit theorem, the acceptance probability of a randomization is asymptotically pa = P (χ2K ≤ a). Therefore, for relatively
large sample size, there usually exist many randomizations satisfying the balance criterion with
a > 0. In practice, we would like to choose the asymptotic acceptance probability to be small, e.g.,
pa = 0.001. However, we do not want pa to be too small, such as accepting only those assignments
with the smallest Mahalanobis distance. Too small pa will result in few randomizations, making
6

the repeated sampling inference intractable, even asymptotically, as well as the randomization tests
powerless (Morgan and Rubin 2012). Furthermore, as illustrated by later examples, the gain from
reducing pa usually decreases as pa becomes smaller.

3.2.

Multiple correlation between potential outcomes and covariates

We define the finite population squared multiple correlation between the potential outcome Y (z)
and the covariates X as R2 (z) for z = 1, 0, and the finite population squared multiple correlation between the individual causal effect and the covariates as R2 (τ ). Note that R2 (1), R2 (0)
and R2 (τ ) are quantities of the finite population, which do not depend on the randomization or
rerandomization scheme. Similar measures also appeared in Cochran (1965) and Rubin (1976).
We further define an R2 -type measure that is a function of the finite population quantities as
well as the proportions of the group sizes:

2

R =

SY2 (1)
r1 Vτ τ

2

R (1) +

SY2 (0)
r0 Vτ τ

R2 (0) −

Sτ2 2
R (τ ).
Vτ τ

When the causal effect is additive, Sτ2 = 0 and SY2 (1) = SY2 (0) , and then R2 = R2 (1) = R2 (0) reduces
to the squared multiple correlation between X and Y (1) or Y (0).
The following proposition states that under CRE R2 is the proportion of the sampling variance
of τ̂Y explained by τ̂X in linear projection.
Proposition 1. The sampling squared multiple correlation between τ̂Y and τ̂X under CRE is R2 ,
which can be equivalently written as

R2 = Corr(τ̂Y , τ̂X ) =

r1−1 SY2 (1)|X + r0−1 SY2 (0)|X − Sτ2|X
r1−1 SY2 (1) + r0−1 SY2 (0) − Sτ2

,

where SY2 (z)|X and Sτ2|X are the finite population variances of the linear projections of the potential
outcomes and individual causal effects on covariates.

3.3.

Asymptotic sampling distribution of τ̂Y under ReM

With rerandomization, we accept the randomizations satisfying the covariate balance criterion,
√
and therefore the sampling distribution of n(τ̂Y − τY ) over rerandomizations is the same as its
7

sampling distribution over a CRE conditional on

√

nτ̂X satisfying the covariate balance criterion.

Although the following proposition holds for rerandomization with more general balance criteria,
we first state it for ReM.
Proposition 2. Under ReM and Condition 1,




√

 n(τ̂Y − τY )

 √
nτ̂X

√

.

nτ̂X ∈ M ∼

 
A
 
B

B ∈ M,

(1)

recalling from earlier that (A, B 0 ) is a random vector following N (0, V ).
Simply stated,

√

n(τ̂Y − τY ) has two parts: the part unrelated to the covariates, which we call

ε0 , and thus unaffected by rerandomization, and the other part related to the covariates, which we
call LK,a , and thus affected by rerandomization. Therefore, the asymptotic distribution of τ̂Y is
a linear combination of two independent random variables: ε0 ∼ N (0, 1) is a standard Gaussian
random variable, and LK,a is a random variable following the distribution of D1 | D 0 D ≤ a, where
D = (D1 , . . . , DK )0 ∼ N (0, IK ).
Theorem 1. Under ReM and Condition 1,
√

n(τ̂Y − τY ) |

√

.

nτ̂X ∈ M ∼


p
√
p
Vτ τ
1 − R2 · ε0 + R2 · LK,a ,

(2)

where ε0 is independent of LK,a .
The coefficients of the linear combination are functions of R2 , which measures the association
between the potential outcomes and the covariates. When R2 = 0, the right hand side of (2) be√
comes a Gaussian random variable, the same as the asymptotic distribution of n(τ̂Y − τY ) under
√
CRE in Section 2.3; when R2 = 1, (2) reduces to Vτ τ · LK,a , a random variable with bounded
√
√
support [− aVτ τ , aVτ τ ]. Importantly, the definition of R2 is based on linear projections but not
linear models of the potential outcomes. Our asymptotic theory is based on the distribution of the
randomization without imposing any modeling assumptions on the potential outcomes. Furthermore, under rerandomization, the asymptotic distribution in (2) has a clear geometric interpretation
as displayed in Figure 1, in which we fix Vτ τ at 1 without loss of generality, θ is the angle between
√
√
n(τ̂Y − τY ) and its projection on the space spanned by nτ̂X , and then R is the cosine of θ.
8

√
Figure 1: Geometry of rerandomization. O is the origin, θ is the angle between n(τ̂
√ Y − τY ) and
√
its projection on nτ̂X . The ellipse around O is the acceptance region√under ReM. 1 − R2 · ε0 is
√
√
the component of n(τ̂Y − τY ) orthogonal to the space of nτ̂X , and R2 · LK,a is the projection
√
onto the space of nτ̂X under ReM.

3.4.

Representation and simulation of the asymptotic distribution under ReM

The asymptotic distribution in (2) involves a random variable LK,a that does not appear in standard
statistical problems. Algebraically, LK,a ∼ D1 | D 0 D ≤ a is the first coordinate of a K dimensional
standard Gaussian vector, subject to the constraint that the squared length of the vector does not
exceed a. This type of truncation of Gaussian distributions is apparently unstudied except for
Tallis (1963) and Morgan and Rubin (2012). Because the standard Gaussian vector is spherically
symmetric (Dempster 1969; Rubin 1976; Fang et al. 1989), it can be written as a product of two
independent random components, a χK random variable and a random vector uniformly distributed
on the (K − 1) dimensional unit sphere. The truncation condition, D 0 D ≤ a, affects only the
first component χK , leaving the second component unchanged. Basic properties of spherically
symmetrical distributions allow us to represent LK,a using some known distributions, which allows
for easy simulation of LK,a .
Let χ2K,a ∼ χ2K | χ2K ≤ a be a truncated χ2 random variable, UK the first coordinate of the
uniform random vector over the (K − 1) dimensional unit sphere, S a random sign taking ±1 with
probability 1/2, and βK ∼ Beta (1/2, (K − 1)/2) a Beta random variable degenerating to a point
mass at 1 when K = 1.

9

Proposition 3. LK,a can be represented as
LK,a ∼ D1 | D 0 D ≤ a ∼ χK,a UK ∼ χK,a S

p
βK ,

(3)

where (χK,a , UK ) are mutually independent, and (χK,a , S, βK ) are mutually independent. LK,a is
symmetric and unimodal around zero, with variance Var(LK,a ) = vK,a = P (χ2K+2 ≤ a)/P (χ2K ≤
a) < 1.
Because both ε0 and LK,a are symmetric and both are unimodal at zero, their linear combination
is also symmetric and unimodal at zero according to Wintner (1936)’s Theorem. The same is true
√
for the asymptotic distribution of n(τ̂Y − τY ) in (2). The representation in (3) allows for easy
√
simulation of LK,a , as well as the asymptotic distribution of n(τ̂Y − τY ) in (2), which is relevant
for statistical inference discussed later.
Without loss of generality, we fix Vτ τ at 1, and consider the distribution of

Q=

p

1 − R 2 · ε0 +

√

R2 · LK,a ,

(4)

which depends on R2 , the dimension of the covariates K, and the asymptotic acceptance probability
of rerandomization pa = P (χ2K ≤ a). We simulate values of Q using independent and identically
distributed (i.i.d) draws from (4). First, we fix K = 10 and pa = 0.001. Figure 2a shows the
probability densities of Q with different values of R2 , which approaches to that of LK,a as R2
increases. Because ε0 is more diffusely distributed than the truncated variable LK,a , the probability
density of Q will concentrate more around 0 with increasing R2 , as shown in Figure 2a.
Second, we fix K = 3 and R2 = 0.6. Figure 2b shows the probability densities of Q with
different values of asymptotic acceptance probability pa ; the CRE corresponds to pa = 1. With
smaller pa , the distribution of Q becomes more concentrated around 0. Asymptotically, using
smaller acceptance probabilities in ReM gives us more precise estimators for the average causal
effect. However, when R2 < 1, which is usually the case in practice, the gain of ReM by decreasing
the threshold a becomes less as a becomes smaller. For example, the density of Q with pa = 0.0001
is almost the same as the one with pa = 0.001 in Figure 2b, and the percentage reduction in variance
of Q achieved by decreasing pa from 0.001 to 0.0001 is only 5.7%.
10

0.6
0.5

pa = 1
pa = 0.5
pa = 0.1
pa = 0.001

−4

−2

0
2
1 − R2 ⋅ ε0+ R2 ⋅ LK,a

0.0

0.1

0.2

0.2

Density
0.4

Density
0.3
0.4

0.6

0.8

R2 = 0
R2 = 0.3
R2 = 0.5
R2 = 0.7
R2 = 0.9

4

−4

−2

0
2
1 − R2 ⋅ ε0+ R2 ⋅ LK,a

4

(b) R2 = 0.6 and K = 3

(a) K = 10 and pa = 0.001

Figure 2: Asymptotic distribution under ReM with Vτ τ fixed at 1

3.5.

Asymptotic unbiasedness, sampling variance and quantile ranges

Theorem 1 characterizes the asymptotic behavior of τ̂Y over ReM, which immediately implies the
following conclusions as extensions of Morgan and Rubin (2012).
First, the asymptotic distribution in (2) is symmetric around 0, implying that τ̂Y is asymptotically unbiased for τY . Let Ea (·) and Vara (·) denote the expectation and covariance matrix (or
variance for scalar cases) of the asymptotic sampling distribution of a sequence of random vectors.
√
√
Corollary 1. Under ReM and Condition 1, Ea { n(τ̂Y − τY ) | nτ̂X ∈ M} = 0.
Morgan and Rubin (2012) gave a counter-example showing that, in an experiment with unequal
treatment group sizes, τ̂Y can be biased for τY over ReM. As conjectured by Morgan and Rubin
(2015), our result suggests that the bias is often small with large samples. Corollary 1 extends
Morgan and Rubin (2012, Theorem 2.1) and ensures the asymptotic unbiasedness of τ̂Y for experiments with any ratio of group sizes. Corollary 1 also implies that any covariate asymptotically has
the same means under treatment and control.
Furthermore, from Proposition 3 and Theorem 1, we can calculate the asymptotic sampling
variances of τ̂X and τ̂Y , and the percentage reductions in asymptotic sampling variances (PRIASV)
under ReM compared to CRE. Recalling that vK,a = P (χ2K+2 ≤ a)/P (χ2K ≤ a), we summarize the
results below.
11

Corollary 2. Under ReM and Condition 1, the asymptotic sampling covariance of τ̂X is

Vara

√

nτ̂X |

√


nτ̂X ∈ M = vK,a Vxx ,

and the PRIASV of any component of τ̂X is 1 − vK,a . The asymptotic sampling variance of τ̂Y is

Vara

√

n(τ̂Y − τY ) |

√


nτ̂X ∈ M = Vτ τ 1 − (1 − vK,a )R2 ,

(5)

and the PRIASV of τ̂Y is (1 − vK,a )R2 .
Note that the asymptotic sampling covariance and sampling variance of τ̂X and τ̂Y are actually
the limits of vK,a Vxx and Vτ τ {1 − (1 − vK,a )R2 } in the sequence of finite populations. However, for
descriptive convenience, we omit these limit signs when discussing the expectation and covariance
of asymptotic sampling distributions. When a is close to 0, or equivalently the asymptotic acceptance probability is small, the asymptotic sampling variance on the right hand side of (5) reduces
to Vτ τ (1 − R2 ), which is identical to the asymptotic sampling variance of the regression adjusted
estimator under CRE discussed in Lin (2013) as an extension of Fisher (1925, 1935). Therefore,
rerandomization does covariate adjustment in the design stage, and regression does covariate adjustment in the analysis stage. Cox (2009) and Morgan and Rubin (2012) discussed related issues.
When the causal effect is additive, R2 is equal to the finite population squared multiple correlation between X and Y (0). Therefore, Corollary 2 is an asymptotic version of Theorem 3.2 in
Morgan and Rubin (2012).
Under ReM, in addition to the sampling variance reduction result concerning τ̂Y in Corollary 2,
we consider the reduction in the length of the (1 − α) quantile range of τ̂Y compared to that under
CRE. We choose the length of the (1 − α) quantile range, because of its connection to constructing
confidence intervals as discussed shortly.
Let zξ be the ξth quantile of a standard Gaussian distribution. Let νξ (R2 , pa , K) be the ξth
quantile of the distribution of Q in (4). Note that νξ (0, pa , K) = zξ . Because pa and K are usually
known by design, we write νξ (R2 , pa , K) as νξ (R2 ) for notational simplicity. Under ReM, the (1−α)

12

quantile range of the asymptotic distribution of

√

n(τ̂Y − τY ) is

h
p
p i
QRα (Vτ τ , R2 ) = να/2 (R2 ) Vτ τ , ν1−α/2 (R2 ) Vτ τ ,

(6)

and the corresponding quantile range under CRE is
h
p i
p
QRα (Vτ τ , 0) = zα/2 Vτ τ , z1−α/2 Vτ τ .

(7)

Theorem 2. Under Condition 1, the length of the (1−α) quantile range of the asymptotic sampling
√
distribution of n(τ̂Y − τY ) under ReM is less than or equal to that under CRE, with the difference
nondecreasing in R2 .

3.6.

Sampling variance estimation and confidence intervals

Asymptotic sampling variance and quantile range for τ̂Y depend on Vτ τ and R2 , which are determined by the finite population covariances among potential outcomes and covariates. To obtain a
sampling variance estimator and to construct an asymptotic confidence interval for τY , we need to
estimate these finite population variances and covariances. Let s2Y (z) , s2Y (z)|X and sY (z),X be the
sample variance of Y , sample variance of linear projection of Y on X, and sample covariance of
Y and X in treatment arm z. We show in the Supplementary Material that under ReM they are
asymptotically unbiased for their population analogues SY2 (z) , SY2 (z)|X and SY (z),X . Therefore, we
can estimate Vτ τ by (Ding et al. 2016)
2 −1
V̂τ τ = r1−1 s2Y (1) + r0−1 s2Y (0) − (sY (1),X − sY (0),X )(SX
) (sX,Y (1) − sX,Y (0) ).

We then estimate R2 by
n
 2 −1
o
R̂2 =V̂τ−1
r1−1 s2Y (1)|X + r0−1 s2Y (0)|X − sY (1),X − sY (0),X SX
sX,Y (1) − sX,Y (0) .
τ

(8)

We set R̂2 to be 0 if the estimator in (8) is negative.
According to (5), we can estimate the sampling variance of τ̂Y by V̂τ τ {1 − (1 − vK,a )R̂2 }/n,
and according to (6), we can construct a large sample (1 − α) confidence interval for τY using
13

√
τ̂Y −QRα (V̂τ τ , R̂2 )/ n. The sampling variance estimator is smaller than Neyman (1923)’s sampling
variance estimator for CRE, and the confidence interval is shorter than Neyman (1923)’s confidence
interval for CRE. Not surprisingly, unless the residual from the linear projection of individual causal
effect on the covariates is constant, the above sampling variance estimator and confidence interval
are both asymptotically conservative, in the sense that the probability limit of variance estimator
is larger than or equal to the actual sampling variance, and the limit of coverage probability of
confidence interval is larger than or equal to (1 − α). Therefore, if we conduct ReM in the design
stage but analyze data as in CRE, the consequential sampling variance estimator and confidence
intervals will be overly conservative. These results are all intuitive, and we present the algebraic
details for the proofs of these results in the Supplementary Material, where the unimodality of LK,a
plays an important role in the conservativeness of confidence intervals. Interestingly, as shown in
the Supplementary Material, we do not need more moment conditions beyond Condition 1 to ensure
the asymptotic properties of the variance estimator and the confidence intervals.
We also conduct simulations in the Supplementary Material with non-additive and additive
causal effects, where the results agree with our theory for ReM.

4.
4.1.

Rerandomization with tiers of covariates
Mahalanobis distance with tiers of covariates criterion

When covariates are thought to have different levels of importance for the outcomes, Morgan and
Rubin (2015) proposed rerandomization using the Mahalanobis distance with differing criteria for
different tiers of covariates. We partition the covariates into T tiers indexed by t = 1, . . . , T with
decreasing importance, with kt covariates in tier t. Let Xi = (Xi [1], . . . , Xi [T ]), where Xi [t] denotes
the covariates in tier t. Define Xi [t] = (Xi [1], . . . , Xi [t]), the covariates in the first t tiers. Following
2
the notation in Morgan and Rubin (2015), we let SX[t−1]
be the finite population covariance matrix

of the covariates in first t − 1 tiers, and SX[t],X[t−1] be the finite population covariance matrix
between X[t] and X[t − 1]. We first apply a block-wise Gram–Schmidt orthogonalization to the
covariates to create the orthogonalized covariates:

Ei [1] = Xi [1],
14


−1
2
Ei [t] = Xi [t] − SX[t],X[t−1] SX[t−1]
Xi [t − 1],

(2 ≤ t ≤ T )

where Ei [t] is the residual of the projection of the covariates Xi [t] in tier t onto the space spanned
by the covariates in previous tiers; Ei = (Ei [1], . . . , Ei [T ]). Let τ̂E[t] be the difference-in-means
2
vector of Ei [t] between treatment and control groups, and SE[t]
the finite population covariance

matrix of Ei [t]. The Mahalanobis distance in tier t is

Mt =

n1 n0 0  2 −1
τ̂
SE[t]
τ̂E[t] ,
n E[t]

and rerandomization using the Mahalanobis distance with tiers of covariates (ReMT) accepts those
treatment assignments with Mt ≤ at , where at ’s are predetermined constants (1 ≤ t ≤ T ). We can
√
show that the criterion depends only on nτ̂X and Vxx . If T = 1, then ReMT is simply ReM. We
√
use T to denote the acceptance region for nτ̂X under ReMT. The theory below extends Morgan
and Rubin (2015) using the concepts from our Section 3.

4.2.

Multiple correlation between potential outcomes and covariates with tiers

Similar to Section 3.2, we define the finite population squared multiple correlation between the
potential outcome Y (z) and the orthogonalized covariates in tier t as ρ2t (z), and the finite population
squared multiple correlation between the individual causal effect and the orthogonalized covariates
in tier t as ρ2t (τ ). We further define an R2 -type measure as the function of these finite population
quantities and the proportions of group sizes:

ρ2t

=

SY2 (1)
r1 Vτ τ

ρ2t (1)

+

SY2 (0)
r0 Vτ τ

ρ2t (0) −

Sτ2 2
ρ (τ ),
Vτ τ t

(1 ≤ t ≤ T )

which under the additive causal effect assumption reduces to ρ2t = ρ2t (1) = ρ2t (0), the squared
multiple correlation between E[t] and Y (1) or Y (0).
Under CRE, ρ2t is the sampling squared multiple correlation between τ̂Y and τ̂E[t] , and can be
equivalently written as

ρ2t

= Corr(τ̂Y , τ̂E[t] ) =

r1−1 SY2 (1)|E[t] + r0−1 SY2 (0)|E[t] − Sτ2|E[t]
r1−1 SY2 (1) + r0−1 SY2 (0) − Sτ2
15

,

(1 ≤ t ≤ T )

where SY2 (z)|E[t] and Sτ2|E[t] are the finite population variances of the projections of the potential
outcomes and individual causal effects on the orthogonalized covariates in tier t. For descriptive
P
simplicity, we introduce ρ2T +1 = 1 − Tt=1 ρ2t = 1 − R2 for later discussion.

4.3.

Asymptotic distribution of τ̂Y

√
0 ) in (1) still holds for ReMT, with region M replaced
The weak convergence of n(τ̂Y − τY , τ̂X
√
by region T . Intuitively, n(τ̂Y − τY ) can be decomposed into (T + 1) parts: the part unrelated
to covariates and the T projections onto the space spanned by the orthogonalized covariates in T
tiers. Due to the construction of the orthogonalized covariates, these (T + 1) parts are orthogonal
to each other and the constraint for balance on the Mahalanobis distance in tier t affects only the
t-th projection.
As earlier, let ε0 ∼ N (0, 1), and extending earlier notation, let Lkt ,at ∼ Dt1 | Dt0 Dt ≤ at , where
Dt = (Dt1 , . . . , Dtkt ) ∼ N (0, Ikt ) for 1 ≤ t ≤ T .
Theorem 3. Under ReMT and Condition 1,
√

n(τ̂Y − τY ) |

√

nτ̂X ∈ T

.

p
∼
Vτ τ

ρT +1 · ε0 +

T
X

!
ρt · Lkt ,at

,

(9)

t=1

where (ε0 , Lk1 ,a1 , . . . , LkT ,aT ) are mutually independent.
Obviously, in (9), ε0 is the part of

√

n(τ̂Y − τY ) that is unrelated to the covariates, and Lkt ,at

is the part related to the orthogonalized covariates Ei [t] in tier t. According to Proposition 3, the
distribution in Theorem 3 involves distributions that are easy to simulate.

4.4.

Asymptotic unbiasedness, sampling variance and quantile ranges

Theorem 3 characterizes the asymptotic behavior of

√

n(τ̂Y − τY ) under ReMT, which extends

Morgan and Rubin (2015) as follows.
First, the asymptotic distribution in (9) is symmetric around 0, implying that τ̂Y is asymptotically unbiased for τY . Therefore, all observed or unobserved covariates have asymptotically
balanced means.
√
√
Corollary 3. Under ReMT and Condition 1, Ea { n(τ̂Y − τY ) | nτ̂X ∈ T } = 0.
16

The asymptotic sampling variance of τ̂X under ReMT has a complicated but conceptually
obvious form, and we give it in the Supplementary Material. Below we present only the PRIASV
of τ̂Y ; the PRIASVs for covariates are special cases of the same corollary because covariates are
formally “outcomes” unaffected by the treatment. Recall the definition of vkt ,at = P (χ2kt +2 ≤
at )/P (χ2kt ≤ at ).
Corollary 4. Under ReMT and Condition 1, the asymptotic sampling variance of τ̂Y is

Vara

√

n(τ̂Y − τY ) |

√

T
X
1−
(1 − vkt ,at )ρ2t

(
nτ̂X ∈ T

=Vτ τ

)
,

(10)

t=1

and the PRIASV of τ̂Y is

PT

t=1 (1

− vkt ,at )ρ2t .

When the causal effect is additive, ρ2t becomes the finite population squared multiple correlation
between E[t] and Y (0). Therefore, Corollary 4 is an asymptotic extension of Morgan and Rubin
(2015, Theorem 4.2). When the thresholds at ’s are close to zero, the asymptotic sampling variance
P
on the right hand side of (10) reduces to Vτ τ (1 − Tt=1 ρ2t ) = Vτ τ (1 − R2 ), which is identical to that
of the regression adjusted estimator under CRE (Lin 2013).
We now compare the quantile range under ReMT to that under CRE. Let νξ (ρ21 , ρ22 , . . . , ρ2T ) be
P
the ξth quantile of ρT +1 ε0 + Tt=1 ρt Lkt ,at . Although νξ (ρ21 , ρ22 , . . . , ρ2T ) depends also on pat and kt
(1 ≤ t ≤ K), we omit them to avoid notational clatter. The (1−α) quantile range of the asymptotic
√
distribution of n(τ̂Y − τY ) under ReMT is
h
p
p i
QRα (Vτ τ , ρ21 , . . . , ρ2T ) = να/2 (ρ21 , . . . , ρ2T ) Vτ τ , ν1−α/2 (ρ21 , . . . , ρ2T ) Vτ τ .

(11)

The stronger the correlation between the outcome and the orthogonalized covariates in tier t,
the more reduction in quantile range we have when using ReMT rather than CRE. The following
theorem is immediate.
Theorem 4. Under Condition 1, the (1 − α) quantile range of the asymptotic distribution of
√
n(τ̂Y − τY ) under ReMT is narrower than, or equal to the one under CRE, and the reduction in
length of the quantile range is nondecreasing in ρ2t for all 1 ≤ t ≤ T .

17

4.5.

Sampling variance estimation and confidence interval

We can estimate Vτ τ and ρ2t (1 ≤ t ≤ T ) in the same way as in ReM, and we estimate ρ2T +1 by
1 − R̂2 . In practice, we set ρ̂2t (1 ≤ t ≤ T ) to 0 when it is negative due to sampling variability, and
standardize their sum to R̂2 . According to (10) and (11), we can estimate the sampling variance
of τ̂Y and (1 − α) confidence intervals for τY by replacing the unknown quantities with their point
estimates. The sampling variance estimator is smaller than Neyman (1923)’s sampling variance
estimator for CRE, and the confidence interval is shorter than Neyman (1923)’s confidence interval
for CRE; both are asymptotically conservative in general, and only when the residual from the linear
projection of individual causal effect on the covariates is constant, are they asymptotically exact.
Therefore, analyzing data from ReMT as from CRE, the resulting sampling variance estimator
and confidence intervals are overly conservative. These intuitive statements appear to require
notationally lengthy proofs, which are relegated to the Supplementary Material. Specifically, the
proof for the conservativeness of confidence intervals utilizes the unimodality of the Lkt ,at ’s.

5.

Rerandomization with more general balance criterion

5.1.

More general balance criterion

As pointed out by Morgan and Rubin (2012), the criterion can be any accept-reject function of
the treatment assignment and covariate balance. We can always use randomization tests for a
sequence of sharp null hypotheses, and thereby construct fiducial confidence intervals by inverting
these randomization tests (e.g., under the additive causal effects assumption). In this section, we
discuss the repeated sampling properties of the difference-in-means estimator, where we consider
√
covariate balance criteria that depend only on nτ̂X and Vxx , including ReM and ReMT as special
√
cases, and write the binary covariate balance indicator function as φ( nτ̂X , Vxx ). Let G denote
the acceptance region for rerandomization with the general covariate balance criterion φ (ReG),
√
i.e., nτ̂X ∈ G = {µ : φ(µ, Vxx ) = 1}.
For technical reasons, we require φ to satisfy the following conditions. First, φ is almost surely
continuous. Second, Var(B | φ(B, Vxx ) = 1), as a function of Vxx with B ∼ N (0, Vxx ), is
continuous for all Vxx > 0. Third, P (φ(B, Vxx ) = 1) > 0, for all Vxx > 0 with B ∼ N (0, Vxx ).
18

Fourth, φ(µ, Vxx ) = φ(−µ, Vxx ) for all µ and Vxx > 0. The first two conditions impose certain
smoothness on φ, and the third condition prevents the acceptance region from being a set of measure
zero. The fourth condition imposes symmetry considerations, because relabeling the treatment and
control units should not change the balance. Both ReM and ReMT satisfy these conditions. Below,
we summarize theory in parallel with Sections 3 and 4.

5.2.

Asymptotic sampling properties

−1 B,
The weak convergence in (1) holds with M replaced by G. The projection of A onto B is Vτ x Vxx
−1 B. Therefore, the asymptotic distribution of τ̂ is
and the residual is ε = A − Vτ x Vxx
Y

√

n(τ̂Y − τY ) |

√

.

−1
nτ̂X ∈ G ∼ ε + Vτ x Vxx
B | B ∈ G,

(12)

where ε ∼ N (0, Vτ τ (1 − R2 )) is independent of B ∼ N (0, Vxx ). The asymptotic distribution in
(12) can be easily simulated.
The symmetry condition φ(µ, Vxx ) = φ(−µ, Vxx ) implies that the distribution in (12) is symmetric around 0, which further implies that τ̂Y is asymptotically unbiased for τY . Viewing covariates
as outcomes unaffected by the treatment, all observed or unobserved covariates asymptotically have
the same means in treatment and control groups.

5.3.

Advice for the investigator

Because Vxx is known in finite population inference, we can choose ReGs that result in better
√
√
covariate balance before the physical experiments. Because Vara { nτ̂X | nτ̂X ∈ G} = Var(B |
B ∈ G), it is important to check whether Var(B | B ∈ G) ≤ Var(B) holds, i.e., Var(B) − Var(B |
B ∈ G) is positive semi-definite. This ensures that ReG using acceptance region G reduces the
sampling covariance matrix of the difference-in-means of the covariates. Otherwise, we need to
change the criterion. If Vxx,φ ≡ Var(B | B ∈ G) ≤ Var(B) = Vxx holds, then we can derive the
PRIASV for all observed covariates under ReG. From the decomposition

Vara

√
√
−1
−1
Vxx,φ Vxx
Vxτ ,
n(τ̂Y − τY ) | nτ̂X ∈ G = Vτ τ (1 − R2 ) + Vτ x Vxx

19

we can derive the PRIASV of

√

n(τ̂Y − τY ) under ReG. However, without imposing further condi-

tions on φ, there is no guarantee that the asymptotic quantile range of τ̂ will be shorter under ReG
than that under CRE. Therefore, in the design stage of the experiment, we recommend choosing
balance criteria that are expected to lead to both variance and quantile range reductions, such as
ReM and ReMT.

5.4.

Sampling variance estimation and confidence interval

We can show that s2Y (z) , s2Y (z)|X and s2Y (z),X are asymptotically unbiased for SY2 (z) , SY2 (z)|X and
SY2 (z),X under ReG. Therefore, we can unbiasedly estimate Vτ τ by V̂τ τ as earlier, and R2 by R̂2 in the
same form as (8). We can then estimate Vτ x = r1−1 SY (1),X + r0−1 SY (0),X and the sampling variance
of τbY by replacing the unknown quantities with their point estimates. We can then estimate the
−1 B | B ∈ G is unimodal, as
asymptotic distribution in (12), as well as its quantile ranges. If Vτ x Vxx

in ReM and ReMT, the final confidence interval is asymptotically conservative. We relegate more
details to the Supplementary Material.

6.

An education example with tiers of covariates

We illustrate our theory using the data from the Student Achievement and Retention Project
(Angrist et al. 2009), a randomized evaluation of academic services and incentives at one of the
satellite campuses of a large Canadian university, involving college freshmen. A treatment group
of 150 students was offered an array of support services and substantial cash awards for meeting
a target first year grade point average (GPA), and a control group of many more (1006) students
received only standard university support services.
To illustrate the benefit of rerandomization, we use the 15 covariates as listed in Table 1, and
exclude students with missing values, resulting in 118 students in the treatment group and 856 in
the control. To make the simulation relevant to the real data, we fix unknown parameters based on
some simple model fitting: We fit a linear regression of the observed first year GPA on the treatment
indicator, all covariates and their interactions, and use the fitted model to generate all potential
outcomes under non-additivity. Note that the generating models for the potential outcomes are
not linear in the covariates themselves. To make the data generating process realistic, we simulate
20

Table 1: Covariates in the Student Achievement and Retention Project. The numbers of covariates
in these three tiers are 1, 4 and 10, and the thresholds are a1 = 0.016, a2 = 1.064 and a3 = 4.865.
Tier
1
2
3

Covariates
high school GPA
whether lives at home, gender, age,
whether rarely puts off studying for tests
whether mother is a college graduate, whether mother is a high school graduate,
mother tongue (English or other), whether plans to work while in school,
whether father is a college graduate, whether father is a high school graduate,
whether never puts off studying for tests, whether wants more than a bachelor degree,
whether intends to finish in 4 years, whether at the first choice school

eight pseudo sets of potential outcomes using the fitted model with different choices for the variance
of the residuals. The error terms for Y (1) and Y (0) are independent, and therefore conditional
on the covariates, the potential outcomes are simulated as uncorrelated, but they have a positive
correlation marginally. The final potential outcomes are all truncated to lie on [0, 4], mimicking
the value of the GPA. We choose different variances of residuals such that the values of R2 for the
eight simulated data sets are located approximately evenly within interval [0, 0.5]. One choice for
the variance of residuals is the one estimated from the fitted linear model, and the corresponding
R2 is about 0.23.
Table 1 partitions the covariates into three tiers with decreasing a priori importance to the
outcome. As suggested by Morgan and Rubin (2015), for tiers with increasing numbers of covariates,
we choose at such that P (χ2kt ≤ at ) = (0.001)1/3 = 0.1 for t = 1, 2, 3. We simulate data under
ReMT, and obtain the confidence intervals based on our asymptotic theory for ReMT and Neyman
(1923)’s results for CRE. Figure 3a shows the empirical coverage probabilities of our and Neyman
(1923)’s confidence intervals, showing that Neyman (1923)’s CRE confidence intervals are highly
conservative. Note that there are 15 covariates and only 118 units in the treatment group, and the
sample size is not extremely large. Despite this, our asymptotic confidence interval works well in
this example.
To evaluate the performance of ReMT compared to CRE, we compare the average length of
Neyman (1923)’s confidence interval under CRE with the confidence interval under ReMT. From
Figure 3b, the percentage reduction in average lengths of the confidence intervals under ReMT
compared to Neyman (1923)’s under CRE is nondecreasing in R2 . We also compare the empirical

21

95% quantile range of τ̂Y under ReMT and CRE, and the percentage reduction in the lengths of
quantile ranges are close to the percentage reduction for average lengths of confidence intervals.
When R2 is close to that of the real data set (i.e. 0.23), the percentage increase in the effective
sample size, that is, the sample size needed in CRE in order for τ̂Y to achieve the same 95% quantile
range under ReMT, is about 24%. When R2 is about twice as large as with the real data (i.e. 0.5),

30

the percentage increase in the effective sample size increases to 80%.
ReMT
Neyman ReMT

●

●

●

percentage
10
15
20

●

●
●

●

●

●

●

●

●
●

●
●
●

0

0.93

Reduction in CI
Reduction in QR

25

●

5

coverage probability
0.95
0.97
0.99

●

0.1

0.2

0.3
R2

0.4

0.5

0.1

(a) Empirical coverage probabilities of our
and Neyman (1923)’s 95% confidence
intervals under ReMT

0.2

0.3
R2

0.4

0.5

(b) Percentage reductions of average lengths
of confidence intervals ( ) and quantile
ranges (4) comparing ReMT with CRE

Figure 3: Eight data sets simulated based on the Student Achievement and Retention Project

7.

Conclusions, Connections and Extensions

Extending Morgan and Rubin (2012, 2015), we show, using analysis and simulations, that rerandomization balances covariates better than complete randomization, and provides a more precise
difference-in-means estimator for the average causal effect. The asymptotic distributions of the
difference-in-means estimator under rerandomization with strigent constraints are close to that of
the regression adjusted estimator under CRE (Lin 2013), implying that rerandomization does the
covariate adjustment in the design stage and avoids outcome modeling. The new asymptotic distributions allow us to construct confidence intervals for the average causal effect, when the classical
Neyman (1923)’s inference for CRE is overly conservative.

22

Appendix
In this Appendix, we provide proofs for the asymptotic distribution of

√

n(τ̂Y − τY ) under ReM and

ReMT, and the representation for random variable LK,a . First, we need the following two lemmas.
Lemma 1. Let LK,a ∼ D1 | D 0 D ≤ a, where D = (D1 , . . . , DK )0 ∼ N (0, IK ). For any K
dimensional unit vector h, we have LK,a ∼ h0 D | D 0 D ≤ a.
Proof of Lemma 1. We construct an orthogonal matrix H whose first row is h0 . Then D ∼ HD.
Therefore, LK,a ∼ D1 | D 0 D ≤ a ∼ h0 D | (HD)0 HD ≤ a ∼ h0 D | D 0 D ≤ a.
Lemma 2. Let UK be the first coordinate of the uniform random vector over the (K − 1) dimensional unit sphere. Let S be a random sign taking ±1 with probability 1/2, βK ∼ Beta (1/2, (K − 1)/2)
be a Beta random variable degenerating to a point mass at 1 when K = 1, and (S, βK ) are mutually
√
independent. Then UK ∼ S βK .
Proof of Lemma 2. When K = 1, it is easy to see Lemma 2 holds. When K ≥ 2, let D =
(D1 , D2 , . . . , DK )0 ∼ N (0, IK ).The standardized normal random vector with unit length is uni√
formly distributed on the unit sphere (Fang et al. 1989), and therefore UK ∼ D1 / D 0 D. Let S be
a random sign independent of D; then D1 ∼ S|D1 |. Thus (D1 , D2 , . . . , DK ) ∼ (S|D1 |, D2 , . . . , DK ),

UK

√
√
∼ D1 / D 0 D ∼ S|D1 |/ D 0 D = S ·

s

D12
PK

2
k=1 Dk

,

PK

2 follow i.i.d. Gamma distributions
Because D12 , . . . , DK
P
2
with shape parameter 1/2 and scale parameter 2, D12 /( K
k=1 Dk ) ∼ Beta(1/2, (K − 1)/2) is inde-

and S is independent of D12 /(

2
k=1 Dk ).

pendent of S. Therefore, Lemma 2 holds.
−1 B, which has variance c2 =
Proof of Theorem 1. The linear projection of A on B is Vτ x Vxx
−1 V
2
Vτ x Vxx
xτ = Vτ τ R . The residual from the linear projection of A on B is

−1
ε = A − Vτ x Vxx
B ∼ N (0, (1 − R2 )Vτ τ ) ∼

−1/2

where ε0 ∼ N (0, 1). Let h0 = Vτ x Vxx

p
Vτ τ (1 − R2 ) · ε0 ,
−1/2

/c be the standardized vector of Vτ x Vxx

23

with unit

−1/2

length, and D = Vxx

B ∼ N (0, IK ) be the standardization of B. Then
−1/2

−1
A = ε + Vτ x Vxx
B = ε + Vτ x Vxx

D = ε + ch0 D.

According to Proposition 2 and Lemma 1,
√

n(τ̂Y − τY ) |

√

.

nτ̂X ∈ M ∼ A | B ∈ M
p

√
p
Vτ τ
1 − R2 · ε0 + R2 · LK,a ,
∼ ε + ch0 D | D 0 D ≤ a ∼ ε + cLK,a ∼

where ε, or ε0 , is independent of LK,a .
Proof of Theorem 3. We use Γ to denote the linear transformation from Xi to Ei , i.e. Ei =
ΓXi , where Γ depends only on Vxx . Correspondingly, let G = ΓB = (G01 , G02 , . . . , G0T )0 be the
block-wise Gram–Schmidt orthogonalization of B, where Gt is a kt dimensional random vector.
Let Dt = Var(Gt )−1/2 Gt ∼ N (0, Ikt ) be the standardization of Gt . The linear projection of A
on Gt is Cov(A, Gt )Var(Gt )−1 Gt , with variance c2t = Cov(A, Gt )Var(Gt )−1 Cov(Gt , A). Because
0 ) has the same covariance matrix as (A, G0 ) = (A, ΓB 0 ), the sampling variance
(τ̂Y , τ̂E0 ) = (τ̂Y , Γτ̂X

of linear projection of τ̂Y on τ̂E[t] is the same as the variance of linear projection of A on Gt . The
former is Vτ τ ρ2t , and the latter is c2t . Therefore, c2t = Vτ τ ρ2t . The residual from the linear regression
of A on B (or equivalently on G) is

−1
ε = A − Vτ x Vxx
B =A−

T
X

Cov(A, Gt )Var(Gt )−1 Gt

t=1

∼ N

0, Vτ τ

1−

T
X

!!
ρ2t

∼ N (0, Vτ τ ρ2T +1 ) ∼

p
Vτ τ · ρT +1 ε0 ,

t=1

where ε0 ∼ N (0, 1). Let h0t = Cov(A, Gt )Var(Gt )−1/2 /ct be the standardized vector of Cov(A, Gt )Var(Gt )−1/2
with unit length. Then A has the following decomposition:

A = ε+

= ε+

T
X
t=1
T
X

Cov(A, Gt )Var(Gt )−1 Gt = ε +

T
X
t=1

ct h0t Dt .

t=1

24

Cov(A, Gt )Var(Gt )−1/2 Dt

Because Proposition 2 holds for ReMT with M replaced by T , and (D1 , . . . , DT ) are mutually
independent,
√

n(τ̂Y − τY ) |

√

.

nτ̂X ∈ T ∼ A | B ∈ T ∼ ε +

T
X


ct h0t Dt | Dt0 Dt ≤ at .

t=1

According to Lemma 1,
√

n(τ̂Y − τY ) |

√

.

nτ̂X ∈ T ∼ ε +

T
X

ct Lkt ,at

p
∼ Vτ τ

ρT +1 ε0 +

T
X

!
ρt Lkt ,at

,

t=1

t=1

where Lkt ,at ∼ Dt1 | Dt0 Dt ≤ at , (ε, Lk1 ,a1 , . . . , LkT ,aT ) are mutually independent, and (ε0 , Lk1 ,a1 , . . . , LkT ,aT )
are mutually independent. Therefore, Theorem 3 holds.
Proof of Proposition 3. Let D = (D1 , . . . , DK )0 ∼ N (0, IK ) and LK,a ∼ D1 | D 0 D ≤ a. First,
we show that LK,a is symmetric and unimodal around 0. Let f (·) be the density of standard
Gaussian distribution. It is easy to show that LK,a ∼ D1 | D 0 D ≤ a is symmetric around 0, and
has density
P
2
2
P (χ2K−1 ≤ a − l2 )
f (l)P ( K
k=2 Dk ≤ a − l )
=
f
(l)
p(l) =
P (D 0 D ≤ a)
P (χ2K ≤ a)
Therefore, LK,a is unimodal, because for any |l1 | < |l2 | <

p(l1 ) =f (l1 )

√

√
√
(− a < l < a).

a,

P (χ2K−1 ≤ a − l22 )
P (χ2K−1 ≤ a − l12 )
≥
f
(l
)
= p(l2 ).
2
P (χ2K ≤ a)
P (χ2K ≤ a)

Second, the variance formula for LK,a follows from Morgan and Rubin (2012, Theorem 3.1).
√
Third, we represent LK,a by known distributions. Let RK = D 0 D be the length of vector D,
and D/RK be the normalized vector of D with unit length. From the property of the multivariate
2 ∼ χ2 , D/R follows the uniform distribution on the K − 1 dimensional
Gaussian distribution, RK
K
K

unit sphere, and they are independent (Fang et al. 1989). Let UK be the first coordinate of D/RK ,
2 ≤ a. Because χ
2
then LK,a ∼ UK RK | RK
K,a ∼ RK | RK ≤ a, and χK,a is independent of UK ,
2 ≤ a ∼ χ
we have LK,a ∼ UK RK | RK
K,a UK . According to Lemma 2, we have LK,a ∼ χK,a UK ∼
√
χK,a S βK .

25

REFERENCES
J. Angrist, D. Lang, and P. Oreopoulos. Incentives and services for college achievement: Evidence
from a randomized trial. American Economic Journal: Applied Economics, 1:136–163, 2009.
S. Athey and G.W. Imbens. The econometrics of randomized experiments. Handbook of Economic
Field Experiments, 1:73–140, 2017.
R. A. Bailey. Restricted randomization. Biometrika, 70:183–198, 1983.
M. Bruhn and D. McKenzie. In pursuit of balance: Randomization in practice in development field
experiments. American Economic Journal: Applied Economics, 1:200–232, 2009.
W. G. Cochran. The planning of observational studies of human populations. Journal of the Royal
Statistical Society. Series A (General), 128:234–266, 1965.
W. G. Cochran. Sampling Techniques. New York: Wiley, 3rd edition, 1977.
William G. Cochran and Gertrude M. Cox. Experimental Designs. John Wiley & Sons, Inc., 2
edition, 1992.
D. R. Cox. Randomization and concomitant variables in the design of experiments. In P. R. Krishnaiah G. Kallianpur and J. K. Ghosh, editors, Statistics and Probability: Essays in Honor of
C. R. Rao, pages 197–202. North-Holland, Amsterdam, 1982.
D. R. Cox. Randomization in the design of experiments. International Statistical Review, 77:
415–429, 2009.
A. P. Dempster. Elements of Continuous Multivariate Analysis. Addison-Wesley, Reading, Massachusetts, 1969.
P. Ding, A. Feller, and L. Miratrix. Decomposing treatment effect variation. arXiv preprint
arXiv:1605.06566, 2016.
K. T. Fang, S. Kotz, and K. W. Ng. Symmetric Multivariate and Related Distributions. Chapman
and Hall/CRC, 1989.

26

S. E. Fienberg and D. V. Hinkley. R. A. Fisher : An Appreciation. New York : Springer-Verlag,
1980.
R. A. Fisher. Statistical Methods for Research Workers. Edinburgh: Oliver and Boyd, 1st edition,
1925.
R. A. Fisher. The arrangement of field experiments. Journal of the Ministry of Agriculture of Great
Britain, 33:503–513, 1926.
R. A. Fisher. The Design of Experiments, 1st Edition. Edinburgh, London: Oliver and Boyd, 1935.
R. Greevy, B. Lu, J. H. Silber, and P. Rosenbaum. Optimal multivariate matching before randomization. Biostatistics, 5:263–275, 2004.
P. M. Grundy and M. J. R. Healy. Restricted randomization and quasi-latin squares. Journal of
the Royal Statistical Society, Series B (Methodological), 12:286–291, 1950.
J. Hájek. Limiting distributions in simple random sampling from a finite population. Publications
of the Mathematics Institute of the Hungarian Academy of Science, 5:361–74, 1960.
B. B. Hansen and J. Bowers. Covariate balance in simple, stratified and clustered comparative
studies. Statistical Science, 23:219–236, 2008.
P. W. Holland. Statistics and causal inference. Journal of the American statistical Association, 81:
945–960, 1986.
G. W. Imbens and D. B. Rubin. Causal Inference for Statistics, Social, and Biomedical Sciences:
An Introduction. Cambridge: Cambridge University Press, 2015.
L. Keele, C. McConnaughy, I. White, P. M. E. M. List, and D. Bailey. Adjusting experimental
data. In Experiments in Political Science Conference. 2009.
A. M. Krieger, D. Azriel, and A. Kapelner. Nearly Random Designs with Greatly Improved Balance.
arXiv preprint arXiv:1612.02315, 2016.
E. L. Lehmann. Fisher, Neyman, and the Creation of Classical Statistics. New York: Springer,
2011.
27

X. Li and P. Ding. General forms of finite population central limit theorems with applications to
causal inference. Journal of the American Statistical Association, in press, 2016.
W. Lin. Agnostic notes on regression adjustments to experimental data: Reexamining Freedman’s
critique. The Annals of Applied Statistics, 7:295–318, 2013.
K. L. Morgan and D. B. Rubin. Rerandomization to improve covariate balance in experiments.
The Annals of Statistics, 40:1263–1282, 2012.
K. L. Morgan and D. B. Rubin. Rerandomization to balance tiers of covariates. Journal of the
American Statistical Association, 110:1412–1421, 2015.
J. Neyman. On the application of probability theory to agricultural experiments. essay on principles
(with discussion). section 9 (translated). reprinted ed. Statistical Science, 5:465–472, 1923.
D. B. Rubin. Multivariate matching methods that are equal percent bias reducing, I: some examples.
Biometrics, 32:109–120, 1976.
D. B. Rubin. Comment to W. R. Shadish, M. H. Clark and P. M. Steiner. Journal of the American
Statistical Association, 103:1350–1353, 2008.
L. J. Savage. The Foundations of Statistical Inference. Methuen and Co. Led., London, 1962.
T. P. Speed. Introduction to Fisher (1926). In S. Kotz and N. L. Johnson, editors, Breakthroughs
in Statistics, pages 71–81. Springer, 1992.
Student. Comparison between balanced and random arrangements of field plots. Biometrika, 29:
363–378, 1938.
G. M. Tallis. Elliptical and radial truncation in normal populations. The Annals of Mathematical
Statistics, 34:940–944, 1963.
A. Wintner. On a class of Fourier transforms. American Journal of Mathematics, 58:45–90, 1936.
F. Yates. Comment to F. J. Anscombe. Journal of the Royal Statistical Society, Series A (General),
111:204–205, 1948.
W. J. Youden. Randomization and experimentation. Technometrics, 14:13–22, 1972.
28

Supplementary Material
Section A1 uses simulations to evaluate the asymptotic approximations for the sampling distributions of τ̂ , as well as the coverage probablilities of 95% confidence intervals for τ under ReM.
√
0 ) under ReG, the asymptotic unbiasedSection A2 shows the weak convergence of n(τ̂Y − τY , τ̂X
ness of τ̂Y and balance in means of all covariates, and the formula for the sampling squared multiple
correlation between τ̂Y and τ̂X under CRE in Proposition 1. Section A3 shows the percentage reductions in asymptotic sampling variances and lengths of quantile ranges under rerandomization.
Section A4 shows the asymptotic conservativeness of sampling variance estimators and confidence
intervals.

A1.

Numerical Examples

We conduct numerical examples where the group sizes are very different and the potential outcomes
are simulated from a nonlinear model with binary values. Let r1 = 0.1, r0 = 0.9 and K = 3, so the
two treatment group sizes are very different. Let β1 = (2, 3, 4)0 and β0 = (0, 1, 1)0 . The covariates
for all units are i.i.d samples from Xk ∼ Bernoulli(0.5) for 1 ≤ k ≤ K, where (X1 , . . . , XK ) are
mutually independent. The binary potential outcomes are i.i.d samples from:

Y (z) = I z + βz0 (X − 0.51K ) + δz ≥ 0 ,

where I(·) is the indicator function.

i.i.d

δ0 , δ1 ∼ N (0, 1), (z = 0, 1),

(A1)

We simulate three data sets with different sample sizes

(1000, 3000, 5000), and the causal effects for these simulated data sets are not additive. For example, when n = 1000, τY = 0.121, SY2 (1) = 0.24, SY2 (0) = 0.25, and Sτ2 = 0.33. Figure A1 shows the
√
histograms of n(τ̂Y − τY ) under both ReM with pa = 0.001 and CRE, based on 105 rerandomizations and 105 complete randomizations, as well as their asymptotic approximations using (2) and
Gaussian distributions. Although the potential outcomes models are not linear, the asymptotic
distributions are close to their corresponding theoretical repeated sampling distributions, and the
asymptotic approximations become better as the sample sizes increase.
In the above numerial example, the causal effects are not additive. Figure A2 shows the empirical
coverage probabilities of the 95% confidence intervals with different sample sizes (1000, 2000, 3000, 4000, 5000)
1

0
2
n(τ^Y − τY)

4

6

0.4
0.2

Density

−6

−4

(a) n = 1000

−2

0
2
n(τ^Y − τY)

4

6

(b) n = 3000

0.0

0.1

0.2
0.1

−2

0.0

−4

histogram−ReM
histogram−CRE
density−ReM
density−CRE

0.3

0.4

Density

0.2

Density

0.1
0.0

−6

histogram−ReM
histogram−CRE
density−ReM
density−CRE

0.3

0.3

0.4

histogram−ReM
histogram−CRE
density−ReM
density−CRE

−6

−4

−2

0
2
n(τ^Y − τY)

4

6

(c) n = 5000

√
Figure A1: Histograms and asymptotic densities of n(τ̂Y − τY ). The grey and white histograms
are the empirical distributions, and the solid and dotted lines are the asymptotic densities under
ReM and CRE.
and treatment and control proportions (r1 , r0 ) = (0.1, 0.9). We also generate other data sets with
additive causal effects, in which the data generating process is the same as (A1) except that Yi (0)
is replaced by Yi (1) − τY . As anticipated, with non-additive causal effects, the empirical coverage probabilities are larger than 95%, but with additive causal effects, the empirical coverage

0.99

probabilities are close to 95%.

●

non−additive
additive
●

●

3000

4000

●

●

0.94

0.95

0.96

0.97

coverage probability

0.98

●

1000

2000

sample size

5000

Figure A2: Empirical coverage probabilities of 95% confidence intervals with (4) and without (
additivity

2

)

A2.

More Details on Weak Convergence

√
0 ) under rerandomization with general
We consider the asymptotic distribution of n(τ̂Y − τY , τ̂X
√
balance criterion φ( nτ̂X , Vxx ) satisfying the conditions in Section 5.1:
(1) φ is almost surely continuous.
(2) Var(B | φ(B, Vxx ) = 1), as a function of Vxx with B ∼ N (0, Vxx ), is continuous for all
Vxx > 0.
(3) P (φ(B, Vxx ) = 1) > 0, for any Vxx > 0 with B ∼ N (0, Vxx ).
(4) φ(µ, Vxx ) = φ(−µ, Vxx ), for all µ and Vxx > 0.
We write the limit of V as




 Vτ τ
V∞ = lim V = lim 
n→∞
n→∞
Vxτ





Vτ x   Vτ τ,∞ Vτ x,∞ 
=
,
Vxx
Vxτ,∞ Vxx,∞

which is assumed to be positive definite. Let G = {µ : φ(µ, Vxx ) = 1} be the acceptance region for
√
nτ̂X , and G∞ = {µ : φ(µ, Vxx,∞ ) = 1} be its limit.
Proposition A1. Under ReG, as n → ∞,


√
 n(τ̂Y − τY )
 √

nτ̂X


√



d
 A∞ 
nτ̂X ∈ G −→ 

B∞

B∞ ∈ G∞ ,

0 ) ∼ N (0, V ), in the sense that, for any continuity set C of (A , B 0 ) | B
where (A∞ , B∞
∞
∞
∞ ∈ G∞ ,
∞

P

√

0
n(τ̂Y − τY , τ̂X
)∈C|

√


0
nτ̂X ∈ G → P (A∞ , B∞
) ∈ C | B∞ ∈ G∞ .

Proof of Proposition A1. According to the finite population central limit theorem,
√

n(τ̂Y − τY ),

√

0
nτ̂X
, Vxx



3

d

−→


0
A∞ , B∞
, Vxx,∞ .

The continuous mapping theorem implies
 d
√
√ 0
√
n(τ̂Y − τY ), nτ̂X
, φ( nτ̂X , Vxx ) −→


0
A∞ , B∞
, φ(B∞ , Vxx,∞ ) .

(A2)

Let I = (0.5, 1.5) ∈ R1 be an open interval. Becasue φ is a 0-1 function, I is a continuity set of
φ(B∞ , Vxx,∞ ), in the sense that P (φ(B∞ , Vxx,∞ ) ∈ ∂I) = 0. According to (A2) and Portmanteau’s
Theorem, as n → ∞,
 √
 √
P φ( nτ̂X , Vxx ) = 1 = P φ( nτ̂X , Vxx ) ∈ I → P {φ(B∞ , Vxx,∞ ) ∈ I} = P {φ(B∞ , Vxx,∞ ) = 1} .
0 ) | φ(B , V
For any continuity set C ∈ RK+1 of (A∞ , B∞
∞
xx,∞ ) = 1, C × I is also a continuity set
0 , φ(B , V
of (A∞ , B∞
∞
xx,∞ )). This is because


0
P (A∞ , B∞
, φ(B∞ , Vxx,∞ )) ∈ ∂(C × I)


0
0
≤ P (A∞ , B∞
, φ(B∞ , Vxx,∞ )) ∈ ∂C × I + P (A∞ , B∞
, φ(B∞ , Vxx,∞ )) ∈ C × ∂I

0
+P (A∞ , B∞
, φ(B∞ , Vxx,∞ )) ∈ ∂C × ∂I

0
= P (A∞ , B∞
) ∈ ∂C, φ(B∞ , Vxx,∞ ) ∈ I

0
= P (A∞ , B∞
) ∈ ∂C | φ(B∞ , Vxx,∞ ) ∈ I · P {φ(B∞ , Vxx,∞ ) ∈ I} = 0,
0 ) | φ(B , V
where the last equality follows from the fact that C is a continuity set of (A∞ , B∞
∞
xx,∞ ) =

1. Thus, according to (A2) and Portmanteau’s Theorem, as n → ∞,
 √

√
√ 0
√
0
n(τ̂Y − τY , τ̂X
) ∈ C, φ( nτ̂X , Vxx ) = 1 = P
n(τ̂Y − τY ), nτ̂X
, φ( nτ̂X , Vxx ) ∈ C × I


0
0
→ P (A∞ , B∞
, φ(B∞ , Vxx,∞ )) ∈ C × I = P (A∞ , B∞
) ∈ C, φ(B∞ , Vxx,∞ ) = 1 .
P

√

0 ) | φ(B , V
Hence for any continuity set C of (A∞ , B∞
∞
xx,∞ ) = 1, as n → ∞,

√
0 ) ∈ C, φ(√nτ̂ , V ) = 1}
√
P { n(τ̂Y − τY , τ̂X
xx
X
√
P
n(τ̂Y −
∈ C | φ( nτ̂X , Vxx ) = 1 =
P {φ( nτ̂X , Vxx ) = 1}
0 ) ∈ C, φ(B , V

P {(A∞ , B∞
∞
xx,∞ ) = 1}
0
→
= P (A∞ , B∞
) ∈ C | φ(B∞ , Vxx,∞ ) = 1 .
P {φ(B∞ , Vxx,∞ ) = 1}
√

0
τY , τ̂X
)

4

Therefore, Proposition A1 holds.
Proposition A1 implies the following corollary, including Proposition 2 as a special case.
Corollary A1. Under ReG,



√
 n(τ̂Y − τY )
 √

nτ̂X

√

nτ̂X

 
. A
∈G ∼  
B

B ∈ G,

where (A, B 0 ) ∼ N (0, V ).
d

0 ) ∼ N (0, V ). As n → ∞, V → V , and then (A, B 0 , V ) −→
Proof of Corollary A1. Let (A∞ , B∞
∞
∞
xx
0 ,V
(A∞ , B∞
xx,∞ ) . The same logic as the proof of Proposition A1 implies

 
A
 
B





d
 A∞ 
φ(B, Vxx ) = 1 −→ 

B∞

φ(B∞ , Vxx,∞ ) = 1

Therefore, according to Proposition A1, Corollary A1 holds.
The following corollary shows the asymptotic distribution of τ̂Y under ReG.
Corollary A2. Under ReG,
√

n(τ̂Y − τY ) |

√

.

−1
nτ̂X ∈ G ∼ ε + Vτ x Vxx
B | B ∈ G,

where ε ∼ N (0, Vτ τ (1 − R2 )) is independent of B ∼ N (0, Vxx ).
Proof of Corollary A2. The residual from the linear projection of A on B is
−1
ε = A − Vτ x Vxx
B ∼ N (0, (1 − R2 )Vτ τ ),

which is independent of B. According to Corollary A1,
√

n(τ̂Y − τY ) |

√

nτ̂X ∈ G

.

∼ A|B∈G

5

−1
∼ ε + Vτ x Vxx
B | B ∈ G.

The following corollary shows the asymptotic unbiasedness of τ̂Y and balance in means of all
covariates, which includes Corollaries 1 and 3 as special cases.
√
√
Corollary A3. Under ReG, Ea { n(τ̂Y − τY ) | nτ̂X ∈ G} = 0.
√
√
Proof of Corollary A3. According to Proposition A1, Ea { n(τ̂Y − τY ) | nτ̂X ∈ G} = E(A∞ |
0 ) ∼ N (0, V ). Because φ satisfies φ(µ, V
B∞ ∈ G∞ ), where (A∞ , B∞
∞
xx,∞ ) = φ(−µ, Vxx,∞ ), we
0 ) ∼ (−A , −B 0 ), we have
know that B∞ ∈ G∞ if and only if −B∞ ∈ G∞ . Using (A∞ , B∞
∞
∞

E(A∞ | B∞ ∈ G∞ ) = E(−A∞ | −B∞ ∈ G∞ ) = E(−A∞ | B∞ ∈ G∞ ) = −E(A∞ | B∞ ∈ G∞ ).
√
√
Thus, E(A∞ | B∞ ∈ G∞ ) = 0, and Ea { n(τ̂Y − τY ) | nτ̂X ∈ G} = 0. Because covariates are
outcomes unaffected by treatment, the difference-in-means of any covariate has asymptotic mean
0.
Proof of Proposition 1. We have
−1
Vτ x Vxx
Vxτ



2 −1
(r1 r0 )−1 SX
r1−1 SX,Y (1) + r0−1 SX,Y (0)

r0 2
r1
2 −1
=
SY (1)|X + SY2 (0)|X + 2SY (1),X SX
SX,Y (0)
r1
r0
n
o

2 −1
SX,Y (0)
= r1−1 SY2 (1)|X + r0−1 SY2 (0)|X − SY2 (1)|X + SY2 (0)|X − 2SY (1),X SX
=

r1−1 SY (1),X + r0−1 SY (0),X

= r1−1 SY2 (1)|X + r0−1 SY2 (0)|X − Sτ2|X .

The sampling squared multiple correlation between τ̂Y and τ̂X under CRE has the following equivalent forms:
−1 V
r1−1 SY2 (1)|X + r0−1 SY2 (0)|X − Sτ2|X
Vτ x Vxx
xτ
=
Corr(τ̂Y , τ̂X ) =
Vτ τ
r1−1 SY2 (1) + r0−1 SY2 (0) − Sτ2

=

SY2 (1)
r1 Vτ τ

2

R (1) +

SY2 (0)
r0 Vτ τ

Therefore, Proposition 1 holds.

6

R2 (0) −

Sτ2 2
R (τ ) = R2 .
Vτ τ

A3.

Improvements Under Rerandomization

A3.1.

Reductions in asymptotic variances

First we investigate the reduction in asymptotic sampling variances under ReM and ReMT, and
2 as the limit of R2 , and ρ2
2
then we consider ReG. We introduce R∞
t,∞ as the limit of ρt (1 ≤ t ≤ T ).
2 and ρ2
The existences of R∞
t,∞ are guaranteed by the convergence of V .

Proof of Corollary 2. Recall that B∞ ∼ N (0, Vxx,∞ ). According to Proposition A1 and the
results for Gaussian covariates (Morgan and Rubin 2012, Theorem 3.1), the asymptotic sampling
√
variance of nτ̂X is

Vara


√
√
nτ̂X | nτ̂X ∈ M = Var(B∞ | B∞ ∈ M∞ ) = vK,a Var(B∞ ) = vK,a Vxx,∞ = lim vK,a Vxx .
n→∞

√
Because Vara ( nτ̂X ) = Var(B∞ ) = Vxx,∞ , we can deduce the PRIASV of τ̂X .
According to Theorem 1, under ReM, the asymptotic sampling variance of τ̂Y is
√



2
2
2
2
nτ̂X ∈ G = Vτ τ,∞ (1 − R∞
)Var(ε0 ) + R∞
Var(LK,a ) = Vτ τ,∞ 1 − R∞
+ R∞
vK,a


2
=Vτ τ,∞ 1 − (1 − vK,a )R∞
= lim Vτ τ 1 − (1 − vK,a )R2 .
Vara

√

n(τ̂Y − τY ) |

n→∞

√
Because Vara { n(τ̂Y − τY )} = Vτ τ,∞ , the PRIASV of τ̂Y is
√
√
Vara { n(τ̂Y − τY ) | nτ̂X ∈ G}
2
√
1−
= (1 − vK,a )R∞
= lim (1 − vK,a )R2 .
n→∞
Vara { n(τ̂Y − τY )}

Proof of Corollary 4. We first derive the asymptotic sampling variance and PRIASV of τ̂Y , and
then derive those of τ̂X . According to Theorem 3, for ReMT, the asymptotic sampling variance of
τ̂Y is

Vara

√

n(τ̂Y − τY ) |

√

(
nτ̂X ∈ T

= Vτ τ,∞

ρ2T +1,∞ Var(ε0 )

+

T
X
t=1

7

)
ρ2t,∞ Var(Lkt ,at )

T
X

(
ρ2T +1,∞

= Vτ τ,∞

+

)
ρ2t,∞ vkt ,at

t=1

(
1−

= Vτ τ,∞

T
X

(1 − vkt ,at )ρ2t,∞

= lim Vτ τ
n→∞

t=1

where the last line follows from

PT +1
t=1

T
X
1−
(1 − vkt ,at )ρ2t

(

)

)
,

t=1

ρ2t,∞ = 1. The PRIASV of τ̂Y is

√
√
T
T
X
Vara { n(τ̂Y − τY ) | nτ̂X ∈ T } X
2
√
(1 − vkt ,at )ρ2t .
=
(1 − vkt ,at )ρt,∞ = lim
1−
n→∞
Vara { n(τ̂Y − τY )}
t=1
t=1
2
Let X[tj ] be the jth covariate in tier t, and Rl,t
be the finite population squared multiple
j

2
correlation between Xi [tj ] and Xi [l] for 1 ≤ l ≤ T , with R0,t
= 0. The PRIASV for the outcome
j

implies that the PRIASV of τ̂X[tj ] is
T
X

lim

n→∞



2
2
− Rl−1,t
(1 − vkl ,al ) Rl,t
j



j

l=1

(
= lim

n→∞

RT2 ,t −

T
X

j



2
2
vkl ,al Rl,t
− Rl−1,t
j

j



)
.

l=1

2 = 1 for l ≥ t, we can further simplify the PRIASV of τ̂
Because Rl,t
X[tj ] as
j

(
lim

n→∞

1−

t
X



2
2
vkl ,al Rl,t
− Rl−1,t
j


j

)
.

l=1

To derive the asymptotic sampling variance of τ̂X , we use the notation introduced in the proof
of Theorem 3. Let Γ∞ be the limit of the linear transformation matrix Γ, and G∞ = Γ∞ B∞ =
(G01,∞ , G02,∞ , . . . , G0T,∞ )0 be the block-wise Gram–Schmidt orthogonalization of B∞ , where Gt,∞ is
a kt dimensional random vector. According to Proposition A1 and the fact that (G1,∞ , G2,∞ , . . . , GT,∞ )
are mutually independent, the asymptotic sampling variance of τ̂X is

Vara

√

nτ̂X |

√



0
−1
nτ̂X ∈ T = Var (B∞ | B∞ ∈ T∞ ) = Var Γ−1
∞ G∞ | Gt,∞ Var(Gt,∞ ) Gt,∞ ≤ at , 1 ≤ t ≤ T

0 −1
= Γ−1
∞ diag {vk1 ,a1 Var(G1,∞ ), . . . , vkT ,aT Var(GT,∞ )} Γ∞
n
o
−1
2
2
= lim Γ−1 diag vk1 ,a1 (r1 r0 )−1 SE[1]
, . . . , vkT ,aT (r1 r0 )−1 SE[T
Γ0
]
n→∞


−1
n2 −1
2
2
= lim
Γ diag vk1 ,a1 SE[1]
, . . . , vkT ,aT SE[T
Γ0
.
]
n→∞ n1 n0

8

According to Proposition A1, for ReG,

Vara
Vara

√

√

nτ̂X |

n(τ̂Y − τY ) |

√


nτ̂X ∈ T = Var(B∞ | B∞ ∈ T ) ≡ Vxx,φ,∞ = lim Vxx,φ ,
n→∞

√

nτ̂X ∈ T

= Var(A∞ | B∞ ∈ T∞ )
−1
−1
= Var A∞ − Vτ x,∞ Vxx,∞
B∞ + Vτ x,∞ Vxx,∞
B∞ | B∞ ∈ T∞



−1
−1
= Var(A∞ − Vτ x,∞ Vxx,∞
B∞ ) + Var Vτ x,∞ Vxx,∞
B∞ | B∞ ∈ T∞
2
−1
−1
= Vτ τ,∞ (1 − R∞
) + Vτ x,∞ Vxx,∞
Vxx,φ,∞ Vxx,∞
Vxτ,∞

−1
−1
= lim Vτ τ (1 − R2 ) + Vτ x Vxx
Vxx,φ Vxx
Vxτ .
n→∞

We can then immediately check whether ReM reduces the sampling covariance matrix of the
difference-in-means of the covariates.

A3.2.

Reductions in quantile ranges in ReM and ReMT

To prove Theorem 2, we need the following two lemmas.
Lemma A3. Let ε0 , η ∼ N (0, 1) be independent. For any a > 0 and c ≥ 0,

P


p
1 − ρ2 · ε0 + ρη ≥ c | η 2 ≤ a

is a decreasing function of ρ for ρ ∈ [0, 1].
Proof of Lemma A3. For any a > 0, let F (·) and f (·) denote the cumulative distribution and probability density of N (0, 1), and let G(·) and g(·) denote the cumulative distribution and probability
density of η | η 2 ≤ a. We have
 Z
p
P
1 − ρ2 · ε0 + ρη ≥ c | η 2 ≤ a =

∞

η≥

P

−∞

Z

∞

(
1−G

=
−∞

!
1 − ρ2 · x 2
| η ≤ a dF (x)
ρ
!)
p
c − 1 − ρ2 · x
dF (x).
ρ

c−

p

Taking the partial derivative with respect to ρ, we have
 Z ∞
∂ p
P
1 − ρ2 · ε0 + ρη ≥ c | η 2 ≤ a =
−g
∂ρ
−∞
9

! x√ 1 − c
p
1 − ρ2 · x − c
1−ρ2
dF (x)
ρ
ρ2



! Z
∞
t + cρ
t + cρ
tρ + c
p
−g (t)
−g (t)
=
=
dF
f
2
ρ(1 − ρ )
ρ(1 − ρ2 )
1 − ρ2
−∞
−∞
!
Z ∞

tρ + c
2 −3/2
g(t)f p
= − 1−ρ
(t + cρ)dt.
1 − ρ2
−∞
Z

∞

tρ + c
p
1 − ρ2

!

ρ
p
dt
1 − ρ2

The integral part in the above formula is
!
!
Z √a
tρ + c
1
tρ + c
g(t)f p
(t + cρ)dt =
f (t)f p
(t + cρ)dt
P (η 2 ≤ a) −√a
1 − ρ2
1 − ρ2
−∞
 2



Z √a
Z √a
2
t
(t + cρ)2
1
e−c /2
(tρ + c)2
exp − −
exp −
(t + cρ)dt =
(t + cρ)dt
2πP (η 2 ≤ a) −√a
2
2(1 − ρ2 )
2πP (η 2 ≤ a) −√a
2(1 − ρ2 )


Z √a+cρ
2
e−c /2
u2
exp −
udu ≥ 0.
2πP (η 2 ≤ a) −√a+cρ
2(1 − ρ2 )

Z

=
=

∞

Therefore, ∂P

p

1 − ρ2 · ε0 + ρη ≥ c | η 2 ≤ a /∂ρ ≤ 0.

Lemma A4. Let ε0 ∼ N (0, 1), LK,a ∼ D1 | D 0 D ≤ a, where D = (D1 , . . . , DK )0 ∼ N (0, IK ), and
(ε0 , LK,a ) are mutually independent. Then, for any a > 0 and c ≥ 0,

P

p


1 − ρ2 · ε0 + ρLK,a ≥ c

is a decreasing function of ρ for ρ ∈ [0, 1].

Proof of Lemma A4. The independence of ε0 and D implies

P


p

p
1 − ρ2 · ε0 + LK,a ≥ c = P
1 − ρ2 · ε0 + ρD1 ≥ c | D 0 D ≤ a .

Assume 0 ≤ ρ1 ≤ ρ2 ≤ 1, and (d2 , . . . , dK ) satisfies

PK

2
k=2 dk

< a. Conditioning on (D2 , . . . , DK ) =

(d2 , . . . , dK ), Lemma A3 implies

P

!
K
q
X
2
2
2
1 − ρ1 · ε0 + ρ1 D1 ≥ c | D1 ≤ a −
Dk , D2 = d2 , . . . , DK = dK
k=2

≥P

!
K
q
X
1 − ρ22 · ε0 + ρ2 D1 ≥ c | D12 ≤ a −
Dk2 , D2 = d2 , . . . , DK = dK .
k=2

10

Taking expection for both sides, we have

q

q
0
0
2
2
P
1 − ρ1 · ε0 + ρ1 D1 ≥ c | D D ≤ a ≥ P
1 − ρ2 · ε0 + ρ2 D1 ≥ c | D D ≤ a .

Therefore, Lemma A4 holds.

Proof of Theorem 2. According to Theorem 1, the lengths of (1 − α) quantile ranges of the
p
p
√
2 ) V
asymptotic distributions of n(τ̂Y −τY ) under ReM and CRE are 2ν1−α/2 (R∞
τ τ,∞ and 2z1−α/2 Vτ τ,∞ ,
2 ) and Lemma A4, we know that ν
2
respectively. According to the definition of ν1−α/2 (R∞
1−α/2 (R∞ )
2 .
is a decreasing function of R∞

To prove Theorem 4, we need the following four lemmas. We first define a random variable to
be SUM if it is symmetric and unimodal around zero.
Lemma A5. Let ζ0 , ζ1 and ζ2 be three jointly independent random variables. If
(1) ζ0 is continuous and SUM, or ζ0 = 0;
(2) ζ1 and ζ2 are symmetric around 0;
(3) P (ζ1 ≥ c) ≤ P (ζ2 ≥ c) for any c > 0,
then P (ζ0 + ζ1 ≥ c) ≤ P (ζ0 + ζ2 ≥ c) for any c > 0.

Proof for Lemma A5. Note that when ζ0 = 0, Lemma A5 holds automatically. We consider only
the case where ζ0 is continuous and SUM. Let Fζ0 (·) be the cumulative distribution function of ζ0 .
For any c > 0,
Z

∞

P (ζ0 + ζ1 ≥ c) =
P (ζ1 ≥ c − x)dFζ0 (x)
−∞
Z c
Z ∞
=
P (ζ1 ≥ c − x)dFζ0 (x) +
P (ζ1 ≥ c − x)dFζ0 (x)
−∞
0

Z

Z

c
∞

P (ζ1 ≥ t)dFζ0 (c − t) +
P (ζ1 ≥ −t)dFζ0 (c + t)
0
Z∞∞
Z ∞
=
P (ζ1 ≥ t)d {−Fζ0 (c − t)} +
P (ζ1 ≤ t)dFζ0 (c + t)
0
0
Z ∞
Z ∞
=
P (ζ1 ≥ t)d {Fζ0 (t − c) − 1} +
{1 − P (ζ1 ≥ t)} dFζ0 (c + t)
=

0

0

11

Z

∞

=

Z

∞

P (ζ1 ≥ t)d {Fζ0 (t − c) − Fζ0 (c + t)}

dFζ0 (c + t) +
0

Z

0

∞

Z

∞

P (ζ1 ≥ t)d {−P (t − c ≤ ζ0 ≤ t + c)} .

dFζ0 (c + t) +

=

0

0

Similarly,
Z

∞

P (ζ0 + ζ2 ≥ c) =

∞

Z

P (ζ2 ≥ t)d {−P (t − c ≤ ζ0 ≤ t + c)} .

dFζ0 (c + t) +
0

0

Because ζ0 is SUM and continuous, −P (t − c ≤ ζ0 ≤ t + c) is a continuous increasing function of t
when t ≥ 0. Because P (ζ1 ≥ t) ≤ P (ζ2 ≥ t) for any t > 0, we have that for all c > 0,
∞

Z
P (ζ0 + ζ1 ≥ c) =
Z0
Z ∞
≤
dFζ0 (c + t) +
0

Z

∞

P (ζ1 ≥ t)d {−P (t − c ≤ ζ0 ≤ t + c)}

dFζ0 (c + t) +
0

∞

P (ζ2 ≥ t)d {−P (t − c ≤ ζ0 ≤ t + c)} = P (ζ0 + ζ2 ≥ c).

0

Lemma A6. [Wintner 1936] If ζ1 and ζ2 are SUM and independent, then ζ1 + ζ2 is also SUM.
+1
+1
be two
and {ρ̃t }Tt=1
Lemma A7. Let ε0 , η1 , η2 , . . . , ηT be (T + 1) i.i.d N (0, 1). Let {ρt }Tt=1
PT +1 2 PT +1 2
nonnegative constant sequences satisfying t=1 ρt = t=1 ρ̃t = 1. If there exists 1 ≤ t0 ≤ T such

that
ρt0 ≥ ρ̃t0 ,

ρt = ρ̃t (t 6= t0 , 1 ≤ t ≤ T ),

ρT +1 ≤ ρ̃T +1 ,

then for any c ≥ 0 and at > 0 (1 ≤ t ≤ T ),

P

ρT +1 ε0 +

T
X

!
ρt ηt ≥ c |

ηt2

≤ at , 1 ≤ t ≤ T

≤P

ρ̃T +1 ε0 +

t=1

T
X

!
ρ̃t ηt ≥ c |

ηt2

≤ at , 1 ≤ t ≤ T

.

t=1

Proof of Lemma A7. Without loss of generality, we assume t0 = 1. Then ρ1 ≥ ρ̃1 and ρ21 + ρ2T +1 =
ρ̃21 + ρ̃2T +1 . According to Lemma A3, for any c ≥ 0,








ρT +1
ρ1
ρ̃T +1
ρ̃1
P q
ε0 + q
η1 ≥ c | η12 ≤ a1  ≤ P  q
ε0 + q
ρ̃1 ≥ c | η12 ≤ a1 
2
2
2
2
2
2
2
2
ρ1 + ρT +1
ρ1 + ρT +1
ρ̃1 + ρ̃T +1
ρ̃1 + ρT +1

12

which implies that, for any c ≥ 0,



P ρT +1 ε0 + ρ1 η1 ≥ c | η12 ≤ a1 ≤ P ρ̃T +1 ε0 + ρ̃1 η1 ≥ c | η12 ≤ a1 .

According to Proposition 3 and Lemma A6,

PT

t=2 ρt


ηt | ηt2 ≤ at is SUM. Lemma A5 implies that

for any c > 0,

P

ρT +1 ε0 +

T
X

!
ρt ηt ≥ c | ηt2 ≤ at , 1 ≤ t ≤ T

≤P

ρ̃T +1 ε0 +

t=1

T
X

!
ρ̃t ηt ≥ c | ηt2 ≤ at , 1 ≤ t ≤ T

.

t=1

Therefore, Lemma A7 holds.
Lemma A8. Let ε0 ∼ N (0, 1), Lkt ,at ∼ Dt1 | Dt0 Dt ≤ at , where Dt = (Dt1 , . . . , Dtkt ) ∼ N (0, Ikt ),
+1
+1
and (ε0 , Lk1 ,a1 , Lk2 ,a2 , . . . , LkT ,aT ) are mutually independent. Let {ρt }Tt=1
and {ρ̃t }Tt=1
be two
PT +1 2 PT +1 2
nonnegative constant sequences satisfying t=1 ρt = t=1 ρ̃t = 1. If there exists 1 ≤ t0 ≤ T such

that
ρt0 ≥ ρ̃t0 ,

ρt = ρ̃t (t 6= t0 , 1 ≤ t ≤ T ),

ρT +1 ≤ ρ̃T +1 ,

then for any c ≥ 0 and at > 0 (1 ≤ t ≤ T ),

P

ρT +1 ε0 +

T
X

!
ρt Lkt ,at ≥ c

≤P

ρ̃T +1 ε0 +

T
X

!
ρ̃t Lkt ,at ≥ c .

t=1

t=1

Proof of Lemma A8. Let Dt = (Dt1 , . . . , Dtkt )0 ∼ N (0, Ikt ), 1 ≤ t ≤ T , and (D1 , . . . , DT , ε0 ) be
mutually independent. Conditioning on Dtj = dtj (1 ≤ t ≤ T, j ≥ 2), according to Lemma A7,

P

≤ P

ρT +1 ε0 +

ρ̃T +1 ε0 +

T
X

ρt Dt1 ≥ c |

2
Dt1

≤ at −

kt
X

t=1

i=2

T
X

kt
X

2
ρ̃t Dt1 ≥ c | Dt1
≤ at −

t=1

!
2
Dti
, Dtj

= dtj , 1 ≤ t ≤ T, j ≥ 2
!

2
Dti
, Dtj = dtj , 1 ≤ t ≤ T, j ≥ 2 .

i=2

Taking expectations of both sides, we have

P

ρT +1 ε0 +

T
X

!
ρt Dt1 ≥ c |

t=1

13

Dt0 Dt

≤ at , 1 ≤ t ≤ T

≤ P

ρ̃T +1 ε0 +

T
X

!
ρ̃t Dt1 ≥ c |

Dt0 Dt

≤ at , 1 ≤ t ≤ T

.

t=1

Therefore, Lemma A8 holds.

Proof of Theorem 4. According to Theorem 3, the lengths of (1 − α) quantile ranges of the
p
√
asymptotic distributions of n(τ̂Y −τY ) under ReMT and CRE are 2ν1−α/2 (ρ21,∞ , ρ22,∞ , . . . , ρ2T,∞ ) Vτ τ,∞
p
and 2z1−α/2 Vτ τ,∞ , respectively. According to the definition of ν1−α/2 (ρ21,∞ , ρ22,∞ , . . . , ρ2T,∞ ) and
Lemma A8, ν1−α/2 (ρ21,∞ , ρ22,∞ , . . . , ρ2T,∞ ) is a decreasing function of ρ2t,∞ , for any 1 ≤ t ≤ T .
Therefore, Theorem 4 holds.

A4.
A4.1.

Conservativeness in Inference
Conservativeness of the sampling variance estimators

The following lemma, which does not require more moment conditions than Condition 1, is useful for obtaining asymptotically conservative estimators for the sampling variances and sampling
distributions.
Lemma A9. Under ReG, if Condition 1 holds, then for any (Ai , Bi ) equal to (Yi (1), Yi (1)), (Yi (0), Yi (0)),
(Yi (1), Xki ), (Yi (0), Xki ) or (Xki , Xli ), we have

sAB (z) − SAB = op (1),

(z = 0, 1)

sAB (z) is the sample covariance between the Ai ’s and Bi ’s under treatment arm z, and SAB is the
finite population covariance between the Ai ’s and Bi ’s.
Proof of Lemma A9. The key is to bound the variance of sAB (z) under ReG. According to the law
of total expectation,
h
i
h
i

√
2 √
2
nτ̂X ∈ G · E {sAB (z) − SAB } | nτ̂X ∈ G
E {sAB (z) − SAB } = P
h
i

√
√
+P
nτ̂X ∈
/ G · E {sAB (z) − SAB }2 | nτ̂X ∈
/G
h
i

√
√
≥P
nτ̂X ∈ G · E {sAB (z) − SAB }2 | nτ̂X ∈ G .
14

Therefore,
h
i
i
−1 h
√
√
E {sAB (z) − SAB }2 | nτ̂X ∈ G ≤ P
nτ̂X ∈ G
E {sAB (z) − SAB }2
=P

√

nτ̂X ∈ G

−1

Var {sAB (z)} ,

(A3)

where (A3) follows from the fact that sAB (z) is unbiased for SAB under CRE (Cochran 1977)
because, under CRE, units receiving treatment arm z is a simple random sample of size nz .
P
P
obs (z) = n−1
Let Āobs (z) = n−1
z
z
i:Zi =z Ai and B̄
i:Zi =z Bi be the averages of Ai ’s and Bi ’s
under treatment arm z. We first decompose the sample covariance sAB (z) as








X
1
nz
(Ai − Ā)(Bi − B̄) − Āobs (z) − Ā B̄ obs (z) − B̄
,
sAB (z) =

nz − 1  nz
i:Zi =z

and then obtain an upper bound of its variance under CRE


1 X



Var{sAB (z)} =
Var
(Ai − Ā)(Bi − B̄) − Āobs (z) − Ā B̄ obs (z) − B̄

 nz
(nz − 1)2
i:Zi =z






n


o
2
X
2nz
Var 1
≤
(Ai − Ā)(Bi − B̄) + Var Āobs (z) − Ā B̄ obs (z) − B̄  ,
 nz

(nz − 1)2
n2z

i:Zi =z

(A4)
which follows from the Cauchy-Schwarz inequality. Below we consider the two terms in (A4)
separately. The first term in (A4) is bounded by:


1 X

Var
(Ai − Ā)(Bi − B̄)
 nz

i:Zi =z

2


n 
n

X
X
1
1
1
1
(Ai − Ā)(Bi − B̄) −
=
−
(Aj − Ā)(Bj − B̄)


nz
n n−1
n
i=1

≤

1 1
nz n − 1

n
X

j=1

(Ai − Ā)2 (Bi − B̄)2

i=1
n

≤

n 1
1 X
· max (Aj − Ā)2 ·
(Bi − B̄)2 .
nz n 1≤j≤n
n−1
i=1

15

(A5)

Because (Āobs (z) − Ā)2 ≤ max1≤j≤n (Aj − Ā)2 , the second term in (A4) is bounded by:
n

o
Var Āobs (z) − Ā B̄ obs (z) − B̄

2 
2 
obs
obs
B̄ (z) − B̄
≤ E
Ā (z) − Ā
h
i
≤ max (Aj − Ā)2 · E {B̄ obs (z) − B̄}2
1≤j≤n
n
o
= max (Aj − Ā)2 · Var B̄ obs (z)
1≤j≤n

=
≤

2

max (Aj − Ā) ·

1≤j≤n



1
1
−
nz
n



n

1 X
(Bi − B̄)2
n−1

1
n 1
· max (Aj − Ā)2 ·
nz n 1≤j≤n
n−1

i=1

n
X

(Bi − B̄)2 .

(A6)

i=1

Therefore, according to (A4)–(A6), we can bound Var{sAB (z)} by:
n

Var{sAB (z)} ≤

n 1
1 X
4n2z
2
·
·
max
(A
−
Ā)
·
(Bi − B̄)2 ,
j
(nz − 1)2 nz n 1≤j≤n
n−1
i=1

√
which converges to zero under Condition 1. Recall that P ( nτ̂X ∈ G) converges to the asymptotic
√
acceptance probability pa > 0. Then, according to (A3), E[{sAB (z) − SAB }2 | nτ̂X ∈ G] = o(1).
By the Markov inequality, under ReG and Condition 1,

sAB (z) − SAB = Op

!
r h
i
2 √
= op (1).
E {sAB (z) − SAB } | nτ̂X ∈ G

Therefore, Lemma A9 holds.
Lemma A10. Under ReG, if Condition 1 holds, then

2
s2Y (z) − SY2 (z) = op (1), sY (z),X − SY (z),X = op (1), s2X,z − SX
= op (1),

(z = 0, 1)

where s2X,z is the sample variance of X under treatment arm z.
Proof of Lemma A10. Lemma A10 follows directly from Lemma A9.
Let Ṽτ τ = Vτ τ + Sτ2 − Sτ2|X ≥ Vτ τ and Ṽτ τ,∞ be the limit of Ṽτ τ . Under ReM, according to

16

Lemma A10, V̂τ τ is asymptotically unbiased for Ṽτ τ , and
V̂τ τ R̂2 = r1−1 s2Y (1)|X + r0−1 s2Y (0)|X − sY (1),X − sY (0),X



2
SX

−1

sX,Y (1) − sX,Y (0)



is asymptotically unbiased for Vτ τ R2 . Thus, the sampling variance estimator is asymptotically
unbiased for

Ṽτ τ − (1 − vK,a )Vτ τ R2 ≥ Vτ τ − (1 − vK,a )Vτ τ R2 → Vara

√

n(τ̂Y − τY ) |

√

nτ̂X ∈ M .

Under ReMT, according to Lemma A10, V̂τ τ is asymptotically unbiased for Ṽτ τ , and V̂τ τ ρ̂2t is
asymptotically unbiased for Vτ τ ρ2t . Thus, the sampling variance estimator is asymptotically unbiased for

Ṽτ τ −

T
X

(1 − vkt ,at )Vτ τ ρ2t ≥ Vτ τ −

t=1

T
X

(1 − vkt ,at )Vτ τ ρ2t → Vara

√

n(τ̂Y − τY ) |

√

nτ̂X ∈ T

.

t=1

Under ReG, according to Lemma A10, V̂τ τ is asymptotically unbiased for Ṽτ τ , V̂τ τ R̂2 is asymptotically unbiased for Vτ τ R2 , and V̂τ x is asymptotically unbiased for Vτ x . Therefore, the sampling
variance estimator is asymptotically unbiased for
−1
−1
Ṽτ τ − Vτ τ R2 + Vτ x Vxx
Vxx,φ Vxx
Vxτ
−1
−1
≥Vτ τ (1 − R2 ) + Vτ x Vxx
Vxx,φ Vxx
Vxτ → Vara

√

n(τ̂Y − τY ) |

√

nτ̂X ∈ G .

Above all, the sampling variance estimators are asymptotically conservative.

A4.2.

Conservativeness of confidence interval

First, we consider ReM. According to Lemma A10,
q

q
V̂τ τ

1−

R̂2


q
q
p
d
2
2
2 ·L
· ε0 + R̂ · LK,a −→ Ṽτ τ,∞ − Vτ τ,∞ R∞ · ε0 + Vτ τ,∞ R∞
K,a .

p
Thus ν1−α/2 (R̂2 ) V̂τ τ is consistent for the (1−α/2)th quantile of the distribution on the right hand
p
2 ) V
side of the above formula, which is larger than or equal to ν1−α/2 (R∞
τ τ,∞ due to Proposition
17

3 and Lemma A5.
Second, we consider ReMT. According to Lemma A10,
q
V̂τ τ

ρ̂T +1 ε0 +

T
X

!
ρ̂t Lkt ,at

d

−→

T q
q
X
Vτ τ,∞ ρ2t,∞ · Lkt ,at .
Ṽτ τ,∞ − Vτ τ,∞ + Vτ τ,∞ ρ2T +1,∞ · ε0 +

t=1

t=1

p
Thus ν1−α/2 (ρ̂21 , . . . , ρ̂2T ) V̂τ τ is consistent for the (1 − α/2)th quantile of the distribution on the
p
right hand side of above formula, which is larger than or equal to ν1−α/2 (ρ21,∞ , . . . , ρ2T,∞ ) Vτ τ,∞
due to Proposition 3, and Lemmas A5 and A6.
Finally, we consider ReG, and construct the confidence interval. Let V̂ε = V̂τ τ (1 − R̂2 ) be the
√
−1 B | B ∈ G, where
variance estimator for ε in (12). Let qξ (λ) be the ξth quantile of λε0 + V̂τ x Vxx
ε0 ∼ N (0, 1) is independent of B ∼ N (0, Vxx ). For any ξ ≥ 0.5, let q̂ξ = max0≤λ≤V̂ε qξ (λ). The

√
√ 
final confidence interval for τY is then τ̂Y − q̂1−α/2 / n, τ̂Y + q̂1−α/2 / n . According to Lemma
A10, for any λ ≥ 0,
√

d

−1
λε0 + V̂τ x Vxx
B | B ∈ G −→

√

−1
λε0 + Vτ x,∞ Vxx,∞
B∞ | B∞ ∈ G∞ ,

(A7)

where B∞ ∼ N (0, Vxx,∞ ). Let ω1−α/2 (λ) be the (1 − α/2)th quantile of the distribution on the
right side of (A7). Then q1−α/2 (λ) is consistent for ω1−α/2 (λ). According to Lemma A10, V̂ is
2 ≥V
2
consistent for Ṽτ τ,∞ − Vτ τ,∞ R∞
τ τ,∞ (1 − R∞ ). Under some regularity conditions,

p

q̂1−α/2 = max q1−α/2 (λ) −→
0≤λ≤V̂ε

max
2
0≤λ≤Ṽτ τ,∞ −Vτ τ,∞ R∞


2
ω1−α/2 (λ) ≥ ω1−α/2 Vτ τ,∞ (1 − R∞
) .

−1 B | B ∈ G is unimodal, q̂
When Vτ x Vxx
1−α/2 = q1−α/2 (V̂ε ) according to Lemma A5.

18

