# Balanced and Robust Randomized Treatment Assignments: The Finite Selection
Model for the
Health Insurance Experiment and Beyond††thanks: We thank John Golden, Angela
Lee, and Bijan Niknam for helpful research assistance and comments. We also
thank participants at Euro-CIM 2023 for their valuable comments. This work was
supported through a grant from the Alfred P. Sloan Foundation (G-2020-13946).

Ambarish Chattopadhyay Stanford Data Science, Stanford University, 450 Jane
Stanford Way Wallenberg, Stanford, CA 94305; email<EMAIL_ADDRESS>Carl N. Morris Department of Statistics, Harvard University, 1 Oxford Street
Cambridge, MA 02138; email<EMAIL_ADDRESS>José R. Zubizarreta
Departments of Health Care Policy, Biostatistics, and Statistics, Harvard
University, 180 Longwood Avenue, Office 307-D, Boston, MA 02115; email:
<EMAIL_ADDRESS>
###### Abstract

The Finite Selection Model (FSM) was developed by Carl Morris in the 1970s for
the design of the RAND Health Insurance Experiment (HIE) (Morris 1979,
Newhouse et al. 1993), one of the largest and most comprehensive social
science experiments conducted in the U.S. The idea behind the FSM is that each
treatment group takes its turns selecting units in a fair and random order to
optimize a common assignment criterion. At each of its turns, a treatment
group selects the available unit that maximally improves the combined quality
of its resulting group of units in terms of the criterion. In the HIE and
beyond, we revisit, formalize, and extend the FSM as a general tool for
experimental design.

Leveraging the idea of D-optimality, we propose and analyze a new selection
criterion in the FSM. The FSM using the D-optimal selection function has no
tuning parameters, is affine invariant, and when appropriate, retrieves
several classical designs such as randomized block and matched-pair designs.
For multi-arm experiments, we propose algorithms to generate a fair and random
selection order of treatments. We demonstrate FSM’s performance in a case
study based on the HIE and in ten randomized studies from the health and
social sciences. On average, the FSM achieves 68% better covariate balance
than complete randomization and 56% better covariate balance than
rerandomization in a typical study. We recommend the FSM be considered in
experimental design for its conceptual simplicity, efficiency, and robustness.

Keywords: Causal inference; Covariate balance; Experimental design; Multi-
valued treatments

###### Contents

  1. 1 Introduction
    1. 1.1 The RAND Health Insurance Experiment
    2. 1.2 Toward balanced, efficient, and robust experimental designs
    3. 1.3 Contribution and outline
  2. 2 Design of the Health Insurance Experiment
  3. 3 Foundations and overview of the FSM
    1. 3.1 Setup and notation
    2. 3.2 Components of the FSM
  4. 4 The D-optimal selection function
  5. 5 Inference under the FSM
  6. 6 The Health Insurance Experiment
    1. 6.1 Data
    2. 6.2 Balance
    3. 6.3 Efficiency
    4. 6.4 Intuition and further explorations
  7. 7 Ten further studies in the health and social sciences
  8. 8 Practical considerations and extensions
    1. 8.1 Multi-group experiments
    2. 8.2 Stratified experiments
    3. 8.3 Sequential experiments
  9. 9 Summary and remarks
  10. Supplementary materials
    1. A Notation and estimands
    2. B Proofs of theoretical results
    3. C Properties of D-optimal selection function
    4. D Affine invariance and covariate balance
      1. D.1 Connection to A-optimality
    5. E Algorithms for constructing an SOM
      1. E.1 The SCOMARS algorithm
      2. E.2 SOM for multi-group experiments
    6. F FSM for stratified experiments
    7. G FSM for sequential experiments
    8. H A simulation study
      1. H.1 Setup
      2. H.2 Balance
      3. H.3 Efficiency
      4. H.4 Additional tables and figures from the simulation study
        1. H.5 Additional tables and figures from the Health Insurance Experiment
          1. H.6 Additional figures from the case studies

##  1 Introduction

###  1.1 The RAND Health Insurance Experiment

In the 1970’s, the challenge of financing and delivering high-quality and
affordable health care to all Americans was at the center of national policy
debate. At the time, two central questions were “How much more medical care
would people use if it is provided free of charge?” and “What are the
consequences of using more medical care on their health?” To address these and
other related questions, an interdisciplinary team of researchers led by
Joseph P. Newhouse at RAND designed and conducted the Health Insurance
Experiment (HIE), a large-scale, multi-year, randomized public policy
experiment developed and completed between 1971 and 1982. To this day, the HIE
is one of the largest and most comprehensive social science experiments ever
conducted in the U.S. Even now, four decades after its completion, evidence
from the HIE is still fundamental to the national discussion on health care
cost sharing and health care reform.

In the HIE, a representative sample of 2,750 families comprising more than
7,700 individuals was chosen from six urban and rural sites across the United
States. At the beginning of the study, participants completed a baseline
survey providing numerous demographic, medical, and socioeconomic
measurements. Families were then assigned to health insurance plans that
varied substantially in their coinsurance rates and out-of-pocket expenditure
maxima, for a total of 13 possible treatment groups. The goal of the study was
to estimate the marginal averages of utilization and health outcomes in each
of the six sites under each plan.

To provide the strongest possible evidence on health utilization and outcomes,
the study had to be randomized. However, achieving balance for numerous
continuous and categorical baseline covariates through randomization is
challenging in experiments with so many treatment groups and different
implementation sites. In the HIE the groups had to be balanced and
representative of the sites. In the health and social sciences, there is an
ever-increasing need for methods for random assignment of units into multiple
treatment groups that are balanced, efficient, and robust.

###  1.2 Toward balanced, efficient, and robust experimental designs

Randomized experiments are considered to be the gold standard for causal
inference, as randomization provides an unequivocal basis for inference and
control. In randomized experiments, the act of randomization ensures balance
on both observed and unobserved covariates on average. However, a given
realization of the random assignment mechanism may produce substantial
imbalances on one or more covariates. This imbalance problem can be
exacerbated in settings like the HIE, where treatments are multi-valued and
many baseline covariates exist, leading to loss in efficiency of the effect
estimates.

A variety of methods have been proposed in the literature to address this
problem, such as blocking (Fisher 1925, Fisher 1935, Cochran and Cox 1957),
optimal pair-matching (Greevy et al. 2004), greedy pair-switching (Krieger et
al. 2019), and designs using mixed-integer programming (Bertsimas et al.
2015). In particular, rerandomization (Morgan and Rubin 2012) has gained
popularity over the last few years and has become commonplace in experiments.
However, rerandomization may not protect against and be robust to chance
imbalances in functions of the covariates that are not explicitly addressed by
the rerandomization criterion (Banerjee et al. 2017), especially in
experiments with multi-valued ($>$2) treatments. Moreover, defining the
rerandomization criterion requires the selection of a tuning parameter
governing the acceptable degree of imbalance, which may be difficult to choose
and require iteration in practice.

To address these and other related challenges, we revisit and extend the
Finite Selection Model (FSM) for experimental design. The original version of
the FSM was proposed and developed by Carl N. Morris in the design of the HIE
(Morris 1979, Newhouse et al. 1993, Morris and Hill 2000). The idea behind the
FSM is that each treatment group takes turns in a fair and random order to
select units from a pool of available units such that, at each stage, each
treatment group selects the unit that maximally improves the combined quality
of its current group of units. The criterion for measuring quality is
flexible. Among other contributions, in this paper we develop a new criterion
based on D-optimality, which does not require tuning parameters.

To illustrate, Figure 1 exhibits the performance of complete randomization,
rerandomization, and the FSM in a version of the HIE data with four treatment
groups and 20 covariates. For rerandomization, we compute the maximum
Mahalanobis distance (across all pairs of treatment groups) based on the 20
covariates and their squares and pairwise products (i.e., all second-order
transformations), and following Lock (2011), accept 0.1% of the assignments
with the smallest covariate distance (see Sections 6.1 and 6.4 for details).
The figure displays the distribution of absolute standardized mean differences
(ASMD; Rosenbaum and Rubin 1985)111The absolute standardized mean difference
for a single covariate $X$ between treatment groups $g$ and $g^{\prime}$ is
$\text{ASMD}(X)={|\bar{X}_{g}-\bar{X}_{g^{\prime}}|}/{\sqrt{(s^{2}_{g}+s^{2}_{g^{\prime}})/{2}}}$,
where $\bar{X}_{g}$ and $s^{2}_{g}$ are the mean and variance of $X$ in
treatment group $g$, respectively. Please see Rosenbaum and Rubin (1985)) for
details. in covariates and the second-order transformations across multiple
realizations of the randomization mechanisms for the three designs. Lower
values of ASMD indicate better balance on the covariates or their
transformations. Better balance can improve the validity and credibility of a
study, and can also translate into increased efficiency and robustness.

Figure 1: Distributions of ASMD for complete randomization, rerandomization,
and the FSM, for 20 baseline covariates in the HIE data. Without tuning
parameters, the FSM handles multiple ($>$2) treatment groups and substantially
improves covariate balance and, thereby, statistical efficiency.

We observe that, as expected, rerandomization outperforms complete
randomization in terms of imbalances on the main covariates and the second-
order transformations. The FSM, however, markedly outperforms both methods for
both types of covariates without requiring tuning parameters. This analysis
reveals that, while rerandomization performs well by common covariate balance
standards (the majority of the ASMD is smaller than 0.1), there is room for
improvement. As we explain in Section 6, in experiments like the HIE, the
space of possible assignments is vast, and the FSM can meaningfully improve
the assignment of units into treatment groups to achieve better balance and
efficiency.

In a nutshell, the FSM does better because it progressively randomizes units
into treatment groups in a controlled manner towards a criterion that is
common to all groups and robust against general outcome models. As we show in
theory and in practice in sections 4, 6, and 7 the FSM is a flexible tool for
random assignment in various settings.

###  1.3 Contribution and outline

In this paper, we revisit, formalize, and extend the FSM for experimental
design. We show that the FSM can be used for balanced, efficient, and robust
random treatment assignments, outperforming common assignment methods on these
three dimensions. In particular, we re-introduce the FSM under the potential
outcomes framework (Neyman 1923, 1990, Rubin 1974). We use the sequentially
controlled Markovian random sampling (SCOMARS, Morris 1983) algorithm to
determine the selection order of treatments for two-group experiments and
extend it to multi-group experiments. We propose a new selection criterion for
treatments based on the idea of D-optimality and discuss its theoretical
properties. Under suitable conditions, we show that the FSM retrieves several
classical experimental designs, such as randomized block and matched-pair
designs. We explain model-based approaches to inference under the FSM and
develop randomization-based alternatives. We analyze the FSM’s performance
empirically and compare it to common assignment methods. Finally, we discuss
potential extensions of the FSM to more complex experimental design settings,
such as stratified experiments and experiments with sequential arrival of
units. In an accompanying paper (Chattopadhyay et al. 2021), we describe how
these methods can be implemented in the new FSM package for R, which is
publicly available on CRAN.

The paper proceeds as follows. In Section 2, we describe the design of the
RAND Health Insurance Experiment, focusing on the assignment of each family to
one of 13 health insurance plans. In Section 3, we present the setup,
notation, and main components of the FSM. In Section 4, we propose a selection
criterion based on D-optimality and analyze its properties. In Section 5, we
discuss inference under the FSM. In Section 6, we evaluate the performance of
the FSM and compare it to standard methods such as complete randomization and
rerandomization using the HIE data. In Section 7, we perform a similar
comparison using the data from ten experimental studies from the health and
social sciences. Finally, in Section 8 we consider extensions of the FSM to
other settings such as multi-group, stratified, and sequential experiments. In
Section 9, we conclude with a summary and remarks. In the Online Supplementary
Materials, we present all the proofs of the propositions and theorems,
extended theoretical results, further empirical results based on a simulation
study, and supplemental experimenal results on the HIE study and the ten case
studies.

##  2 Design of the Health Insurance Experiment

In the HIE, families were assigned to different health insurance plans using
the original version of the FSM. Initially, assignments were made in each of
the six HIE sites to 12 or 13 fee-for-service plans with varying combinations
of coinsurance (cost sharing) rates and income-related deductibles.
Coinsurance plans consisted of $0\%$ (free care), $25\%$, $50\%$, or $95\%$
coinsurance rates, plus a plan with mixed coinsurance rates, and an individual
deductible plan. Within the cost sharing plans, families were further assigned
to different out-of-pocket maxima where the out-of-pocket expenditures were
capped at 5%, 10%, or 15% of family income, with an annual maximum of $1,000
(Brook et al. 2006). To ensure that the resulting treatment groups were
balanced relative to the population of each site, the FSM considered a discard
group of study non-participants as an additional treatment group.

Listed in chronological order of study initiation, the following sites were
tracked for several years: Dayton, OH; Seattle, WA; Fitchburg, MA; Franklin
County, MA; Charleston, SC; and Georgetown County, SC. The FSM was used,
independently in each of the sites, to make random assignments to improve
balance on up to 22 family-level baseline covariates across treatment groups.
In each of the first two sites, the FSM was used multiple times for separate
independent subsets of families to maintain baseline data schedules. In
addition to estimating the overall marginal effects of health insurance plan
design on healthcare utilization and outcomes, the HIE team also sought to
understand how the experimental results were affected by particular design
choices, e.g., longer versus shorter enrollment duration, receiving versus not
receiving participation incentives, higher versus lower interviewing
frequency. To this end, four additional sub-experiments were conducted, and
the FSM was used to randomize families to the sub-treatment groups.

##  3 Foundations and overview of the FSM

###  3.1 Setup and notation

Consider a sample of $N$ units indexed by $i=1,...,N$. Each of these units is
to be assigned into one of $G$ treatment groups labeled by $g$, with
$g=1,...,G$. Write $n_{g}$ for the pre-specified size of group $g$. Denote
$Z_{i}\in\\{1,2,...,G\\}$ as the assigned treatment group label of unit $i$
and $\bm{Z}=(Z_{1},...,Z_{N})^{\top}$ as the vector of treatment group labels.
Following the potential outcomes framework for causal inference (Neyman 1923,
1990; Rubin 1974), each unit $i$ has a potential outcome under each treatment
$g$, $Y_{i}(g)$, but only one of these outcomes is observed:
$Y^{\text{obs}}_{i}=\sum_{g=1}^{G}\mathbbm{1}(Z_{i}=g)Y_{i}(g)$. Denote
$\bm{Y}(g)=(Y_{1}(g),...,Y_{N}(g))^{\top}$ as the vector of potential outcomes
under treatment $g$. Each unit has a vector of $K$ observed covariates,
$\boldsymbol{X}_{i}$. We write $(\underline{\bm{X}}_{\text{full}})_{N\times
k}$ for the matrix of observed covariates, and $\bar{\bm{X}}_{\text{full}}$
and $\underline{\bm{S}}_{\text{full}}$ for the mean vector and covariance
matrix of these covariates in the full sample, respectively. Denote
$(\underline{\tilde{\bm{X}}}_{\text{full}})_{N\times(k+1)}$ as the design
matrix in the full sample.222The design matrix includes a column of all 1’s
(for the intercept) and $k$ columns of covariates. We assume that
$\underline{\tilde{\bm{X}}}_{\text{full}}$ has full column rank. In Table A1
of the Online Supplementary Materials we provide a list of the notation used
in this paper.

Based on this notation, $Y_{i}(g^{\prime})-Y_{i}(g^{\prime\prime})$ is the
causal effect of treatment $g^{\prime}$ relative to treatment
$g^{\prime\prime}$ for unit $i$. We are interested in estimating the sample
average treatment effect
$\text{SATE}_{g^{\prime},g^{\prime\prime}}=\frac{1}{N}\sum_{i=1}^{N}\\{Y_{i}(g^{\prime})-Y_{i}(g^{\prime\prime})\\}$
and the population average treatment effect
$\text{PATE}_{g^{\prime},g^{\prime\prime}}=\mathbb{E}\\{Y_{i}(g^{\prime})-Y_{i}(g^{\prime\prime})\\}$.
For this, we will randomly assign the units into treatment groups using the
FSM.

###  3.2 Components of the FSM

In the FSM, the $G$ treatment groups take turns selecting units in a random
but controlled order while optimizing a common criterion. This is accomplished
by the two components of the FSM, namely, the selection order matrix and the
selection function.

  1. 1.

Selection order matrix (SOM): An SOM is a matrix that determines the order in
which the treatment groups select the units. Typically, an SOM has two
columns; the first specifies the stages of selection (from $1$ to $N$), and
the second specifies the treatment group that selects first at that stage.

  2. 2.

Selection function: A selection function is a function that determines which
unit gets selected by the choosing treatment group at each stage. Typically, a
selection function is based on an optimality criterion that is common to all
treatment groups.

A good SOM guarantees that the selection of units is fair, so that no single
treatment group selects all the units of a given type, and random, so that
both observed and unobserved covariates are balanced in expectation and there
is a basis for inference. A good selection function will produce efficient and
robust inferences under a wide class of possible outcome functions.

To illustrate, Table 1(b)(a) presents an example data set with 12 observations
and one covariate, age. We consider assigning these 12 units into two groups
of equal sizes using the FSM. Table 1(b)(b) shows an example of an SOM in this
setting. The SOM determines the order in which each treatment selects a unit
at each stage. In the example, treatment group 2 selects first in stage 1,
treatment group 1 selects in stage 2, and so on. Treatment groups select units
based on the selection function.

Table 1: (a) Example data set; (b) selection order matrix and an assignment
using the FSM.

(a) Data set

Index | Age
---|---
1 | 24
2 | 30
3 | 34
4 | 36
5 | 40
6 | 41
7 | 45
8 | 46
9 | 50
10 | 54
11 | 56
12 | 60
Mean | 43

(b) Selection order matrix and assignment

Selection order matrix | Unit selected
---|---
Stage | Treatment | Index | Age
1 | 2 | 1 | 24
2 | 1 | 12 | 60
3 | 1 | 2 | 30
4 | 2 | 11 | 56
5 | 1 | 3 | 34
6 | 2 | 10 | 54
7 | 1 | 9 | 50
8 | 2 | 4 | 36
9 | 1 | 5 | 40
10 | 2 | 8 | 46
11 | 2 | 6 | 41
12 | 1 | 7 | 45

In general, it is crucial that the order of selection is random, but that no
group chooses in a disproportionate manner. For two treatment groups of
arbitrary sizes, this can be accomplished by means of the Sequentially
Controlled Markovian Random Sampling (SCOMARS) algorithm (Morris 1983). In the
FSM, SCOMARS specifies the probability of a treatment group selecting at stage
$r$ ($r\in\\{1,2,...,N\\}$), conditional on the number of selections made by
that group up to stage $r-1$. See the Online Supplementary Materials for a
formal description of the algorithm. SCOMARS satisfies the sequentially
controlled condition (Morris 1983), which requires the deviation of the
observed number of selections made by a treatment group up to stage $r$ from
its expectation to be strictly less than one. Intuitively, this condition
ensures that throughout the selection process, no treatment group departs too
much from its expected fair share of choices. Moreover, SCOMARS is Markovian
because for each group, the probability of selection at stage $r$ depends
solely on the number of selections made up to stage $r-1$. For two groups of
equal sizes (as in the example in Table 1(b)), generating an SOM under SCOMARS
boils down to successively generating $N/2$ independent random permutations of
the treatment labels $(1,2)$. In Section 8.1 and in the Online Supplementary
Materials, we describe this and other extensions of SCOMARS to multi-group
experiments. Unless otherwise specified, in the rest of the paper, we will use
SCOMARS to generate the SOM for experiments with two treatment groups.

The selection function gives a value to each of the units available for
selection at each stage. This value depends on the characteristics of each
available unit in addition to those already assigned to the treatment group
that selects next. In principle, any criterion can be used in the selection
function. For example, if the selection function is constant, then the
treatment group selects a unit randomly from the available pool.
Alternatively, the selection function can compute the contribution of each
unit to a measure of the accuracy of the estimator. In this spirit, we propose
the D-optimal selection function, which, at each stage, minimizes the
generalized variance of the estimated regression coefficients in a linear
potential outcome model (see Section 4 for details).

To build intuition, in Table 1(b)(b) we discuss the special case of $k=1$
covariate. With the D-optimal selection function, the choosing group, in its
first choice, selects the unit whose covariate value is farthest from the
full-sample mean of the covariate; and in the subsequent choices, selects the
unit whose covariate value is farthest from its current mean of the covariate.
In the example in Table 1(b), treatment $2$ selects unit $1$ with age $24$,
the farthest age from the full-sample mean $43$. In the next stage, treatment
$1$ selects unit $12$ with age $60$, the farthest age from $43$.333Notice that
for treatment 1’s first selection, the mean of age remains 43 (i.e., the full-
sample mean of age) and is not recalculated based on the 11 unselected units.
Next, treatment $1$ selects unit $2$ with age $30$, the farthest age from its
current mean age $60$. The process continues until all the 12 units are
selected.

In general, with multivariate data, the FSM first selects the units that are
farthest from the full-sample mean of the covariates and successively
approaches this target, ultimately selecting the units that are closest to it.
In the FSM, the SOM produces balance out of an optimality criterion that is
common to all the treatment groups. This is crucial so that all the choosers
know the same, and as they choose, they produce groups that are balanced and
equally robust against the unknown outcome model.

Another important feature of the FSM is that, in addition to several treatment
groups, it can accommodate a discard group of unassigned units. This is
important, for example, in settings where the number of available units for
assignment is greater than the total number of units that can feasibly be
assigned (e.g., because of budgetary constraints). This feature of the FSM was
used in the HIE to secure the representativeness of the treatment groups
relative to the target populations.

##  4 The D-optimal selection function

Here, we formalize the D-optimal selection function and provide an equivalent,
closed-form characterization that explains how this criterion governs the
selection of units at each stage. Without loss of generality, assume that
treatment 1 selects at stage $r$, $r\in\\{1,2,...,N\\}$. Let
$\tilde{n}_{r-1}$, $\bar{\bm{X}}_{r-1}$, $\underline{\bm{S}}_{r-1}$, and
$\underline{\tilde{\bm{X}}}_{r-1}$ be the number, mean vector, covariance
matrix, and the design matrix of the units selected after the $(r-1)$th stage
by treatment 1, respectively.

To define the selection function, we consider a linear potential outcome model
of $Y_{i}(1)$ on $\bm{X}_{i}$, i.e.,
$Y_{i}(1)=\bm{\beta}^{\top}(1,\bm{X}_{i}^{\top})^{\top}+\eta_{i}$, where
$\eta_{i}$ is an error term satisfying
$\mathbb{E}\\{\eta_{i}|\bm{X}_{i}\\}=0$.444More generally, one can consider a
linear model of $Y_{i}(1)$ on a vector of basis functions $\bm{B}(\bm{X}_{i})$
of the covariates. Denote $\mathcal{R}_{r-1}$ as the set of unselected units
after stage $r-1$. For unit $i\in\mathcal{R}_{r-1}$, let
$\underline{\tilde{\bm{X}}}_{r,i}$ be the resulting design matrix in treatment
group 1 if unit $i$ is selected. When
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
invertible, the D-optimal selection function selects unit
$i^{\prime}\in\mathcal{R}_{r-1}$, where
$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}\det(\underline{\tilde{\bm{X}}}^{\top}_{r,i}\underline{\tilde{\bm{X}}}_{r,i})$.
In other words, at the $r$th stage, the D-optimal selection function chooses
the unit in $\mathcal{R}_{r-1}$ that optimally decreases the generalized
variance of the estimated regression coefficients of the fitted linear model
in treatment 1. Ties in the values of the generalized variances are resolved
randomly. When
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
not invertible, we define the D-optimal selection function by using a form of
Ridge augmentation (see Lemma A1 in the Online Supplementary Materials). The
following theorem provides an equivalent characterization of the D-optimal
selection function that elucidates the selection made by the choosing
treatment group at each stage.

######  Theorem 4.1.

Assume treatment 1 chooses at stage $r$. Then the D-optimal selection function
chooses unit $i^{\prime}$ such that

$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1})^{\top}(\underline{\bm{S}}^{*}_{r-1})^{-1}(\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1}),$

where

$\bar{\bm{X}}^{*}_{r-1}=\begin{cases}\bar{\bm{X}}_{\text{full}}&\text{if
$\tilde{n}_{r-1}=0$}\\\
\frac{\bar{\bm{X}}_{r-1}+\epsilon\bar{\bm{X}}_{\text{full}}}{1+\epsilon}&\text{if
$\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
not invertible}\\\ \bar{\bm{X}}_{r-1}&\text{if $\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
invertible}\end{cases}$

and

$\underline{\bm{S}}^{*}_{r-1}=\begin{cases}\underline{\bm{S}}_{\text{full}}&\text{if
$\tilde{n}_{r-1}=0$}\\\
(\frac{1}{\tilde{n}_{r-1}}\underline{\bm{X}}_{r-1}^{\top}\underline{\bm{X}}_{r-1}+\frac{\epsilon}{N}\underline{\bm{X}}_{\text{full}}^{\top}\underline{\bm{X}}_{\text{full}})-(1+\epsilon)\bar{\bm{X}}^{*}_{r-1}\bar{\bm{X}}^{*\top}_{r-1}&\text{if
$\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
not invertible}\\\ \underline{\bm{S}}_{r-1}&\text{if $\tilde{n}_{r-1}\geq 1$
and $\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$
is invertible.}\end{cases}$

Theorem 4.1 shows that at every stage, the D-optimal selection function
selects the unit among the remaining pool of available units whose covariate
vector maximizes a type of Mahalanobis distance. In its first choice,
treatment 1 maximizes the Mahalanobis distance from the covariate distribution
in the full sample (in particular, from $\bar{\bm{X}}_{\text{full}}$), thereby
choosing the most outlying unit available in the full sample. For the
subsequent stages where
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
not invertible, treatment 1 maximizes the Mahalanobis distance from a mixture
covariate distribution between treatment group 1 and the full sample, where
$\epsilon$ determines the mixing rate. Finally, the latter selections by
treatment 1 maximize the Mahalanobis distance from the covariate distribution
in treatment group 1. Therefore, with every selection, treatment 1 maximizes
the overall separation of the covariates from its current mean, which
increases the efficiency of the estimated regression coefficients.

By definition, the D-optimal selection function improves the accuracy of the
fitted linear model in each treatment group by sequentially minimizing the
generalized variance of the estimated regression coefficients. With the
D-optimal selection function, we can also establish several additional
desirable properties of the FSM. In particular, leveraging the connection
between D-optimality and Mahalanobis distance, we can show that FSM with the
D-optimal selection function is affine invariant, i.e., the selections of
units by the treatment groups remain unchanged even if the covariates are
transformed linearly. See Section C in the Online Supplementary Materials for
a proof. An implication of this property is that the FSM is invariant with
respect to changes in the location and scale of the covariates.

The FSM with the D-optimal selection function is appealing also because it can
encompass several classical designs, such as randomized blocked and matched-
pair designs. Theorem 4.2 formalizes this result. In the traditional
randomized block design (RBD), the units are grouped into blocks of size $G$
according to a categorical, blocking variable, and each treatment is randomly
applied to exactly one unit within each block (see, e.g., Cox and Reid 2000,
Section 3.4). Here we consider a more general version of an RBD where the
blocks are of size $c\times G$ (where $c$ is a fixed positive integer) and
each treatment is applied to $c$ units within each block. This is a special
case of a stratified randomized experiment with strata of equal size and equal
allocation among treatments per stratum. In a matched-pair design with $G=2$
treatments, similar units are grouped into pairs, and each treatment is
randomly applied to one unit within each pair. This is also a special case of
a stratified randomized experiment with equal allocation per strata, where the
size of each stratum equals two.

######  Theorem 4.2.

  1. (a)

Consider $N=cBG$ units belonging to $B$ blocks of equal size that are to be
randomly assigned into $G$ treatment groups of equal size, where $c$ is a
fixed positive integer. Then, if the linear model in the FSM consists of an
intercept and indicators of any $B-1$ levels of the blocking variable, the FSM
with the D-optimal selection function produces the same assignment as an RBD.

  2. (b)

Consider $N/2$ identical pairs of units in terms of baseline covariates
$\bm{X}_{i}$ that are to be assigned into $G=2$ treatment groups of equal
size. Assume $\bm{X}_{i}$ is drawn from a continuous distribution. Then, if
the linear model in the FSM consists of the intercept and the covariates
$\bm{X}_{i}$, then the FSM almost surely produces the same assignment
mechanism as a matched-pair design.

In the first case, Theorem 4.2(a) states that, by including the levels of a
blocking variable as regressors, the FSM with the D-optimal selection function
automatically blocks on that variable. Thus, the FSM retrieves an RBD without
explicitly performing separate randomizations within each block. In the second
case, Theorem 4.2(b) states that, by including the covariates as regressors,
the FSM with the D-optimal selection function produces the same assignment as
a matched-pair experiment, without explicitly performing separate
randomizations in each pair. This phenomenon is particularly useful when the
sample consists of near-identical twins but that are difficult to identify a
priori due to multiple covariates.

##  5 Inference under the FSM

Using the FSM we can make model- and randomization-based inferences. Both
modes of inference are feasible for any selection function and any randomized
SOM. In model-based inference, the sample is typically assumed to be drawn
randomly from some superpopulation, and inference for the PATE is done by
modeling the observed outcome distribution conditional on the treatment
indicators and the covariates. For instance, let the potential outcome model
under treatment $g$ be
$Y_{i}(g)=\bm{\beta}^{\top}_{g}\bm{B}(\bm{X}_{i})+\epsilon_{ig}$, where
$\bm{B}(\bm{X}_{i})=(B_{1}(\bm{X}_{i}),...,B_{b}(\bm{X}_{i}))^{\top}$ is a
vector of $b$ basis functions of the covariates, and $\epsilon_{ig}$,
$i\in\\{1,2,...,N\\}$ are mutually independent errors, independent of the
covariates. Under this model, $\text{PATE}_{g^{\prime},g^{\prime\prime}}$ can
be unbiasedly estimated by
$\widehat{\text{PATE}}_{g^{\prime},g^{\prime\prime}}=\hat{\bm{\beta}}^{\top}_{g^{\prime}}\overline{\bm{B}(\bm{X})}-\hat{\bm{\beta}}^{\top}_{g^{\prime\prime}}\overline{\bm{B}(\bm{X})}$,
where $\overline{\bm{B}(\bm{X})}=\frac{1}{N}\sum_{i=1}^{N}\bm{B}(\bm{X}_{i})$
and $\hat{\bm{\beta}}_{g}$ is the OLS estimator of $\bm{\beta}_{g}$ obtained
by fitting a linear regression of $Y^{\text{obs}}_{i}$ on $\bm{B}(\bm{X}_{i})$
in treatment group $g=g^{\prime},g^{\prime\prime}$. We call this the
regression imputation estimator of
$\text{PATE}_{g^{\prime},g^{\prime\prime}}$. The standard error of this
estimator and the corresponding confidence interval for
$\text{PATE}_{g^{\prime},g^{\prime\prime}}$ can be obtained using standard OLS
theory. We note that, in model-based inference, the standard errors and
confidence intervals do not take into account the randomness stemming from the
assignment mechanism.

In randomization-based inference, the potential outcomes and the covariates
are typically considered fixed and the assignment mechanism is the only source
of randomness (see Chapter 2 of Rosenbaum 2002 and chapters 5–7 of Imbens and
Rubin 2015 for overviews). Inference for causal effects can be done via exact
randomization tests for sharp null hypotheses on unit-level causal effects
(Fisher 1935), or via estimation under Neyman’s repeated sampling approach
(Neyman 1923, 1990). Under the FSM, randomization tests for sharp null
hypotheses can be performed by approximating the distribution of the test
statistic through repeated realizations of the FSM. To illustrate, consider
testing the sharp null hypothesis of zero unit-level causal effects, i.e.,
$H_{0}:Y_{i}(2)-Y_{i}(1)=0$ for all $i$, at level $\alpha$ using the FSM.
While any choice of test statistic preserves the validity of the test, a
common choice is the absolute difference-in-means statistic
$|\frac{1}{n_{2}}\sum_{i:Z_{i}=2}Y^{\text{obs}}_{i}-\frac{1}{n_{1}}\sum_{i:Z_{i}=1}Y^{\text{obs}}_{i}|=|\frac{1}{n_{2}}\sum_{i:Z_{i}=2}Y_{i}(2)-\frac{1}{n_{1}}\sum_{i:Z_{i}=1}Y_{i}(1)|=:T\\{\bm{Z},\bm{Y}(1),\bm{Y}(2)\\}$.
Large values of $T\\{\bm{Z},\bm{Y}(1),\bm{Y}(2)\\}$ are considered evidence
against $H_{0}$. Under $H_{0}$, $Y_{i}(2)=Y_{i}(1)=Y^{\text{obs}}_{i}$ and the
vectors of potential outcomes $\bm{Y}(1)$ and $\bm{Y}(2)$ are known and fixed.
The $p$-value of the test is given by
$p=P_{H_{0}}[T\\{\bm{Z},\bm{Y}(1),\bm{Y}(2)\\}\geq t_{\text{obs}}]$, where
$t_{\text{obs}}$ is the value of the test statistic for the observed
realization of $\bm{Z}$ under the FSM. We can compute this $p$-value by Monte
Carlo approximation, i.e., we generate independent vectors of assignments
$\bm{Z}^{(m)}=(Z^{(m)}_{1},...,Z^{(m)}_{N})^{\top}$, $m\in\\{1,2,...,M\\}$
using the FSM and approximate the $p$-value as
$\hat{p}=\frac{1}{M}\sum_{m=1}^{M}\mathbbm{1}\big{[}T\\{\bm{Z}^{(m)},\bm{Y}(1),\bm{Y}(2)\\}\geq
t_{\text{obs}}\big{]}$. We reject $H_{0}$ at level $\alpha$ if
$\hat{p}\leq\alpha$.

Similar tests can be applied for more general sharp hypotheses of treatment
effects (e.g., dilated and tobit effects; Rosenbaum 2002, 2010). We can invert
these tests to obtain a confidence interval for the hypothesized effect
(Rosenbaum 2002, Section 2.6.1). Moreover, we can get a point estimate of the
effect by solving a Hodges-Lehmann estimating equation corresponding to these
tests (Rosenbaum 2002, Section 2.7.2). Finally, under Neyman’s approach, we
can estimate the sample average treatment effect
$\text{SATE}_{g^{\prime},g^{\prime\prime}}$ by the difference-in-means
statistic. In particular, for groups of equal size, this difference-in-means
statistic is unbiased for $\text{SATE}_{g^{\prime},g^{\prime\prime}}$ under
the FSM (see Proposition A1 for a proof).

##  6 The Health Insurance Experiment

###  6.1 Data

We evaluate the performance of the FSM relative to other common treatment
assignment approaches using the baseline data of the HIE. To this end, we
consider a version of the HIE data presented in Aron-Dine et al. (2013). This
dataset comprises the six cost-sharing plans described in Section 2. To make
the group sizes more homogeneous, we combine the groups with $25\%$, $50\%$,
and mixed coinsurance plans. Thus, in our analysis, we have $G=4$ treatment
groups corresponding to $g=1$, “free care” ($n_{1}=564$); $g=2$, “$25\%,50\%$,
or mixed coinsurance” ($n_{2}=456$); $g=3$, “$95\%$ coinsurance”
($n_{3}=372$); and $g=4$, “individual deductible” ($n_{4}=495$). In total,
there are $N=n_{1}+...+n_{4}=1,887$ families. We assign all $N$ families to
the four treatment groups (i.e., without a discard group of non-participants).
In this version of the HIE data, we pool the data across five of the six
sites, and we randomly assign all the families to the four treatment groups.
Due to loss of data, the Dayton site is excluded from this analysis.

We consider $k=20$ family-level baseline covariates, where $X_{1},...,X_{5}$
are scaled non-binary covariates, $X_{6},...,X_{14}$ are binary covariates,
and $X_{15},...,X_{20}$ are binary covariates indicating missing data (see
Table A7 for a description of each baseline covariate). Using this data, we
compare complete randomization (CRD), rerandomization (RR), and the FSM in
terms of balance and efficiency. For the FSM, we generate the SOM by first
using SCOMARS on the combined groups $\\{1,2\\}$ and $\\{3,4\\}$, and then
using SCOMARS again to split each combined group into its component groups.
For the FSM, we also use the D-optimal selection function based on a linear
potential outcome model on the main covariates. The assignments under the FSM
are generated using the open source R package FSM available on CRAN. For
rerandomization, we consider two balance criteria, one based on the Wilks’
lambda statistic (RR Wilks; Lock 2011, Section 5.2) and the other based on the
maximum pairwise Mahalanobis distance between any two treatment groups (RR
Mahalanobis; Morgan and Rubin 2012). The balance criteria for both RR Wilks
and RR Mahalanobis are based on all the main covariates and the squares and
pairwise products of the scaled (non-binary) covariates. Finally, for both
rerandomization methods, we use an acceptance rate of 0.001 (Lock 2011). We
draw 400 independent assignments for each approach. The results under RR Wilks
and RR Mahalanobis are roughly the same (see Section H.5 in the Online
Supplementary Materials), and hence, for conciseness, here we only discuss the
results for RR Mahalanobis. The runtime of each of these assignments was
approximately 78 seconds with RR Mahalanobis and 28 seconds with the FSM on a
Windows 64-bit laptop computer with an Intel(R) Core i7 processor.

###  6.2 Balance

Figures 2(a) and 2(b) display the distributions of ASMD across randomizations
for the main covariates and their second-order transformations (squares and
pairwise products). RR balances the main covariates and the second-order terms
better than CRD. However, in both cases, the FSM improves considerably over
CRD and RR. In fact, with the FSM, the average imbalance is less than half
(0.02) of those under CRD and RR. Also, with both CRD and RR, it is common to
see imbalances greater than 0.1 ASMD, whereas such extreme imbalances are non-
existent with the FSM.

(a) Main covariates

(b) Second order terms

(c) Frobenius norm

Figure 2: Distributions of absolute standardized mean differences (ASMD) of
the main covariates (panel (a)) and their squares and pairwise products (panel
(b)) across randomizations. For each plot, the legend presents the average
ASMD across simulations for each method. Panel (c) shows the distributions of
discrepancies between the correlation matrices of the covariates in treatment
groups 1 and 2, as measured by the Frobenius norm,
$||\underline{\bm{R}}_{1}-\underline{\bm{R}}_{2}||_{F}$. In terms of the main
covariates, second-order transformations, and correlation matrices, the FSM
substantially outperforms CRD and RR.

A related question is how well the methods balance all second-order features
of the joint distribution of the covariates. Figures 2(c) and A4 provide an
answer to this question in the boxplots of the discrepancies between
correlation matrices across randomizations. As a measure of discrepancy, we
consider the Frobenius norm of the difference between correlation matrices in
two groups, i.e.,
$||\underline{\bm{R}}_{g}-\underline{\bm{R}}_{g^{\prime}}||_{F}$, where
$\underline{\bm{R}}_{g}$ is the sample correlation matrix in group $g$ and
$||\cdot||_{F}$ is the Frobenius norm.555The Frobenius norm of a matrix is the
square root of the sum of squares of all its elements. Smaller values of
$||\underline{\bm{R}}_{g}-\underline{\bm{R}}_{g^{\prime}}||_{F}$ indicate
better balance on the correlation matrix of the covariates between the groups
$g$ and $g^{\prime}$. As in the aforementioned second-order transformations,
we see a similar performance between complete randomization and
rerandomization, which is considerably improved by the FSM with a median about
three times smaller.

###  6.3 Efficiency

In this section, we evaluate the estimation accuracy of the methods under
model- and randomization-based approaches to inference. The main differences
between the model- and randomization-based standard errors is that in the
model-based approach, the variance calculation does not explicitly take into
account the variability arising through the randomization distribution,
whereas in the randomization-based approach it does. For illustration, here we
consider estimating the average treatment effect of treatment 3 relative to
treatment 2, i.e., $\text{SATE}_{3,2}$ and $\text{PATE}_{3,2}$. The results
for the average treatment effects with other pairs of treatment groups are
similar.

Under the model-based approach, we consider two potential outcome models, one
that is linear on the main covariates (Model A1), and another that is linear
on the main covariates and the second-order transformations of the scaled
covariates (Model A2). The results are summarized in Table 4. While the
performance of the three methods is similar under Model A1, under Model A2
there are substantial differences, with the FSM outperforming both complete
randomization and rerandomization. In fact, under Model A2, there is a 14-15%
reduction in the average standard error, and a 53-64% reduction in the maximum
standard error, with the FSM.

Table 2: Average and maximum model-based standard errors relative to the FSM
across randomizations. Under Model A1 (linear model on the covariates), the
FSM is slightly more efficient than RR and CRD. Under Model A2 (linear model
on the covariates and their second-order transformations), the FSM is
considerably more efficient than CRD and RR.

(a) Model A1

| Designs
---|---
CRD | RR Mahalanobis | FSM
Average SE | 1.02 | 1.01 | 1.00
Maximum SE | 1.04 | 1.02 | 1.00

(b) Model A2

| Designs
---|---
CRD | RR Mahalanobis | FSM
Average SE | 1.15 | 1.14 | 1.00
Maximum SE | 1.64 | 1.53 | 1.00

Under the randomization-based approach, we consider the generative models
$Y(3)=10+2X_{1}+3X_{2}+0.5X_{3}+0.3X_{4}+\eta$ (Model B1) and
$Y(3)=10+2X_{1}+2X_{2}X_{3}-X_{4}X_{5}+\eta$ (Model B2) where $Y(3)=Y(2)$ and
$\eta\sim\mathcal{N}(0,1.5^{2})$. Here, both the generative models satisfy the
sharp-null hypothesis of zero treatment effect for every unit and hence,
$\text{SATE}_{3,2}=0$. Under each design, $\text{SATE}_{3,2}$ is estimated
using the standard difference-in-means estimator and the corresponding
randomization-based SE is obtained by generating 400 randomizations and
computing the standard deviation of the estimator across these 400
randomizations. The results are summarized in Table 3(b). See Appendix H.5 for
similar comparisons under a set of different generative models of the
potential outcome. In terms of efficiency, we see again a clear advantage of
the FSM. Under both Model B1 and Model B2, the average standard errors of
complete randomization and rerandomization are more than twice of those under
the FSM.

Table 3: Randomization-based standard errors relative to the FSM. The standard
error for the FSM is 0.11 under Model B1 (linear model on the covariates) and
0.64 under Model B2 (linear model on the covariates and their second-order
transformations). Under both models, the FSM is considerably more efficient
than both CRD and RR.

(a) Model B1

| Designs
---|---
CRD | RR Mahalanobis | FSM
SE | 2.47 | 2.08 | 1

(b) Model B2

| Designs
---|---
CRD | RR Mahalanobis | FSM
SE | 2.63 | 2.25 | 1

###  6.4 Intuition and further explorations

Our analysis illustrates some important differences between the FSM, CRD, and
RR. With respect to RR, these differences pertain to the specification, role,
and implementation of the assignment criterion. First, regarding the
specification of the criterion, while RR uses the Mahalanobis distance, the
FSM uses the D-optimality criterion, which, coupled with a suitable SOM, leads
to robust assignments under a more general class of potential outcome models.

Second, regarding the role of this criterion, while RR essentially constrains
the allowable treatment assignments, the FSM seeks to optimize them toward the
criterion. In essence, while RR solves a feasibility problem by resampling,
the FSM aims to solve a maximization problem by step-wise assignment.
Furthermore, the feasibility problem solved by RR depends on the balance
threshold, which can be difficult to select in practice. While a very high
threshold can accept assignments with poor covariate balance, a very low one
can be computationally onerous.

Third, regarding the implementation of the criterion, while RR assigns all
units in one step and then discards imbalanced assignments, the FSM assigns
units in multiple steps (one at a time) in a random but optimal fashion
determined by the selection order and the selection criterion. This difference
is crucial because in experiments like the HIE with several treatment groups
and many covariates, the space of possible treatment assignments is vast. As
shown in our analyses, optimally selecting among these assignments in a step-
wise manner can make a substantial improvement in terms of balance,
efficiency, computational time, and, ultimately, in the use of scarce
resources available for experimentation. 666Figures 1 and 2 show that,
although RR does well under common balance standards (the mean differences are
systematically lower than the typical threshold of 0.1 ASMD), there is room to
select better (more balanced) random treatment assignments, which is achieved
by the FSM.

To better see this, we asked how we would need to modify RR to achieve
comparable performance to the FSM? Using the HIE data, we approximated the
randomization distribution of the imbalance criterion of RR (i.e., the maximum
Mahalanobis distance $M$ across all pairs of treatment groups) by generating
random assignments for 100 hours. See Table 4 for a summary of the results.
The table displays summary statistics of the distribution of $M$ under CRD,
RR, and the FSM. As shown in Table 4, the highest (worst-case) value of $M$
under the FSM is smaller than the smallest (best-case) value of $M$ under CRD
and RR. Importantly, even if we set the RR acceptance rate to 0.0000001 (i.e.,
1 over 10 million), we still have imbalances higher than the worst-case
imbalance of the FSM. In sum, even with an acceptance rate as low as
0.0000001, RR did not perform as well as the FSM, despite taking 100 hours on
average to generate a single assignment, as opposed to the 30 seconds of
running time of the FSM.

Design | Minimum | 1st Quartile | Median | Mean | 3rd Quartile | Maximum
---|---|---|---|---|---|---
CRD | 18.5 | 39.5 | 43.9 | 44.4 | 48.7 | 96.1
RR (0.001) | 18.5 | 25.4 | 26.2 | 25.9 | 26.7 | 27.1
FSM | 2.8 | 4.7 | 5.3 | 5.4 | 6.0 | 10.6

Table 4: Distribution of the maximum pairwise Mahalanobis distance across
groups ($M$). For CRD, we obtain this distribution by generating over 10
million random assignments for 100 hours. For RR (0.001), we obtain this
distribution using 0.1% of all these assignments with the smallest values of
$M$. For the FSM, we obtain this distribution using the 400 random assignments
from Section 6.1.

##  7 Ten further studies in the health and social sciences

In addition to the previous study, we evaluate the performance of the FSM in
ten randomized studies from the health and social sciences. These ten studies
are labelled (1) Crepon, which evaluates the impact of a microcredit program
in rural Morocco on assets, profits, and consumption (Crépon et al. 2015); (2)
Angrist, which evaluates the impact of cash incentives on certification rates
among low-achievers in Israel (Angrist and Lavy 1999); (3) Finkelstein, which
evaluates the impact of the Camden Coalition of Healthcare Providers’
Hotspotting program on hospital readmission rates among patients with high use
of healthcare services (Finkelstein et al. 2020); (4) Durocher, which
evaluates the impact of intravenous infusion versus intramusculur oxytocin on
postpartum blood loss and hemmorhage rates (Durocher et al. 2019); (5)
Lalonde, which evaluates the impact of Nationally Supported Work program on
earnings (LaLonde 1986); (6) Karlan, which evaluates the impact of loans with
an indemnity component on demand for credit and investment decisions of
farmers (Karlan et al. 2014); (7) Dupas, which evaluates the impact of
different cost provisions for allocating dilute-chlorine water treatment
solution on chlorine residuals in households’ stored water (Dupas et al.
2016); (8) Blattman, which evaluates the impact of industrial job offers and
entrepreneurial programs on health, income and other measures (Blattman and
Dercon 2018); (9) Ambler, which evaluates the impact of offering Salvadoran
migrant maching funds for educational remittances on educational investments
and other outcomes Ambler et al. (2015); (10) Wantchekon, which evaluates the
impact of townhall meeting based on programmatic, nonclientelist platforms on
clientelism, voter turnout, and vote shares (Fujiwara and Wantchekon 2013).
Table 5 provides details on the design parameters considered in these studies.

Study | Design parameters | Main covariates | Second-order transformations
---|---|---|---
$N$ | $G$ | $(n_{1},...,n_{G})$ |  $k$ | CRD | RR | FSM | $\frac{\text{CRD}}{\text{FSM}}$ |  $\frac{\text{RR}}{\text{FSM}}$ | CRD | RR | FSM | $\frac{\text{CRD}}{\text{FSM}}$ | $\frac{\text{RR}}{\text{FSM}}$
Crepon | 4465 | 2 | (2266, 2199) |  33 | 0.024 | 0.018 | 0.015 | 1.6 |  1.2 | 0.024 | 0.023 | 0.018 | 1.3 | 1.3
Angrist | 3821 | 2 | (1910,1911) |  20 | 0.025 | 0.014 | 0.002 | 12.5 |  7.0 | 0.026 | 0.023 | 0.003 | 8.7 | 7.7
Finkelstein | 782 | 2 | (389,393) |  10 | 0.062 | 0.020 | 0.010 | 6.2 |  2.0 | 0.059 | 0.048 | 0.013 | 4.5 | 3.7
Durocher | 480 | 2 | (239,241) |  12 | 0.072 | 0.031 | 0.017 | 4.2 |  1.8 | 0.073 | 0.068 | 0.022 | 3.3 | 3.1
Lalonde | 445 | 2 | (222,223) |  10 | 0.083 | 0.044 | 0.014 | 5.9 |  3.1 | 0.077 | 0.070 | 0.019 | 4.1 | 3.7
Karlan | 169 | 2 | (84, 85) |  16 | 0.124 | 0.059 | 0.053 | 2.3 |  1.1 | 0.123 | 0.119 | 0.060 | 2.1 | 2.0
Dupas | 1118 | 3 | (351, 382, 385) |  11 | 0.059 | 0.018 | 0.010 | 5.9 |  1.8 | 0.058 | 0.044 | 0.017 | 3.4 | 2.6
Blattman | 947 | 3 | (358,304,285) |  34 | 0.064 | 0.048 | 0.026 | 2.5 |  1.8 | 0.065 | 0.064 | 0.036 | 1.8 | 1.8
Ambler | 991 | 4 | (360, 211, 203, 217) |  16 | 0.073 | 0.053 | 0.015 | 4.9 |  3.5 | 0.073 | 0.071 | 0.017 | 4.3 | 4.2
Wantchekon | 24 | 2 | (12, 12) |  10 | 0.334 | 0.170 | 0.245 | 1.4 |  0.7 | 0.333 | 0.289 | 0.237 | 1.4 | 1.2
Average |  |  |  |  | 0.092 | 0.048 | 0.041 |  |  | 0.091 | 0.082 | 0.044 |  |
Average* |  |  |  |  | 0.065 | 0.034 | 0.018 |  |  | 0.064 | 0.059 | 0.023 |  |

Table 5: Design parameters and balance results for ten case studies in the
health and social sciences. The second average denoted with an asterisk (∗)
excludes the Wantchekon study because
$\underline{\tilde{\bm{X}}}^{\top}_{r}\underline{\tilde{\bm{X}}}_{r}$ matrix
is non-invertible for the first $r=22$ selections.

For each study, we generate 100 assignments of complete randomization (CRD),
Rerandomization with Mahalanobis distance (based on the main covariates) and
0.001 acceptance rate (RR), and the FSM (based on the main covariates). The
mean ASMD of the main covariates and their squares and interactions under each
method are presented in Table 5. See figures A7 and A8 in the Online
Supplementary Materials for plots of the distributions of these imbalances,
alongside the Frobenius norms of
$\underline{\bm{R}}_{1}-\underline{\bm{R}}_{2}$.777Groups 1 and 2 are chosen
haphazardly as a typical pair of groups. The results for the other pairs of
groups are similar.

Table 5 shows that for each study, CRD achieves a similar mean balance on the
main covariates and their squares and interactions. RR improves balance over
CRD considerably for the main covariates, but only mildly for the squares and
interactions. By contrast, for almost all the studies, the FSM substantially
improves balance over CRD and RR in terms of both the main covariates and
their transformations. The only exception is the Wantchekon study, where the
group sizes barely exceed the number of covariates $k=10$. The FSM is not
designed for settings like this, where the number of covariates is greater
than or close to the minimum treatment group size. In such settings, the
matrix $\underline{\tilde{\bm{X}}}^{\top}_{r}\underline{\tilde{\bm{X}}}_{r}$
is non-invertible for almost every selection stage of the FSM, and therefore,
the D-optimal selection function in the FSM relies on ridge augmentation to
feasibly select units (see Section 4), producing suboptimal selections.

Across the ten studies, the ASMD on the main covariates are 55%
($=\frac{0.092-0.041}{0.092}$) and 15% ($=\frac{0.048-0.041}{0.048}$) lower on
average with the FSM than CRD and RR, respectively. If we exclude Wantchekon,
then these percent reductions in ASMD are amplified to 72% and 47%. Similarly,
across the ten studies, the ASMD on the squares and interactions of the
covariates with the FSM are about 50% smaller than both CRD and RR, and
without Wantchekon, they are at least 60% smaller. In fact, FSM has better
balance on both the main covariates and their second-order transformations
over CRD and RR uniformly across the first nine studies (as shown by the
$\frac{\text{CRD}}{\text{FSM}}$ and $\frac{\text{RR}}{\text{FSM}}$ columns).
For each study, the relative improvement in balance under the FSM over RR is
larger for the second-order transformations than for the main covariates. In
particular, for half of the ten studies, the mean ASMD of the second-order
transformations under RR are at least three times larger than those under the
FSM, implying substantial improvement in balance on these transformations
under the FSM.

Overall, averaging the ASMD of the main covariates and their second-order
transformations across the first nine studies, we see that the FSM achieves
68% better covariate balance than complete randomization and 56% better
covariate balance than rerandomization in a typical study. Across these
studies, the FSM’s performance relative to CRD and RR is consistent with those
in the HIE study in Section 6.888Notably, the relative performances of the
methods in the HIE study are comparable to those of the Ambler study, which
involves roughly half the sample size of the HIE study and similar values of
the other design parameters. For instance, the average ASMD of the main
covariates under CRD, RR, and the FSM in the Ambler study are roughly
$\sqrt{2}$ times those in the HIE study, where $\sqrt{2}$ is the factor that
corrects for the difference in sample size. A similar pattern to the HIE study
is also noted in the plots of
$||\underline{\bm{R}}_{1}-\underline{\bm{R}}_{2}||_{F}$ in figures A7 and A8
in the Online Supplementary Materials, where for most studies, the worst
(least balanced) assignment among all the draws of the FSM has a better
balance on the correlation matrices than the best (most balanced) assignment
among all the draws of CRD and RR. As discussed in Section 6.3, under both
model- and randomization-based approaches to inference, better balance
directly translates to more efficient estimates of treatment effects.

Therefore, similar to the HIE case study, these results show that across a
range of randomized experiments, the FSM is a flexible and robust approach to
randomization.

##  8 Practical considerations and extensions

###  8.1 Multi-group experiments

As discussed, the FSM can readily handle experiments with multiple treatment
groups. In so doing, the key methodological consideration is the choice of an
SOM. As in two-group experiments, we would like to generate an SOM that is
randomized and sequentially controlled, so that at every stage of the random
selection process, the number of selections made by each treatment group up to
that stage is close to its fair share. Constructing a sequentially controlled
SOM for multi-group experiments with arbitrary group sizes is an open problem.
However, such constructions are possible for several practically relevant
configurations of the group sizes, namely (a) groups of equal size, (b) groups
having one of two distinct sizes, and (c) groups of more than two distinct
sizes such that when combined by groups of equal size they have the same total
size. In the Online Supplementary Materials, we provide algorithms to
construct an SOM for all three configurations and prove that the resulting SOM
is sequentially controlled. In practice, for more general group size
configurations, one strategy to generate an SOM is to first identify one of
these three configurations that is structurally similar to the configuration
at hand, and then use the corresponding SOM-generating algorithm. The
resulting algorithm may not always be sequentially controlled, but is still
likely to produce a well-controlled randomized selection order.

###  8.2 Stratified experiments

In stratified experiments, units are grouped into two or more strata, and
within each stratum, units are randomly assigned to treatment. Here we propose
a family of extensions of the FSM to such settings. Typically, in stratified
experiments the treatment group sizes within each stratum are pre-specified by
the investigator. The main challenge arises when the treatment group sizes
differ across strata. To address this challenge, we construct an augmented SOM
with information of the treatment group that selects at each stage and the
stratum that it selects from. This construction guarantees that each treatment
group is assigned the pre-specified number of units in each stratum. In the
Online Supplementary Materials, we discuss two approaches to construct such an
SOM. At a high level, one approach generates a separate SOM for each stratum,
while the other approach uses SCOMARS to determine the order of stratum labels
for each treatment.

###  8.3 Sequential experiments

Sequential experiments are experiments where units progressively become
available for random assignment, possibly in batches of varying sizes. Here we
describe extensions of the FSM to such settings. The simplest approach is to
run an independent FSM for each new batch of available units. However, in
general, this approach fails to account for accrued covariate imbalances
between the treatment groups. To address this issue, we propose an alternative
approach that considers the new batch as a continuation of the previous one.
More specifically, for each unit in the new batch we evaluate the value of the
D-optimal selection function using all the units already assigned to the
selecting treatment group. In this way, this approach tends to remove accrued
covariate imbalances. See the Online Supplementary Materials for technical
details. In sequential experiments, the balance and efficiency of the FSM
assignments tend to increase with the batch sizes.

##  9 Summary and remarks

We revisited, formalized, and extended the FSM for experimental design. We
proposed a new selection function based on D-optimality that requires no
tuning parameters. We showed that, equipped with this selection function, the
FSM has a number of appealing properties. First, the FSM is affine invariant
and hence, it self-standardizes covariates with possibly different units of
measurements. Second, the FSM produces randomized block designs without
explicitly randomizing in each block. Third, the FSM also produces matched-
pair designs without explicitly constructing the matched pairs beforehand and
randomizing within each pair. We described how both model-based and
randomization-based inference on treatment effects can be conducted using the
FSM. For a range of practically relevant configurations of group sizes in
multi-group experiments, we proposed new algorithms to generate a fair and
random selection order of treatments under the FSM. We also discussed
potential extensions of the FSM to stratified and sequential experiments. In a
case study on the RAND Health Insurance Experiment, and ten additional
randomized studies from the health and social sciences, we showed that the FSM
is a robust approach to randomization, exhibiting better performance than
complete randomization and rerandomization in terms of balance and efficiency.

While there are settings where complete randomization may perform better than
the FSM in terms of efficiency, such settings are less common and involve
jagged, i.e., highly non-smooth, potential outcome models. In settings where
these models are reasonably smooth, the FSM is expected to perform well.
Overall, through our extensive explorations with real and simulated
experimental data, the FSM has consistently stood out as a robust design that
can handle multiple treatment groups and a fairly large number of categorical
and continuous covariates without requiring tuning parameters and nor
coarsening covariates. We recommend giving strong considerations to the FSM in
experimental design for its conceptual simplicity, practicality, balance, and
robustness.

## References

  * Ambler et al. (2015) Ambler, K., Aycinena, D., and Yang, D. (2015), “Channeling remittances to education: A field experiment among migrants from El Salvador,” American Economic Journal: Applied Economics, 7, 207–32.
  * Angrist and Lavy (1999) Angrist, J. D. and Lavy, V. (1999), “Using Maimonides’ rule to estimate the effect of class size on scholastic achievement,” The Quarterly Journal of Economics, 114, 533–575.
  * Aron-Dine et al. (2013) Aron-Dine, A., Einav, L., and Finkelstein, A. (2013), “The RAND health insurance experiment, three decades later,” Journal of Economic Perspectives, 27, 197–222.
  * Banerjee et al. (2017) Banerjee, A. V., Chassang, S., and Snowberg, E. (2017), “Decision theoretic approaches to experiment design and external validity,” in Handbook of Economic Field Experiments, Elsevier, vol. 1, pp. 141–174.
  * Bertsimas et al. (2015) Bertsimas, D., Johnson, M., and Kallus, N. (2015), “The power of optimization over randomization in designing experiments involving small samples,” Operations Research, 63, 868–876.
  * Blattman and Dercon (2018) Blattman, C. and Dercon, S. (2018), “The impacts of industrial and entrepreneurial work on income and health: Experimental evidence from Ethiopia,” American Economic Journal: Applied Economics, 10, 1–38.
  * Brook et al. (2006) Brook, R. H., Keeler, E. B., Lohr, K. N., Newhouse, J. P., Ware, J. E., Rogers, W. H., Davies, A. R., Sherbourne, C. D., Goldberg, G. A., Camp, P., et al. (2006), “The health insurance experiment: a classic RAND study speaks to the current health care reform debate,” Santa Monica, CA: RAND Corporation.
  * Chattopadhyay et al. (2021) Chattopadhyay, A., Morris, C. N., and Zubizarreta, J. R. (2021), “Randomized and Balanced Allocation of Units into Treatment Groups Using the Finite Selection Model for R,” arXiv preprint arXiv:2105.02393.
  * Cochran and Cox (1957) Cochran, W. and Cox, G. (1957), Experimental Designs, John Wiley & Sons New York.
  * Cox and Reid (2000) Cox, D. R. and Reid, N. (2000), The Theory of the Design of Experiments, CRC Press.
  * Crépon et al. (2015) Crépon, B., Devoto, F., Duflo, E., and Parienté, W. (2015), “Estimating the impact of microcredit on those who take it up: Evidence from a randomized experiment in Morocco,” American Economic Journal: Applied Economics, 7, 123–50.
  * Dupas et al. (2016) Dupas, P., Hoffmann, V., Kremer, M., and Zwane, A. P. (2016), “Targeting health subsidies through a nonprice mechanism: A randomized controlled trial in Kenya,” Science, 353, 889–895.
  * Durocher et al. (2019) Durocher, J., Dzuba, I. G., Carroli, G., Morales, E. M., Aguirre, J. D., Martin, R., Esquivel, J., Carroli, B., and Winikoff, B. (2019), “Does route matter? Impact of route of oxytocin administration on postpartum bleeding: A double-blind, randomized controlled trial,” PloS one, 14, e0222981.
  * Finkelstein et al. (2020) Finkelstein, A., Zhou, A., Taubman, S., and Doyle, J. (2020), “Health care hotspotting—a randomized, controlled trial,” New England Journal of Medicine, 382, 152–162.
  * Fisher (1925) Fisher, R. A. (1925), “Statistical methods for research workers, 13e,” London: Oliver and Loyd, Ltd, 99–101.
  * Fisher (1935) — (1935), The Design of Experiments, London: Oliver & Boyd.
  * Fujiwara and Wantchekon (2013) Fujiwara, T. and Wantchekon, L. (2013), “Can informed public deliberation overcome clientelism? Experimental evidence from Benin,” American Economic Journal: Applied Economics, 5, 241–55.
  * Greevy et al. (2004) Greevy, R., Lu, B., Silber, J. H., and Rosenbaum, P. R. (2004), “Optimal multivariate matching before randomization,” Biostatistics, 5, 263–275.
  * 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.
  * Imbens and Rubin (2015) Imbens, G. W. and Rubin, D. B. (2015), Causal Inference in Statistics, Social, and Biomedical Sciences, Cambridge University Press.
  * Karlan et al. (2014) Karlan, D., Osei, R., Osei-Akoto, I., and Udry, C. (2014), “Agricultural decisions after relaxing credit and risk constraints,” The Quarterly Journal of Economics, 129, 597–652.
  * Krieger et al. (2019) Krieger, A. M., Azriel, D., and Kapelner, A. (2019), “Nearly random designs with greatly improved balance,” Biometrika, 106, 695–701.
  * LaLonde (1986) LaLonde, R. J. (1986), “Evaluating the econometric evaluations of training programs with experimental data,” The American Economic Review, 604–620.
  * Lock (2011) Lock, K. F. (2011), Rerandomization to Improve Covariate Balance in Randomized Experiments, Harvard University.
  * 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.
  * Morris (1979) Morris, C. (1979), “A finite selection model for experimental design of the health insurance study,” Journal of Econometrics, 11, 43–61.
  * Morris (1983) — (1983), “Sequentially controlled Markovian random sampling (SCOMARS),” Institute of Mathematical Statistics Bulletin, 12, 237.
  * Morris and Hill (2000) Morris, C. and Hill, J. (2000), “The health insurance experiment: design using the finite selection model,” Public Policy and Statistics: Case Studies from RAND, Springer Science & Business Media, 29–53.
  * Newhouse et al. (1993) Newhouse, J. P. et al. (1993), Free for All?, Harvard University Press.
  * Neyman (1923, 1990) Neyman, J. (1923, 1990), “On the application of probability theory to agricultural experiments,” Statistical Science, 5, 463–480.
  * Rosenbaum (2002) Rosenbaum, P. R. (2002), Observational Studies, Springer.
  * Rosenbaum (2010) — (2010), “Design sensitivity and efficiency in observational studies,” Journal of the American Statistical Association, 105, 692–702.
  * Rosenbaum and Rubin (1985) Rosenbaum, P. R. and Rubin, D. B. (1985), “Constructing a control group using multivariate matched sampling methods that incorporate the propensity score,” The American Statistician, 39, 33–38.
  * Rubin (1974) Rubin, D. B. (1974), “Estimating causal effects of treatments in randomized and nonrandomized studies.” Journal of Educational Psychology, 66, 688.
  * Stuart (2010) Stuart, E. A. (2010), “Matching methods for causal inference: a review and a look forward,” Statistical Science, 25, 1–21.

## Supplementary materials

###  A Notation and estimands

Table A1: Notation $N$ | $\triangleq$ |  Full sample size
---|---|---
$i$ | $\triangleq$ |  Index of unit, $i=1,...,N$
$G$ | $\triangleq$ |  Number of treatments
$g$ | $\triangleq$ |  Index of treatment group, $g=1,2,...,G$
$n_{g}$ | $\triangleq$ |  Size of treatment group $g$
$k$ | $\triangleq$ |  Number of baseline covariates
$\boldsymbol{X}_{i}$ | $\triangleq$ |  Observed vector of baseline covariates of unit $i$
$\underline{\bm{X}}_{\text{full}}$ | $\triangleq$ |  $N\times k$ matrix of covariates in the full sample
$\underline{\tilde{\bm{X}}}_{\text{full}}$ | $\triangleq$ |  $N\times k+1$ design matrix in the full sample
$\bar{\bm{X}}_{\text{full}}$ | $\triangleq$ |  $k\times 1$ vector of means of the baseline covariates in the full sample
$\underline{\bm{S}}_{\text{full}}$ | $\triangleq$ |  $k\times k$ covariance matrix of the baseline covariates in the full sample
$Y_{i}(g)$ | $\triangleq$ |  Potential outcome of unit $i$ under treatment $g$
$\bm{Y}(g)$ | $\triangleq$ |  Vector of potential outcomes under treatment $g$, $(Y_{1}(g),...,Y_{N}(g))^{\top}$
$Z_{i}$ | $\triangleq$ |  Treatment assignment indicator of unit $i$, $Z_{i}\in\\{1,2,...,G\\}$
$\bm{Z}$ | $\triangleq$ |  Vector of treatment assignment indicators, $(Z_{1},...,Z_{N})^{\top}$
$Y^{\text{obs}}_{i}$ | $\triangleq$ |  Observed outcome of unit $i$, $Y^{\text{obs}}_{i}=\sum_{g=1}^{G}\mathbbm{1}(Z_{i}=g)Y_{i}(g)$
Table A2: Estimands $Y_{i}(g^{\prime})-Y_{i}(g^{\prime\prime})$ | $\triangleq$ |  Unit level causal effect of treatment $g^{\prime}$ relative to treatment $g^{\prime\prime}$ for unit $i$; $g^{\prime},g^{\prime\prime}\in\\{1,2,...,G\\}$
---|---|---
$\text{SATE}_{g^{\prime},g^{\prime\prime}}$ | $\triangleq$ |  $\frac{1}{N}\sum_{i=1}^{N}\\{Y_{i}(g^{\prime})-Y_{i}(g^{\prime\prime})\\}$, the Sample Average Treatment Effect of treatment $g^{\prime}$ relative to treatment $g^{\prime\prime}$
$\text{PATE}_{g^{\prime},g^{\prime\prime}}$ | $\triangleq$ |  $\mathbb{E}\\{Y_{i}(g^{\prime})-Y_{i}(g^{\prime\prime})\\}$, the Population Average Treatment Effect of treatment $g^{\prime}$ relative to treatment $g^{\prime\prime}$

###  B Proofs of theoretical results

######  Lemma A1.

Let treatment 1 be the choosing group at the $r$th stage. Also, let
$\underline{\tilde{\bm{X}}}_{r-1}$ be the $\tilde{n}_{r-1}\times(k+1)$ design
matrix in treatment group 1 after the $(r-1)$th stage, where
$\tilde{n}_{r-1}\geq 1$ and
$\text{rank}(\underline{\tilde{\bm{X}}}_{r-1})=k+1$. The D-optimal selection
function chooses unit $i^{\prime}$ with covariate vector
$\bm{X}_{i^{\prime}}\in\mathbb{R}^{k}$, where

$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$ (A1)

###### Proof.

We follow the notations outlined in Section 4. At the $r$th stage, D-optimal
selection function selects unit $i^{\prime}\in\mathcal{R}_{r-1}$, where
$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}\det(\underline{\tilde{\bm{X}}}^{\top}_{r,i}\underline{\tilde{\bm{X}}}_{r,i})$.
Now, for $i\in\mathcal{R}_{r-1}$,

$\displaystyle\det(\underline{\tilde{\bm{X}}}^{\top}_{r,i}\underline{\tilde{\bm{X}}}_{r,i})$
$\displaystyle=\det\Big{\\{}\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}+\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)(1,\bm{X}^{\top}_{i})\Big{\\}}$ (A2)
$\displaystyle=\det(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})\det\Big{\\{}\underline{\bm{I}}+(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-\frac{1}{2}}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-\frac{1}{2}}\Big{\\}}$
(A3)
$\displaystyle=\det(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})\Big{\\{}1+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)\Big{\\}},$ (A4)

where the final equality holds since for two matrices
$\underline{\bm{A}}_{m\times n}$ and $\underline{\bm{B}}_{n\times m}$,
$\det(\underline{\bm{I}}_{m}+\underline{\bm{A}}\underline{\bm{B}})=\det(\underline{\bm{I}}_{n}+\underline{\bm{B}}\underline{\bm{A}})$.
Equation A4 implies that the selected unit $i^{\prime}$ maximizes
$(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$. This completes the proof.

∎

#### Proof of Theorem 4.1

###### Proof.

We use the notations in Section 3.1 and Table A1. We first consider the case
where $\tilde{n}_{r-1}=0$. The selected unit $i^{\prime}$ satisfies,

$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right).$ (A5)

Now, denoting $\bm{e}_{1}=(1,0,...,0)$ as the $k\times 1$ first standard unit
vector, we have

$\displaystyle(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$
$\displaystyle=(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}1\\\
\bar{\bm{X}}_{\text{full}}\end{smallmatrix}\right)+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}0\\\
\bm{X}_{i}-\bar{\bm{X}}_{\text{full}}\end{smallmatrix}\right)$
$\displaystyle=(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}}\frac{\bm{e}_{1}}{N}+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}0\\\
\bm{X}_{i}-\bar{\bm{X}}_{\text{full}}\end{smallmatrix}\right)$
$\displaystyle=\frac{1}{N}+\\{0,(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}})^{\top}\\}(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}0\\\
\bm{X}_{i}-\bar{\bm{X}}_{\text{full}}\end{smallmatrix}\right)$
$\displaystyle=\frac{1}{N}+\frac{1}{N}(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}})^{\top}(\underline{\bm{S}}_{\text{full}})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}}).$
(A6)

Here the last equality holds since, by the formula for the inverse of a
partitioned matrix,
$(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}=\left(\begin{smallmatrix}\underline{\bm{B}}_{11}&\underline{\bm{B}}_{12}\\\
\underline{\bm{B}}_{21}&\underline{\bm{B}}_{22}\\\ \end{smallmatrix}\right)$,
where
$\underline{\bm{B}}^{-1}_{22}=\underline{\bm{X}}^{\top}_{\text{full}}\underline{\bm{X}}_{\text{full}}-N\bar{\bm{X}}_{\text{full}}\bar{\bm{X}}^{\top}_{\text{full}}=N\underline{\bm{S}}_{\text{full}}.$
This completes the proof of the $\tilde{n}_{r-1}=0$ case. The proof for the
case where $\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
invertible follows similar steps and hence is omitted.

We now consider the case where $\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
not invertible. We denote
$\bar{\bm{X}}^{*}_{r-1}=\frac{\bar{\bm{X}}_{r-1}+\epsilon\bar{\bm{X}}_{\text{full}}}{1+\epsilon}$
and
$\underline{\bm{S}}^{*}_{r-1}=(\frac{1}{\tilde{n}_{r-1}}\underline{\bm{X}}_{r-1}^{\top}\underline{\bm{X}}_{r-1}+\frac{\epsilon}{N}\underline{\bm{X}}_{\text{full}}^{\top}\underline{\bm{X}}_{\text{full}})-(1+\epsilon)\bar{\bm{X}}^{*}_{r-1}\bar{\bm{X}}^{*\top}_{r-1}$.
The selected unit $i^{\prime}$ satisfies,

$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(1,\bm{X}^{\top}_{i})\Big{(}\frac{1}{\tilde{n}_{r-1}}\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}+\frac{\epsilon}{N}\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}}\Big{)}^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$ (A7)

Denoting
$\underline{\tilde{\bm{G}}}=\left(\begin{smallmatrix}\sqrt{\frac{1}{\tilde{n}_{r-1}}}\underline{\bm{X}}_{r-1}\\\
\sqrt{\frac{\epsilon}{N}}\underline{\bm{X}}_{\text{full}}\end{smallmatrix}\right)$,
we have

$\displaystyle(1,\bm{X}^{\top}_{i})\Big{(}\frac{1}{\tilde{n}_{r-1}}\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}+\frac{\epsilon}{N}\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}}\Big{)}^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$
$\displaystyle=(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{G}}}^{\top}\underline{\tilde{\bm{G}}})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$
$\displaystyle=(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{G}}}^{\top}\underline{\tilde{\bm{G}}})^{-1}\left(\begin{smallmatrix}0\\\
\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1}\end{smallmatrix}\right)+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{G}}}^{\top}\underline{\tilde{\bm{G}}})^{-1}\left(\begin{smallmatrix}1\\\
\bar{\bm{X}}^{*}_{r-1}\end{smallmatrix}\right)$
$\displaystyle=(0,(\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1})^{\top})(\underline{\tilde{\bm{G}}}^{\top}\underline{\tilde{\bm{G}}})^{-1}\left(\begin{smallmatrix}0\\\
\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1}\end{smallmatrix}\right)+\frac{1}{1+\epsilon}$
$\displaystyle=(\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1})^{\top}(\underline{\bm{S}}^{*}_{r-1})^{-1}(\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1})+\frac{1}{1+\epsilon}.$
(A8)

Here, the third equality holds since $\left(\begin{smallmatrix}1\\\
\bar{\bm{X}}^{*}_{r-1}\end{smallmatrix}\right)=\frac{1}{1+\epsilon}\underline{\tilde{\bm{G}}}^{\top}\underline{\tilde{\bm{G}}}\bm{e}_{1}$
and the fourth equality holds since
$(\underline{\tilde{\bm{G}}}^{\top}\underline{\tilde{\bm{G}}})^{-1}=\left(\begin{smallmatrix}\underline{\bm{B}}_{11}&\underline{\bm{B}}_{12}\\\
\underline{\bm{B}}_{21}&\underline{\bm{B}}_{22}\\\ \end{smallmatrix}\right)$,
where
$\underline{\bm{B}}^{-1}_{22}=(\frac{1}{\tilde{n}_{r-1}}\underline{\bm{X}}_{r-1}^{\top}\underline{\bm{X}}_{r-1}+\frac{\epsilon}{N}\underline{\bm{X}}_{\text{full}}^{\top}\underline{\bm{X}}_{\text{full}})-(1+\epsilon)\bar{\bm{X}}^{*}_{r-1}\bar{\bm{X}}^{*\top}_{r-1}=\underline{\bm{S}}^{*}_{r-1}.$
This completes the proof. ∎

#### Proof of Theorem 4.2

###### Proof.

(a) We first consider the setting of a standard block design where $N=BG$
(i.e., $c=1$). The blocks are labelled $1,2,...,B$. Here, the SOM is
constructed by stacking $B$ independent random permutations of the ‘chunk’
$(1,2,...,G)$. We will show that the choices made by the treatment groups in
the FSM follow the assignment mechanism of an RBD.

Consider the first randomized chunk of the SOM, i.e., a random permutation of
$(1,2,...,G)$. At the first stage of this randomized chunk, the choosing
treatment group aims to maximize
$(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$. Note that we can write
$\underline{\tilde{\bm{X}}}_{\text{full}}$ as
$\underline{\tilde{\bm{X}}}_{\text{full}}=\left(\begin{smallmatrix}\underline{\bm{D}}\\\
\vdots\\\ \underline{\bm{D}}\end{smallmatrix}\right)$, where
$\underline{\bm{D}}_{B\times B}=\left(\begin{smallmatrix}1&1&0&...&0&0\\\
1&0&1&...&0&0\\\ \vdots&\vdots&\vdots&...&\vdots&\vdots\\\ 1&0&0&...&0&1\\\
1&0&0&...&0&0\\\ \end{smallmatrix}\right)$. Now, consider a transformation of
the rows of the design matrix given by
$\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}=(\underline{\bm{D}}^{\top})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$. The transformed design matrix is
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{\text{full}}=\underline{\tilde{\bm{X}}}_{\text{full}}\underline{\bm{D}}^{-1}=\left(\begin{smallmatrix}\underline{\bm{I}}_{B}\\\
\vdots\\\ \underline{\bm{I}}_{B}\end{smallmatrix}\right)$. We note that the
$\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$s nothing but standard
unit vectors. Now,

$\displaystyle(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$
$\displaystyle=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{\text{full}}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{\text{full}})^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}.$
(A9)

Therefore, the selection function remains the same under the above
transformation. Now,
$\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{\text{full}}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{\text{full}})^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}=\frac{1}{G}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}=\frac{1}{G}$
for all $i$, which essentially implies that the choosing group has no
preference among the units for selection and hence chooses any one of the $N$
units randomly. Similarly, at the subsequent stages of this randomized chunk,
the corresponding choosing groups select one of the remaining units randomly.

Next, we consider the second randomized chunk of the SOM. Without loss of
generality, suppose treatment 1 gets to choose first in this chunk. Also,
without loss of generality, suppose that in its first choice, treatment 1 had
selected a unit from block 1. We claim that in this selection, treatment 1
will choose one of the remaining units randomly from any block other than
block 1, which respects the assignment mechanism of an RBD.

To prove the claim, we first consider the objective function at this stage.
Treatment 1 aims to maximize
$(1,\bm{X}^{\top}_{i})\Big{(}\frac{1}{\tilde{n}_{r-1}}\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}+\frac{\epsilon}{N}\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}}\Big{)}^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$. Here, we denote the current stage by $r$.
Using the same transformation as in the case of the first chunk, we can write
the objective function as
$\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\Big{(}\frac{1}{\tilde{n}_{r-1}}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}+\frac{\epsilon}{N}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{\text{full}}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{\text{full}}\Big{)}^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$,
where
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}=\underline{\tilde{\bm{X}}}_{r-1}\underline{\bm{D}}^{-1}$.
Since
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{\text{full}}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{\text{full}}=G\underline{\bm{I}}_{B}$,
it is equivalent to maximize

$\displaystyle\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\Big{(}\underline{\bm{I}}_{b}+\frac{B}{\tilde{n}_{r-1}\epsilon
G}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}\Big{)}^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$
$\displaystyle=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\Big{(}\underline{\bm{I}}_{b}+\delta\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}\Big{)}^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$
(A10)
$\displaystyle=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\Big{\\{}\underline{\bm{I}}_{b}-\delta\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}(\underline{\bm{I}}_{\tilde{n}_{r-1}}+\delta\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1})^{-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}\Big{\\}}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}.$
(A11)

Here, $\delta=\frac{B}{\tilde{n}_{r-1}\epsilon G}$. The final equality holds
by the Woodbury matrix identity. Now, in this case,
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}=(1,0,...,0)$
(since treatment 1 has only selected one unit from block 1 up to this stage).
So, the objective function in Equation A11 equals
$1-\frac{\delta}{1+\delta}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1})\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}$.
Since $\delta>0$, it is equivalent to minimize
$\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1})\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\left(\begin{smallmatrix}1&\bm{0}^{\top}_{1\times(B-1)}\\\
\bm{0}_{(B-1)\times
1}&\underline{\bm{0}}_{(B-1)\times(B-1)}\end{smallmatrix}\right)\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$,
which takes the value $0$ for a unit in any block other than block 1 and $1$
for a unit in block $1$. This proves the claim for treatment 1. Moreover, by
similar reasoning, the claim holds for all the other treatment groups in this
randomized chunk.

Next, we consider a general randomized chunk of the SOM. Once again, without
loss of generality, suppose treatment 1 gets to choose first in this chunk.
Also, for simplicity of exposition and without loss of generality, suppose
treatment 1 has already selected from blocks $1,2,...,b$, implying that
$\tilde{n}_{r-1}=b$ and
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}=\left(\begin{smallmatrix}\underline{\bm{I}}_{b}&\underline{\bm{0}}_{b\times(B-b)}\end{smallmatrix}\right)$.
This form of
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}$, along
with Equation A11 implies that it is equivalent to minimize
$\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1})\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\left(\begin{smallmatrix}\underline{\bm{I}}_{b}&\underline{\bm{0}}_{b\times(B-b)}\\\
\underline{\bm{0}}^{\top}_{(B-b)\times
b}&\underline{\bm{0}}_{(B-b)\times(B-b)}\end{smallmatrix}\right)\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$,
which is minimized for any unit $i$ belonging to the blocks $b+1,...,B$. This
shows that at this stage, treatment 1 randomly chooses a unit from a block
other than the blocks it has already chosen from. By similar reasoning, at
subsequent stages of this randomized chunk, the choosing group follows the
same selection strategy for their own group. This completes the proof of the
theorem for the setting of a standard block design.

We now prove the theorem for the general block design setting with $N=cBG$,
$c>1$. The proof strategy is exactly the same as the $c=1$ setting. Here the
SOM is generated by randomly permuting the chunk $(1,2,...,G)$ $B\times c$
times. Once the selections are completed for the the first $B$ chunks, the
resulting assignment resembles that of a standard RBD (by the previous proof),
where each treatment group randomly chooses exactly one unit from each block.
For the $(B+1)$th chunk, suppose, without loss of generality, that treatment 1
gets to choose first. At this stage (denoted by stage $r$), treatment 1 tries
to maximize,

$\displaystyle(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)$
$\displaystyle=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1})^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$
$\displaystyle=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}=1,$
(A12)

where the penultimate equality holds since
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}=\underline{\bm{I}}_{B}$.
Thus, similar to the first randomized chunk in the setting of $c=1$, treatment
1 (and the other treatments) randomly chooses one of the available units.

Finally, we consider a general chunk. Without loss of generality, suppose
treatment 1 gets to choose first in this chunk. We can write the corresponding
transformed design matrix
$\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}$ as

$\displaystyle\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1}=\left(\begin{smallmatrix}&\underline{\bm{I}}_{B}\\\
&\underline{\bm{I}}_{B}\\\ &\vdots\\\ &\underline{\bm{I}}_{B}\\\
\underline{\bm{I}}_{b}&&\underline{\bm{0}}_{b\times(B-b)}\end{smallmatrix}\right).$
(A13)

Here, without loss of generality, we have assumed that treatment 1 has chosen
$c_{0}+b$ times from the first $b$ blocks and $c_{0}$ times from the remaining
blocks, where $c_{0}<c$. This implies that treatment 1 aims to maximize.

$\displaystyle\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}(\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}^{\top}_{r-1}\underline{\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}}_{r-1})^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i}$
$\displaystyle=\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}^{\top}_{i}\Big{\\{}c_{0}\underline{\bm{I}}_{B}+\left(\begin{smallmatrix}\underline{\bm{I}}_{b}\\\
\underline{\bm{0}}^{\top}_{(B-b)\times
b}\end{smallmatrix}\right)\left(\begin{smallmatrix}\underline{\bm{I}}_{b}&\underline{\bm{0}}_{b\times(B-b)}\end{smallmatrix}\right)\Big{\\}}^{-1}\tilde{\raisebox{0.0pt}[0.85pt]{$\tilde{\bm{X}}$}}_{i},$
(A14)

which has the same form as the objective function in Equation A10 in the $c=1$
setting. Thus, following similar arguments as in the $c=1$ setting, we
conclude that at this stage, treatment 1 selects a unit randomly from blocks
$b+1,...,B$, which conforms to the assignment mechanism of an RBD. Also, at
subsequent stages of the randomized chunk, the choosing group follows the same
selection strategy for their own group. This completes the proof of the
theorem.

(b) With two groups of equal sizes, the SOM consists of successive random
permutations of the ‘chunk’ $(1,2)$. By Theorem 4.1, for the first pair of
stages of selection, the objective function (to maximize) is given by

$\displaystyle(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}})^{\top}(\underline{\bm{S}}_{\text{full}})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}}).$
(A15)

Under the assumption of identical twins and continuous data generating
distributions, with probability 1, there are exactly two units (one being a
twin of the other), whose common covariate value $\bm{X}^{(1)}$ (say)
maximizes the objective function in Equation A15. Therefore, the choosing
group at the first stage selects one of these two identical twins randomly,
and in the next stage, the other treatment selects the remaining twin. This
respects the assignment mechanism of a matched-pair design.

Consider the next pairs of stages. The objective function of the choosing
treatment group is given by:

$\displaystyle(\bm{X}_{i}-\frac{1}{1+\epsilon}\bm{X}^{(1)})^{\top}\Big{\\{}\bm{X}^{(1)}\bm{X}^{(1)\top}+\frac{\epsilon}{N}\underline{\bm{X}}^{\top}_{\text{full}}\underline{\bm{X}}_{\text{full}}-(1+\epsilon)\bm{X}^{(1)}\bm{X}^{(1)\top}\Big{\\}}^{-1}(\bm{X}_{i}-\frac{1}{1+\epsilon}\bm{X}^{(1)})$
(A16)

Similar to the previous case, here also we have (with probability 1) exactly
two units, one being a twin of the other, whose common covariate value
$\bm{X}^{(2)}$ maximizes the objective function in Equation A16. Thus, the
choosing group at the first stage of this pair selects one of these two twins
randomly, and in the next stage, the other treatment chooses the remaining
twin. Proceeding in this manner, it follows that, at the end of the selection
process, each treatment group ends up selecting one twin randomly from
$\frac{N}{2}$ identical twins, which is equivalent to a matched-pair design.
This completes the proof.

∎

#### Proof of Proposition A3

###### Proof.

The A-optimal selection function aims to minimize

$\displaystyle\text{trace
}\Big{\\{}\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r,i}\underline{\tilde{\bm{X}}}_{r,i})^{-1}\Big{\\}}$
(A17)
$\displaystyle=\text{trace}\Big{[}\underline{\bm{T}}\\{\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}+\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)(1,\bm{X}^{\top}_{i})\\}^{-1}\Big{]}$
$\displaystyle=\text{trace}\Big{\\{}\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}-\underline{\bm{T}}\frac{(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}}{1+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}\Big{\\}}$
$\displaystyle=\text{trace}\\{\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\\}-\text{trace}\Big{\\{}\underline{\bm{T}}\frac{(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}}{1+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}\Big{\\}}$
$\displaystyle=\text{trace}\\{\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\\}-\text{trace}\Big{\\{}\frac{(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}{1+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}\Big{\\}}$ (A18)

∎

Here the second equality holds due to the Sherman-Morrison-Woodbury formula,
the third and fourth equality hold due to the linearity and cyclicality of
$\text{trace}(\cdot)$, respectively. Equation A18 shows that it is equivalent
to maximize
$\text{trace}\Big{\\{}\frac{(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}{1+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}\Big{\\}}$. This completes the proof.

######  Proposition A1.

For $n_{1}=n_{2}=...=n_{G}$, the difference-in-means statistic between
treatment $g^{\prime}$ and $g^{\prime\prime}$ is unbiased for
$\text{SATE}_{g^{\prime},g^{\prime\prime}}$.

###### Proof.

With equal-sized groups, by symmetry, every unit has an equal chance of
belonging to one of the $G$ treatment groups. That is,
$P(Z_{i}=g)=\frac{1}{G}$ for all $g\in\\{1,2,...,G\\}$. Therefore,

$\displaystyle\mathbb{E}\Big{\\{}\frac{1}{n_{g}}\sum_{i:Z_{i}=g}Y^{\text{obs}}_{i}|\bm{Y}(g)\Big{\\}}$
$\displaystyle=\mathbb{E}\Big{\\{}\frac{G}{N}\sum_{i=1}^{N}\mathbbm{1}(Z_{i}=g)Y_{i}(g)|\bm{Y}(g)\Big{\\}}$
$\displaystyle=\frac{G}{N}\sum_{i=1}^{N}P(Z_{i}=g)Y_{i}(g)$
$\displaystyle=\frac{1}{N}\sum_{i=1}^{N}Y_{i}(g).$ (A19)

The proposition follows from Equation A19 applied to treatment groups
$g^{\prime}$ and $g^{\prime\prime}$. ∎

###  C Properties of D-optimal selection function

###  D Affine invariance and covariate balance

######  Theorem A2.

  1. (a)

The FSM with the D-optimal selection function is invariant under affine
transformations of the covariate vector.

  2. (b)

For continuous, symmetrically distributed covariates and two groups of equal
size, the FSM with the D-optimal selection function almost surely produces
exact mean-balance on all even transformations of the centered covariate
vector.

#### Proof of Theorem A2

###### Proof.

(a) We consider the case where $\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
invertible. The proofs for the other two cases are similar. By Theorem 4.1, in
this case, the chosen unit $i^{\prime}$ satisfies,

$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(\bm{X}_{i}-\bar{\bm{X}}_{r-1})^{\top}(\underline{\bm{S}}_{r-1})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{r-1}).$
(A20)

Consider an affine transformation of the covariate $\bm{X}$ given by
$\bm{U}=\underline{\bm{A}}\bm{X}+\bm{b}$, where $\underline{\bm{A}}$ is a
$k\times k$ invertible matrix and $\bm{b}$ is a vector of dimension k. Let the
corresponding values of $\bar{\bm{X}}_{r-1}$ and $\underline{\bm{S}}_{r-1}$ be
$\bar{\bm{U}}_{r-1}$ and $\underline{\bm{S}}_{U,r-1}$, respectively. We
observe that,

$\displaystyle(\bm{U}_{i}-\bar{\bm{U}}_{r-1})^{\top}(\underline{\bm{S}}_{U,r-1})^{-1}(\bm{U}_{i}-\bar{\bm{U}}_{r-1})$
$\displaystyle=\\{\underline{\bm{A}}(\bm{X}_{i}-\bar{\bm{X}}_{r-1})\\}^{\top}(\underline{\bm{A}}\underline{\bm{S}}_{r-1}\underline{\bm{A}}^{\top})^{-1}\underline{\bm{A}}(\bm{X}_{i}-\bar{\bm{X}}_{r-1})$
$\displaystyle=(\bm{X}_{i}-\bar{\bm{X}}_{r-1})^{\top}(\underline{\bm{S}}_{r-1})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{r-1}).$
(A21)

This shows that the D-optimal selection function remains unchanged under
affine transformations and hence, FSM with the D-optimal selection function is
affine invariant.

(b) The in-sample symmetry of the data essentially implies that if $\bm{X}$
belongs to the sample, then $-\bm{X}$ also belongs to the sample. Moreover, by
the assumption of a continuous data generating distribution, with probability
1, the covariate values are different up to reflection. Now, consider an even
transformation $g(\cdot)$, i.e., $g(-\bm{X})=g(\bm{X})$. With two groups of
equal sizes, the SOM consists of successive random permutations of the ‘chunk’
$(1,2)$. By Theorem 4.1, for the first pair of stages of selection, the
objective function (to maximize) is given by

$\displaystyle(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}})^{\top}(\underline{\bm{S}}_{\text{full}})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{\text{full}})$
$\displaystyle=\bm{X}^{\top}_{i}(\underline{\bm{S}}_{\text{full}})^{-1}\bm{X}_{i}.$
(A22)

It follows that, if a unit in the sample with covariate $\bm{X}^{(1)}$
maximizes the objective function in Equation A22, then so does the unit with
covariate $-\bm{X}^{(1)}$. Moreover, due to the continuous data generating
distribution, with probability 1, these are the only two units that maximize
this objective function. Therefore, if treatment 1 selects the unit with
covariate $\bm{X}^{(1)}$, treatment 2 selects the unit with covariate
$-\bm{X}^{(1)}$, and vice-versa. This preserves exact balance on $g(\bm{X})$.

Now, consider the next pair of stages. Without loss of generality, suppose
treatment 1 had chosen a unit with covariate $\bm{X}^{(1)}$ and treatment 2
had chosen a unit with covariate $-\bm{X}^{(1)}$ in their respective previous
choices. Also, without loss of generality, assume that in this pair of stages,
treatment 1 gets to choose first. By Theorem 4.1, treatment 1 aims to
maximize,

$\displaystyle(\bm{X}_{i}-\bar{\bm{X}}^{*}_{r-1})^{\top}(\underline{\bm{S}}^{*}_{r-1})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{r-1}^{*})$
$\displaystyle=\\{\bm{X}_{i}-\frac{1}{1+\epsilon}\bm{X}^{(1)}\\}^{\top}\Big{\\{}\bm{X}^{(1)}\bm{X}^{(1)\top}+\frac{\epsilon}{N}\underline{\bm{X}}^{\top}_{\text{full}}\underline{\bm{X}}_{\text{full}}-(1+\epsilon)\bm{X}^{(1)}\bm{X}^{(1)\top}\Big{\\}}^{-1}\\{\bm{X}_{i}-\frac{1}{1+\epsilon}\bm{X}^{(1)}\\}.$
(A23)

Also, during treatment 2’s turn in this pair of stages, it tries to maximize

$\displaystyle(\bm{X}_{i}+\frac{1}{1+\epsilon}\bm{X}^{(1)})^{\top}\Big{\\{}\bm{X}^{(1)}\bm{X}^{(1)\top}+\frac{\epsilon}{N}\underline{\bm{X}}^{\top}_{\text{full}}\underline{\bm{X}}_{\text{full}}-(1+\epsilon)\bm{X}^{(1)}\bm{X}^{(1)\top}\Big{\\}}^{-1}(\bm{X}_{i}+\frac{1}{1+\epsilon}\bm{X}^{(1)})$
$\displaystyle=\\{(-\bm{X}_{i})-\frac{1}{1+\epsilon}\bm{X}^{(1)}\\}^{\top}\Big{\\{}\bm{X}^{(1)}\bm{X}^{(1)\top}+\frac{\epsilon}{N}\underline{\bm{X}}^{\top}_{\text{full}}\underline{\bm{X}}_{\text{full}}-(1+\epsilon)\bm{X}^{(1)}\bm{X}^{(1)\top}\Big{\\}}^{-1}\\{(-\bm{X}_{i})-\frac{1}{1+\epsilon}\bm{X}^{(1)}\\}.$
(A24)

Equations A23 and A24 imply that if treatment 1 chooses a unit with covariate
value $\bm{X}^{(2)}$, then with probability 1, treatment 2 chooses the unit
with covariate value $-\bm{X}^{(2)}$, and vice versa. This shows that, at the
end of the second pair of stages in the SOM, exact mean balance on $g(\bm{X})$
is preserved. Proceeding in this manner it follows that, at the end of the
selection process, with probability 1, both the treatment groups will have
exact balance on $g(\bm{X})$. This completes the proof.

∎

It follows from Theorem A2(a) that, for any SOM, the choices made by each
treatment group remain unchanged even if the covariate vectors are transformed
via an affine transformation (e.g., changing the units of measurement of the
covariates). Therefore, the FSM with the D-optimal selection function self-
standardizes the covariates. In addition, if the covariate vector is
symmetrically distributed in the sample, then by Theorem A2(b), the FSM
exactly balances even transformations such as the second, fourth order
moments, and the pairwise products of the centerd covariates. An implication
of Theorem A2(b) is that, for covariates drawn from symmetric continuous
distributions (such as the Normal, t, and Laplace distributions), the FSM
tends to balance all these transformations due to the approximate symmetry of
the covariates in the sample. The choice of the D-optimal selection function
is thus robust in the sense that it allows the FSM to balance a family of
transformations of the covariate vector by design, without explicitly
including them in the assumed linear model nor requiring the specification of
tuning parameters.

####  D.1 Connection to A-optimality

The original FSM used a criterion based on A-optimality as the selection
function (see Morris 1979). In this section, we compare the A-and D-optimal
selection functions. The A-optimal selection function requires prespecifying a
policy matrix $\underline{\bm{P}}_{p\times(k+1)}$ and a corresponding vector
of policy weights $\bm{w}_{p\times 1}$. Here, $\underline{\bm{P}}$ transforms
the original vector of regression coefficients to a vector of $p$ linear
combinations that are of policy interest, and $\bm{w}$ assigns weights to each
combination according to their importance. If treatment 1 gets to choose at
the $r$th stage, then this criterion selects the unit that minimizes the
resulting
$\text{trace}\Big{\\{}\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r,i}\underline{\tilde{\bm{X}}}_{r,i})^{-1}\Big{\\}}$,
where
$\underline{\bm{T}}=\underline{\bm{P}}^{\top}\text{diag}(\bm{w})\underline{\bm{P}}$.
Proposition A3 shows an equivalent characterization of the A-optimal selection
function.

######  Proposition A3.

Let treatment 1 be the choosing group at the $r$th stage. Assume that
$\tilde{n}_{r-1}\geq 1$ and
$\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1}$ is
invertible. The A-optimal selection function chooses unit $i^{\prime}$ with
covariate vector $\bm{X}_{i^{\prime}}\in\mathbb{R}^{k}$, where
$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}\frac{(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\underline{\bm{T}}(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}{1+(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)}.$

The A-optimality criterion provides a family of selection functions depending
on $\underline{\bm{P}}$ and $\bm{w}$. For some choices of $\underline{\bm{P}}$
and $\bm{w}$, the selection function is not affine invariant , e.g.,
$\underline{\bm{P}}=\bm{I}$ and $\bm{w}=(1,1,...,1)^{\top}$, while for other
choices it is, e.g.,
$\underline{\bm{P}}=\underline{\tilde{\bm{X}}}_{\text{full}}$ and
$\bm{w}=(1,1,...,1)^{\top}$. In particular, the A-optimal selection function
with $\underline{\bm{P}}=\underline{\tilde{\bm{X}}}_{\text{full}}$ and
$\bm{w}=(1,1,...,1)^{\top}$ is closely related to the D-optimal selection
function. To see this, consider a case where in the selection process, the
design matrices in each treatment group scale similarly relative to the design
matrix in the full sample, i.e.,
$(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}=c_{r}(\underline{\tilde{\bm{X}}}^{\top}_{\text{full}}\underline{\tilde{\bm{X}}}_{\text{full}})^{-1}$
for some constant $c_{r}>0$. In this case, the A-optimal selection function
chooses unit $i^{\prime}$ such that
$i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(1,\bm{X}^{\top}_{i})(\underline{\tilde{\bm{X}}}^{\top}_{r-1}\underline{\tilde{\bm{X}}}_{r-1})^{-1}\left(\begin{smallmatrix}1\\\
\bm{X}_{i}\end{smallmatrix}\right)\iff
i^{\prime}\in\arg\\!\max\limits_{i\in\mathcal{R}_{r-1}}(\bm{X}_{i}-\bar{\bm{X}}_{r-1})^{\top}(\underline{\bm{S}}_{r-1})^{-1}(\bm{X}_{i}-\bar{\bm{X}}_{r-1})$,
which is equivalent to the D-optimal selection function. Hence, in this case,
the FSM under the D-optimal and A-optimal selection functions make similar
choices of units.

###  E Algorithms for constructing an SOM

####  E.1 The SCOMARS algorithm

Consider a setting with $G=2$ treatment groups of arbitrary sizes $n_{1}$ and
$n_{2}$. Let $W_{r}$ be the binary indicator for selection of group 1 stage
$r$, $r\in\\{1,2,…,N\\}$, with $p_{r}:=P(W_{r}=1)$ being the marginal
probability of selection at stage $r$. Write $S_{r}:=\sum_{j=1}^{r}W_{j}$ and
$F_{r}:=\mathbb{E}(S_{r})=\sum_{j=1}^{r}p_{j}$. A treatment assignment is
sequentially controlled if $|S_{r}-F_{r}|<1$ for all $r\in\\{1,2,…,N\\}$.

The SCOMARS algorithm proceeds as follows:

  * •

Stage 1, $P(W_{1}=1)=p_{1}$.

  * •

Stage $r\geq 2$,
$P(W_{r}=1|S_{r-1}=s_{r-1})=P\Big{\\{}U\leq\frac{p_{r}-\max(0,s_{r-1}-F_{r-1})}{1-|s_{r-1}-F_{r-1}|}\Big{\\}}$,
where $U\sim\text{Unif}(0,1)$.

This algorithm satisfies the sequentially controlled condition,
$|S_{r}-F_{r}|<1$ for all $r\in\\{1,2,…,N\\}$ (Morris 1983). It is Markovian
because the probability of selection at stage $r$ depends solely on stage
$r-1$.

####  E.2 SOM for multi-group experiments

We first define the randomized chunk algorithm for generating an SOM for
multi-group experiments with equal group sizes.

######  Definition 1 (Randomized chunk algorithm).

Suppose $n_{1}=n_{2}=...=n_{G}$. The randomized chunk algorithm generates an
SOM by generating and stacking $\frac{N}{G}$ independent random permutations
of the ‘chunk’ $(1,2,...,G)$.

For example, with $N=12$, $g=3$, $n_{1}=n_{2}=n_{3}=4$, one instance of an SOM
generated using randomized chunk is
$(\underbrace{2,1,3},\underbrace{1,2,3},\underbrace{2,1,3},\underbrace{2,3,1})^{\top}$.

The following proposition shows that the randomized chunk algorithm is
sequentially controlled.

######  Proposition A4.

For $G\geq 2$ and $n_{1}=n_{2}=...=n_{G}$, the randomized chunk algorithm
satisfies $|S_{ig}-F_{ig}|\leq\frac{G-1}{G}<1$ for all $g\in\\{1,2,...,G\\}$.

###### Proof.

Let $S_{ig}$ and $F_{ig}$ be the same as defined in Section 8.1
($i\in\\{1,2,...,N\\}$, $g\in\\{1,2,...,G\\}$). For equal sized treatment
groups, $F_{ig}=\frac{i}{G}$. Now, without loss of generality, it suffices to
show that $|S_{i1}-F_{i1}|\leq\frac{G-1}{G}$ for all $i\in\\{1,2,...,N\\}$.
Consider the first chunk in the SOM, which is a random permutation of
$(1,2,...,G)$. If treatment 1 appears in position $i^{*}\in\\{1,2,...,G\\}$
the permutation ($j\in\\{1,2,...,G\\}$), then

$|S_{i1}-F_{i1}|=\begin{cases}\frac{i}{G}&\text{if
$i\in\\{1,...,i^{*}-1\\}$}\\\ 1-\frac{i}{G}&\text{if
$i\in\\{i^{*},...,G\\}$.}\end{cases}$ (A25)

In each case, $|S_{i1}-F_{i1}|\leq\frac{G-1}{G}$ for all
$i\in\\{1,2,...,G\\}$. Moreover, since $|S_{G1}-F_{G1}|=0$, the SOM restarts
itself after the first chunk. Hence, we can conclude that
$|S_{i1}-F_{i1}|\leq\frac{G-1}{G}$ for all $i\in\\{1,2,...,G\\}$. This
completes the proof.

∎

Below we describe two algorithms to generate an SOM for multi-group
experiments and show that they are sequentially controlled. The key idea in
these algorithms is the formation of ‘supergroups’, i.e., combination of one
or more treatment groups. For example, with $g=3$,
$n_{1}=10,n_{2}=20,n_{3}=30$, one can consider two supergroups, namely
$\\{1,2\\}$ of size $10+20=30$ and $\\{3\\}$ of size $30$.

######  Theorem A5.

For $1\leq G_{1}\leq G-1$, let $n_{1}=n_{2}=...=n_{G_{1}}\neq n^{(1)}$, and
$n_{G_{1}+1}=n_{G_{1}+2}=...=n_{G}=n^{(2)}$, where $n^{(1)}\neq n^{(2)}$.
Consider the following three-stage algorithm.

  1. 1.

Run SCOMARS with supergroups $\\{1,...,G_{1}\\}$ and $\\{G_{1}+1,...,G\\}$ to
generate an SOM at the supergroup level.

  2. 2.

Consider the locations of the SOM in step 1 where supergroup
$\\{1,...,G_{1}\\}$ chooses. Then, use randomized chunk to obtain the
selection orders at the levels of the original groups in those locations.

  3. 3.

Consider the locations of the SOM in step 1 where supergroup
$\\{G_{1}+1,...,G\\}$ chooses. Then, use randomized chunk to obtain the
selection orders at the levels of the original groups in those locations.

The above SOM generating algorithm is sequentially controlled.

We first prove a special case of Theorem A5, given below in Lemma A2

######  Lemma A2.

The algorithm in Theorem A5 is sequntially controlled for the special case of
$G_{1}=1$.

###### Proof.

The first step of the algorithm in Theorem A5 runs SCOMARS with treatment
group 1 and the supergroup $\\{2,3,...,G\\}$. Thus, the first step itself
determines the locations of the SOM where treatment 1 gets to choose. Since
SCOMARS is sequentially controlled, we immediately have $|S_{i1}-F_{i1}|<1$
for all $i\in\\{1,2,...,N\\}$.

It remains to show that for $g\in\\{2,3,...,G\\}$, $|S_{ig}-F_{ig}|<1$ for all
$i\in\\{1,2,...,N\\}$. By symmetry, it suffices to show this for $g=2$. Now,
the randomized chunk algorithm on the supergroup $\\{2,3,...,G\\}$ determines
the locations of the SOM where treatment 2 gets to choose. We will prove the
result by first mapping this SOM to an SOM where treatment 1 is absent, and
then by using the sequential controlled property of randomized chunk.

Let us first denote $1\leq r_{1}<r_{2}<...<r_{n_{1}-1}<r_{n_{1}}\leq N$ as the
stages or locations of the SOM where treatment 1 gets to choose. We consider
the following cases,

(i) Case-1: $i\in\\{1,2,...,r_{1}-1\\}$. In this case, by stage $i$, treatment
1 has not made any choices. Now,

$\displaystyle|S_{i2}-F_{i2}|$ $\displaystyle=|S_{i2}-\frac{in^{(2)}}{N}|$
$\displaystyle\leq|S_{i2}-\frac{i}{G-1}|+|\frac{i}{G-1}-\frac{in^{(2)}}{N}|$
$\displaystyle\leq\frac{G-2}{G-1}+i\frac{n_{1}}{N(G-1)}$
$\displaystyle<\frac{G-2}{G-1}+\frac{1}{G-1}=1.$ (A26)

Here the first inequality holds due to triangle inequality. To see that second
inequality, consider a new experiment with treatment groups $\\{2,...,G\\}$ of
size $n^{(2)}$ each and an SOM generated by randomized chunk as in the second
step of the algorithm in Theorem A5. Let $\tilde{S}_{i2}$ be the number of
selections made by treatment $2$ up to stage $i$ in this new experiment and
$\tilde{F}_{i2}=\frac{i}{G-1}$ be its expectation. By Proposition A4,
$|\tilde{S}_{i2}-\tilde{F}_{i2}|\leq\frac{G-2}{G-1}$. Now,
$|S_{i1}-\frac{i}{G-1}|=|\tilde{S}_{i1}-\frac{i}{G-1}|$, which gives us the
second inequality. Finally, the last inequality holds since
$\frac{in_{1}}{N}=F_{i1}<1$.

(ii) Case-2: $i\in\\{r_{t},r_{t}+1,....,r_{t+1}-1\\}$ for some
$t\in\\{1,2,...,n_{1}-1\\}$. In this case, by stage $i$, treatment 1 has made
exactly $t$ choices. Now,

$\displaystyle|S_{i2}-F_{i2}|$ $\displaystyle=|S_{i2}-\frac{in^{(2)}}{N}|$
$\displaystyle\leq|S_{i2}-\frac{i-t}{G-1}|+|\frac{i-t}{G-1}-\frac{in^{(2)}}{N}|$
$\displaystyle\leq\frac{G-2}{G-1}+\frac{1}{G-1}|t-\frac{in_{1}}{N}|$
$\displaystyle<\frac{G-2}{G-1}+\frac{1}{G-1}=1.$ (A27)

Here, the first inequality is due to triangle inequality. To see the second
inequality, we again consider the new experiment described in Case-1. Notice
that,
$|S_{i2}-\frac{i-t}{G-1}|=|\tilde{S}_{(i-t)2}-\tilde{F}_{(i-t)2}|\leq\frac{G-2}{G-1}$,
where the last inequality holds by Proposition A4. Finally, the final
inequality in Equation A27 holds since
$|t-\frac{in_{1}}{N}|=|S_{i1}-F_{i1}|<1$. This completes the proof of the
lemma. ∎

We now prove Theorem A5.

###### Proof.

We first show that, for $g\in\\{1,2,...,G_{1}\\}$,

$|S_{ig}-F_{ig}|<1\quad\forall i\in\\{1,2,...,N\\}.$ (A28)

To show this, we consider steps 1 and 2 of the algorithm as these two steps
are sufficient to determine the location of treatments $1,...,G_{1}$ in the
SOM. We note that, steps 1 and 2 generate an SOM for an experiment with
$G_{1}+1$ treatment groups, namely supergroup $\\{G_{1}+1,...,G\\}$ (of size
$(G-G_{1})n^{(2)}$) and groups $1,2,...,G_{1}$ (each of size $n^{(1)}$). Thus,
by Lemma A2, it follows that Equation A28 holds for $g\in\\{1,2,...,G_{1}\\}$.

To show that Equation A28 holds for $g\in\\{G_{1}+1,...,G\\}$, we first notice
that steps 2 and 3 of the algorithm are completely independent and hence can
be performed in any order. Therefore, by changing the order of steps 2 and 3
and applying the same argument as before, we get that Equation A28 holds for
$g\in\\{G_{1}+1,...,G\\}$. This completes the proof of the theorem. ∎

######  Theorem A6.

Let $G_{1},...,G_{m}$ be such that $1\leq G_{j}\leq G-1$ for all
$j\in\\{1,2,...,m\\}$ and $G_{1}+G_{2}+...+G_{m}=G$. Moreover, for
$j\in\\{1,2,...,m\\}$, let $n^{(j)}$ be the group size of $G_{j}$ many
treatment groups, with $n^{(1)}G_{1}=n^{(2)}G_{2}=...=n^{(m)}G_{m}$. Denote
the collection of $G_{j}$ treatment groups with group sizes $n^{(j)}$ as
supergroup $\mathcal{G}_{j}$. Consider the following multi-stage algorithm.

  1. 1.

Run randomized chunk on supergroups
$\mathcal{G}_{1},\mathcal{G}_{2},...,\mathcal{G}_{m}$ to generate an SOM at
the supergroup level.

  2. 2.

For $j\in\\{1,2,...,m\\}$, consider the locations of the SOM in step 1 where
supergroup $\mathcal{G}_{j}$ chooses. Then, use randomized chunk to obtain the
selection orders at the levels of the original groups in those locations.

The above SOM generating algorithm is sequentially controlled.

To prove this theorem, we first use the following Lemma.

######  Lemma A3.

Let $n_{1}=n_{2}=...=n_{G}=n$. Consider the following SOM generating
algorithm.

  1. 1.

Consider the supergroups $\\{1\\}$ (of size $n$) and $\\{2,3,...,G\\}$ (of
size $(G-1)n$). Generate an SOM at the superpopulation level using SCOMARS.

  2. 2.

Consider the locations of the SOM in step 1 where supergroup $\\{2,3,...,G\\}$
chooses. Then, use randomized chunk to obtain the selection orders at the
levels of the original groups in those locations.

This algorithm is equivalent to the randomized chunk algorithm.

Below we prove this lemma.

###### Proof.

To show that the algorithm is equivalent to randomized chunk, we have to show
that it generates a random permutation of $(1,2,...,G)$ for the first $G$
stages, a fresh random permutation of $(1,2,...,G)$ for the next $G$ stages,
and so on. Since the locations of groups $\\{2,...,G\\}$ are chosen using
randomized chunk, it thus suffices to show that, treatment 1 gets to choose
once (in a random location) in the first $G$ stages, once in the next $G$
stages, and so on.

We use the notation as in Section E.1. Now, suppose among the first $G$
stages, treatment 1 gets to choose at stage $r^{*}$ first. Notice that $r^{*}$
cannot be greater than $G$ as

$\displaystyle
P(W_{G}=1|S_{G-1}=0)=P\Big{\\{}U\leq\frac{\frac{1}{G}-\max(0,0-F_{G-1})}{1-|0-F_{G-1}|}\Big{\\}}=P\Big{\\{}U\leq\frac{1}{G-(G-1)}\Big{\\}}=1.$
(A29)

Now, for $r\in\\{1,2,...,r^{*}-1\\}$ we have,

$\displaystyle
P(W_{r}=1|S_{r-1}=0)=P\Big{\\{}U\leq\frac{\frac{1}{G}-\max(0,0-F_{r-1})}{1-|0-F_{r-1}|}\Big{\\}}=P\Big{\\{}U\leq\frac{1}{G-(r-1)}\Big{\\}}=\frac{1}{G-(r-1)}.$
(A30)

For $r^{*}+1\leq r\leq G$,

$\displaystyle
P(W_{r}=1|S_{r-1}=1)=P\Big{\\{}U\leq\frac{p_{r}-\max(0,1-F_{r-1})}{1-|1-F_{r-1}|}\Big{\\}}=P\Big{\\{}U\leq\frac{\frac{1}{G}-1+\frac{r-1}{G}}{\frac{r-1}{G}}\Big{\\}}=0.$
(A31)

Finally,

$\displaystyle
P(W_{G+1}=1|S_{G}=1)=P\Big{\\{}U\leq\frac{p_{G+1}-\max(0,1-F_{G})}{1-|1-F_{G}|}\Big{\\}}=P\Big{(}U\leq\frac{1}{G}\Big{)}=\frac{1}{G}.$
(A32)

Therefore, by Equation A31, if treatment 1 selects at the $r^{*}$th stage, it
never selects again $2,3,...,G$. Also, by Equation A30, before the $r^{*}$th
stage, the conditional probabilities of treatment 1 selecting are same as what
it would have been under random permutation of the group labels. Finally, by
Equation A32, the process restarts itself at the $(G+1)$th stage, which is
equivalent to starting a fresh new random permutation of the group labels.
This completes the proof of the lemma. ∎

We now prove Theorem A6.

###### Proof.

By the symmetry of the problem, it suffices to show that $|S_{i1}-F_{i1}|<1$
for all $i\in\\{1,2,...,N\\}$. Without loss of generality, we assume that
$\mathcal{G}_{1}=\\{1,2,...,G_{1}\\}$, which implies that treatment 1 belongs
to supergroup $\mathcal{G}_{1}$. Now, it suffices to focus on the following to
steps of the algorithm:

  1. 1.

Run randomized chunk on supergroups
$\mathcal{G}_{1},\mathcal{G}_{2},...,\mathcal{G}_{m}$ to generate an SOM at
the supergroup level.

  2. 2.

Consider the locations of the SOM in step 1 where supergroup $\mathcal{G}_{1}$
chooses. Then, use randomized chunk to obtain the selection orders at the
levels of the original groups in those locations.

This is because, these two steps completely determine the locations of
treatment 1 in the SOM. By Lemma A3, these two steps can be equivalently
performed as follows.

  1. 1.

Consider the supergroups $\mathcal{G}_{1}$ (of size $n^{(1)}G_{1}$) and
$\\{\mathcal{G}_{2},...,\mathcal{G}_{m}\\}$ (of size $(m-1)n^{(1)}G_{1}$).
Generate an SOM at this supergroup level using SCOMARS.

  2. 2.

Consider the locations of the SOM in step 1 where supergroup
$\\{\mathcal{G}_{2},...,\mathcal{G}_{m}\\}$ chooses. Then, use randomized
chunk to obtain the selection orders at the levels of $\mathcal{G}_{j}$ in
those locations.

  3. 3.

Consider the locations of the SOM in step 1 where supergroup $\mathcal{G}_{1}$
chooses. Then, use randomized chunk to obtain the selection orders at the
levels of the original groups in those locations.

We note that this above algorithm is exactly equivalent to the SOM generating
algorithm in Theorem A5 for an experiment with $G_{1}+m-1$ treatment groups,
namely, $1,2,...,G_{1},\mathcal{G}_{2},\mathcal{G}_{3},...,\mathcal{G}_{m}$.
Thus, by Theorem A5, we have $|S_{i1}-F_{i1}|<1$ for all
$i\in\\{1,2,...,N\\}$. ∎

###  F FSM for stratified experiments

In this section, we discuss two potential approaches to use an FSM for
stratified experiments. We consider stratified experiments where the treatment
group sizes within each stratum are set by the investigator beforehand. To
accommodate the FSM to such experiments, we again need to carefully construct
an SOM. In particular, we append the SOM with an additional column of stratum
labels, indicating which stratum the treatment group selects from at each
stage of the selection process. This column of stratum labels is specified in
such a way that the resulting SOM satisfies the group size requirements within
each stratum.

Conceptually, the most straightforward approach is to generate a separate SOM
for each stratum. This is equivalent to setting the column of stratum labels
as
$(\underbrace{1,...,1}_{m_{1}},\underbrace{2,...,2}_{m_{2}},...,\underbrace{S,....,S}_{m_{S}})^{\top}$,
where $S$ is the number of strata and $m_{s}$ is the size of $s$th stratum,
$s\in\\{1,2,...,S\\}$. This approach is easy to implement and can be useful
if, e.g., data on each stratum is available at different stages of the
experiment, akin to a sequential experiment. However, in this approach, the
treatment groups only get to explore the covariate space of a single stratum
for a number of successive stages of selection and hence may not make the most
efficient choices. We address this issue with an alternative approach. For
ease of exposition, we consider two strata: 1 and 2. Let $n_{1g}$ and $n_{2g}$
be the (fixed) sizes of treatment group $g\in\\{1,2,...,G\\}$ in strata 1 and
2, respectively, where $n_{1g}+n_{2g}=n_{g}$. In this approach, we first
generate a usual SOM with group sizes $n_{1},...,n_{G}$. For
$g\in\\{1,2,...,G\\}$, we then select the order of the strata that treatment
$g$ chooses from by running a SCOMARS algorithm with group sizes $n_{1g}$ and
$n_{2g}$. By allowing the treatment groups to select units from different
strata in a balanced manner, this approach mimics the unstratified FSM where
the covariate space of the entire sample is explored for choosing units. Also,
by design, this approach satisfies the size requirement of each treatment
group within each stratum.

###  G FSM for sequential experiments

In this section, we describe our approach to using the FSM for sequential
experiments. Suppose treatment 1 gets to choose at the first stage of
selection for the new batch. Let $\underline{\tilde{\bm{X}}}_{\text{old}}$ be
the design matrix based on units already assigned to treatment 1. Also, for
each unit $i$ in the new batch, let
$\underline{\tilde{\bm{X}}}_{\text{new},i}:=\left(\begin{smallmatrix}\underline{\tilde{\bm{X}}}_{\text{old}}\\\
(1,\bm{X}^{\top}_{i})\end{smallmatrix}\right)$ be the resulting design matrix
in treatment group 1 if unit $i$ is selected. Treatment 1 selects the unit
that maximizes
$\det(\underline{\tilde{\bm{X}}}^{\top}_{\text{new},i}\underline{\tilde{\bm{X}}}_{\text{new},i})$.
In other words, we use the design matrix based on all the units already
assigned to the choosing treatment group to evaluate the D-optimal selection
function for each unit in the new batch, and select the unit that maximizes
the selection function. By carrying over the existing design matrix to the new
batch, this approach tends to correct for any existing covariate imbalances.

###  H A simulation study

####  H.1 Setup

We now compare the performance of the FSM to complete randomization and
rerandomization in a simulation study. Here, $N=120$, $G=2$, $n_{1}=n_{2}=60$,
and $k=6$. The covariates are generated following the design of Hainmueller
(2012):

$\left(\begin{smallmatrix}X_{1}\\\ X_{2}\\\
X_{3}\end{smallmatrix}\right)\sim\mathcal{N}_{3}\left\\{\left(\begin{smallmatrix}0\\\
0\\\ 0\end{smallmatrix}\right),\left(\begin{smallmatrix}2&1&-1\\\ 1&1&-0.5\\\
-1&-0.5&1\end{smallmatrix}\right)\right\\},\hskip
2.84544ptX_{4}\sim\text{Unif}(-3,3),\hskip
2.84544ptX_{5}\sim\chi^{2}_{1},\hskip
2.84544ptX_{6}\sim\text{Bernoulli}(0.5).$ (A33)

In this design, $X_{4}$, $X_{5}$, and $X_{6}$ are mutually independent and
separately independent of $(X_{1},X_{2},X_{3})^{\top}$. We draw a sample of
120 units once from the data generating mechanism in (A33). Conditional on
this sample, we compare four different assignment methods, namely a completely
randomized design (CRD), rerandomization with 0.01 acceptance rate (RR 0.01),
rerandomization with 0.001 acceptance rate (RR 0.001), and the FSM. Both RR
0.01 and RR 0.001 use as rerandomization criteria the Mahalanobis distance
between the two treatment groups on the original covariates. The FSM uses a
linear potential outcome model on the original covariates and the D-optimal
selection function. For each design we draw 800 independent assignments. The