{'fraction_non_alphanumeric': 0.03885260808464908, 'fraction_numerical': 0.016202364222534512, 'mean_word_length': 4.4175250037319005, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Online Balanced Experimental Design

e consider the experimental design problem in an online environment, an important practical task for reducing the variance of estimates in randomized experiments which allows for greater precision, and in turn, improved decision making. In this work, we present algorithms that build on recent advances in online discrepancy minimization which accommodate both arbitrary treatment probabilities and multiple treatments. The proposed algorithms are computational efficient, minimize covariate imbalance, and include randomization which enables robustness to misspecification. We provide worst case bounds on the expected mean squared error of the causal estimate and show that the proposed estimator is no worse than an implicit ridge regression, which are within a logarithmic factor of the best known results for offline experimental design. We conclude with a detailed simulation study showing favorable results relative to complete randomization as well as to offline methods for experimental design with time complexities exceeding our algorithm.


Introduction
Randomized experimentation is a fundamental tool for obtaining counterfactual estimates. The efficacy of randomization comes from a very simple intuition-by randomly assigning treatment status dependence between observed (and unobserved) treatment and pre-treatment covariates necessarily tends to zero as a function of the number of units. In the context of experimentation, this independence condition on observed covariates, commonly known as balance [Imai et al., 2008, Imai andRatkovic, 2014], reduces the variance of estimates of the average treatment effect [Greevy et al., 2004, Higgins et al., 2016, Kallus, 2018, Li et al., 2018, Harshaw et al., 2020. Under the appropriate conditions such corrections can result in large increases in effective sample size, allowing for the detection of the small effects which are commonplace in many large scale studies and industrial applications [Dimmery et al., 2019, Azevedo et al., 2020. These contexts rely heavily on experimentation for decision-making, so reduced variance directly translates into more reliable decisions [Kohavi et al., 2012].
Traditional experimental design like blocking [Greevy et al., 2004, Higgins et al., 2016 or even the novel Gram-Schmidt Walk design [Harshaw et al., 2020] require more than one pass over the sample and their sample complexity is greater than O(n). Even algorithms which admit sequential assignment such as Moore and Moore [2017] suffer from the fact that the algorithm is not linear time and, thus, respondents late in the experiment may take substantially longer to receive a treatment assignment [Cavaille, 2018]. Our work is motivated by this setting to provide a lineartime, single-pass (i.e. sequential) algorithm for balancing experimental design. Our focus is on linear measures of balance (often of particular interest to applied researchers). This provides a new avenue through which experimenters can ensure that their experiments optimize the information they gain from costly samples.
We start by presenting four desiderata for effective, practical online experimental design. First, a method must be computationally efficient. In a review of existing methods for online assignment in the case of survey experiments, Cavaille [2018] finds that existing methods become slow to unusable as increasingly more respondents are included in a study. This resulting speed is fundamentally incompatible with effective administration of an experiment. In short, high latency will cause disproportionately large dropoff in an experiment, which may completely nullify the gains from arXiv:2203.02025v1 [stat.ME] 3 Mar 2022 using more sophisticated experimental design. Any algorithm with greater than linear time complexity exacerbate this problem: higher latency for later subjects than earlier subjects will tend to cause non-random sample attrition, as later subjects (who may be different than earlier respondents) will be more likely to drop out.
Second, an experimental design must reduce covariate imbalance to be effective. This is the entire justification for using methods more sophisticated than Bernoulli randomization, so if the design is not able to improve on balance, then there will be no subsequent reduction in variance and therefore no compelling reason to use it.
Third, the design must incorporate randomization. Harshaw et al. [2020] provides an extensive discussion on the inherent tradeoff between robustness and balance within experimental design. A design which solely optimizes for balance will tend to operate on a knife-edge of accidental bias [Efron, 1971]: the potential bias from an adversarially chosen confounder. With higher accidental bias than Bernoulli randomization, if units do not arrive precisely i.i.d., then the entire design may be compromised. Given that experimental settings are prized specifically for their unbiasedness, this could completely undermine any gains from improved balance. Strong theoretical guarantees on robustness are extremely important in this setting in order to ensure that inferences rest squarely on the design. If assumptions about sampling procedures or data generating processes are necessary to ensure the reliability of estimation, then inferences are not based solely on properties of the experiment, but rather on factors outside the control of the experimenter [Aronow et al., 2021].
Fourth, units which do not show up in the sample should not be included in the balancing. An offline algorithm used in an online environment would fail this condition, because to use it would require to balance the entire population of units expected to show up to an experiment. Units within that population who fail to show would nevertheless be assigned a treatment. Depending on the distribution of units who actually show up in the sample, there are no longer any guarantees that balance will obtain.
Our approach satisfies all four of these desiderata. We propose an online method for covariate balancing requiring linear time and space which provably provides variance reduction. Building upon recent work on discrepancy minimization-the self-balancing walk [Alweiss et al., 2020] and kernel thinning [Dwivedi and Mackey, 2021]-we provide an algorithm whose L 2 -imbalance matches the best known online discrepancy algorithm. Where δ is a failure probability, performing this optimization online results in a log(n/δ) cost in convergence of the average treatment effect over the offline algorithm of discrepancy minimization by Harshaw et al. [2020]. A conjectured lower bound for online discrepancy minimization would further reduce this cost by a square root. We also provide an extension to multiple treatments and non-uniform treatment probabilities.
The rest of the paper is organized as follows. Section 2 with a discussion of related work to position our contribution in the literature. Section 3 defines notation and formally introduces the problem. Section 4 provides our proposed algorithms and methods. Section 5 provides a detailed simulation study of the behavior of the proposed algorithms.

Related Work
There are two common approaches for achieving improved covariate balance in experiments. The first, and most common especially within industrial settings, approach is to perform a post-hoc regression adjustment which includes pre-treatment covariates [Deng et al., 2013, Lin, 2013. The second approach is to consider covariate balance during the design phase of the experiment, i.e., explicitly optimizing treatment assignment in order to minimize imbalance between treatment groups [Greevy et al., 2004, Higgins et al., 2016, Kallus, 2018.
Post-hoc stratification may be seen as asymptotically equivalent in terms of variance reduction to its analogous prestratified design as shown by Miratrix et al. [2013]. Miratrix et al.'s (2013) analysis is limited by two factors: it assumes a fixed number of stratification cells (that do not grow with sample size) and it is conditioned on the poststratification estimator being defined (e.g. treatment and control units within each stratification cell). These limitations may weaken it's asymptotic equivalence argument. A key limitation of post-hoc adjustment approaches is that the desire for simplicity and scalability implies that practitioners typically adjust for only a linear function of the pre-treatment covariates. Indeed, in the common "CUPED" approach, adjustment is performed solely on a linear function of a single pre-treatment outcome measurement [Deng et al., 2013]. Second, many common approaches for constructing stratification cells (i.e. clustering algorithms) may be computationally infeasible in practice for industrial applications when the number of simultaneous experiments and the number of outcome variables of interest are large. Third, without sample splitting (or when naively applied) advanced machine-learning based methods for adjustment may slip in assumptions of correct-specification of the outcome model, or have confidence intervals with poor coverage properties. While cross-fitting may ameliorate some of these problems in larger samples, sample splitting may prove too high a cost when sample sizes are low.
A key shortcoming of design-based covariate balance is the lack of computationally efficient algorithms which provide theoretical guarantees over worst case behavior. Blocking [Fisher, 1935] partitions variables into non-overlapping sets and performs complete randomization within each partition ("blocks"). Higgins et al. [2016] introduced a computationally approximation of blocking which runs in O(n log(n)) time. Kallus [2018], Bertsimas et al. [2015] propose an optimization based approach, Kallus [2018] additionally considers a partially random approach using semi-definite programming. Zhou et al. [2018] provide a method combining batch-based sequential experimentation with rerandomization to achieve balance, but which is not computationally feasible in moderate to large sample sizes. Perhaps the closest to the current work is Harshaw et al. [2020] which propose a balancing design using the Gram-Schmidt walk, an offline method for (linear) discrepancy minimization. Current state of the art for balancing treatment assignment requires polynomial running time and generally requires knowing all of the covariate vectors prior to determining assignment [Higgins et al., 2016, Harshaw et al., 2020, Arbour et al., 2021. As we discuss in section 3, this is a non-starter for online treatment assignment. In this setting, subjects must be allocated as they arrive; it does no good to know how you should have assigned a user at the end of the experiment; you need to know when that subject arrives. It's crucial that when a subject in an experiment arrives they be swiftly allocated to a unit. Especially in an online environment, high latency will lead to attrition, which may counteract any potential gains from greater efficiency. Moreover, the users who attrit may be the very subjects of interest [Munger et al., 2021]. By inducing differential attrition based on patience, the sample in the experiment may differ greatly from the population of interest on unobserved characteristics that make it difficult to extrapolate to a population-level effect [Egami and Hartman, 2020].
There is also a variety of methods aimed at sequential, online assignment in experiments. The seminal work in this literature is Efron [1971] which introduced an online variant of complete randomization which aims to ensure that a pre-specified marginal treatment probability is met without introducing too much accidental bias. Smith [1984] provides a generalization of the Efron [1971] approach which extends gracefully to multiple treatments. There are a variety of online balanced coin designs which seek to reduce covariate imbalance [e.g. Baldi Antognini and Zagoraiou, 2011, Moore and Moore, 2017]. Moore and Moore [2017] is based around Mahalanobis distance. As such, it has polynomial time-complexity at each arrival time. In addition to inefficiency, the theoretical worst-case behavior of this algorithm has not been resolved, even in the stochastic setting. Theoretical guarantees of this sort are paramount in the design setting, as practitioners need to know the credibility of their inferences (and how they may differ from simple Bernoulli randomization).

Background and Problem Setting
We first fix notation. Random variables will be denoted in upper case, with sets in bold. The problem setting, which we refer to as experimental treatment allocation, is as follows. We assume that we observe 1, . . . , n i.i.d. observations of X ∈ R n×d : the covariates 1 . The experimenter is asked to assign a treatment assignment, A ∈ {1, −1} (we will later loosen this to multiple discrete treatment values). We will refer to the assignments of A as treatment and control, respectively. Each unit is imbued with potential outcomes for each treatment, the value of the outcome if that unit had been assigned to the given group: y(1) for treatment and y(−1) for control. After assignment we observe only the potential outcome corresponding to the realized treatment assignment, Y . We assume that the outcomes are not available until the conclusion of the experiment. At the end of the experiment we are interested in measuring the sample average treatment effect (SATE) between any two treatments, k and k with the difference in means estimator: where p(A i ) denotes the probability of assigning treatment A i to instance i. Note that this is simply the differencein-means rather than the more general Horwitz-Thompson estimator [Horvitz and Thompson, 1952], as the treatment probability is marginal rather than conditional. More sophisticated estimators are usable in this setting [e.g. Tsiatis et al., 2008, Aronow andMiddleton, 2013], but we will focus on the simplest as we optimize design as is commonplace for studying design [Kallus, 2017, Harshaw et al., 2020.
If propensity scores are constant, the estimator of the SATE given by equation 1 will be unbiased and consistent for its oracle counterpart, the difference of potential outcomes of the k and k treatments. This SATE is our estimand of interest, as estimated by equation 1.
We will maintain the following assumptions: The problem of experimental allocation is to observe covariate vectors and assign A to units so as to achieve desirable properties of the SATE (for instance, to minimize variance). In the most general setting where no assumptions are placed the relationship between the covariates and outcome complete randomization-randomly drawing assignments without respect to background covariates-is known to be minimax optimal [Kallus, 2018].

Robustness in sequential design
Experiments are prized for their ability to provide unbiased estimates of causal effects with relatively mild assumptions. These assumptions, in fact, typically flow from the design of the experiment rather than more difficult assumptions about the data used in the course of analysis [Sekhon, 2009, Aronow et al., 2021.
In the study of vector balancing, there are three main sampling schemes of interest, listed in order of how adversarial they are: 1. Stochastic arrivals. 2. Oblivious adversarial arrivals. 3. Fully adversarial arrivals.
In the stochastic arrivals setting, units are sampled i.i.d. from some fixed (possibly infinite) population. As such, a given covariate vector is just as likely to arrive early in the sequence as late. In the oblivious adversarial setting, the adversary knows the process which will be used to assign units to groups, but cannot condition specifically on those assignments in making its decisions about the order of arrival of units. This implies the following conditional independence: Assumption 2 (Oblivious Adversary).
This stands in contrast to the fully adversarial case, in which the adversary is able to condition its decision in each period on the full set of past assignments in addition to the history of the covariate vector. We will focus on the oblivious adversary in this paper.

Online Assignment of Treatments
Weighted Online Discrepancy Minimization: We consider a variant of the online Komlós problem [Spencer, 1977], where vectors x 1 , . . . x n arrive one by one and must be immediately assigned a weighted sign of either −2q or 2(1 − q), for 0 < q < 1, such that the weighed discrepancy n i=1 η i x i ∞ , where η i is the weighted sign given to x i , is minimized. Notice that when q = 1/2, the signs become ±1.
Algorithm 1, takes i = 1, . . . , n unit vectors in sequentially and assigns them to a treatment and control, represented by the value of η i . The procedure, an extension of recent work in online discrepancy minimization [Alweiss et al., 2020, Dwivedi andMackey, 2021], assigns treatment with probability proportional to the inner product between a running sum of the signed prior observations. The algorithm and analysis differs from prior work for discrepancy in two aspects which are necessary for use in experimentation. One is to give a ridge regression guarantee, we need to characterize the random vectors output by the algorithm in terms of the projection matrix P . Harshaw et al. [2020] do that for an offline discrepancy minimization algorithm, and here we do it for an online version. The other way we differ is that our algorithm is a slight generalization of the one in Dwivedi and Mackey [2021] in that we have a parameter q. This allows for the case of imbalanced treatment assignments. A straightforward adoption of the analysis in Dwivedi and Mackey [2021] to this case results in a worse dependence on 1/q. Therefore, we derive a sub-exponential concentration bound and get a better dependence on 1/q. Notation: Let P i , i ∈ [n] be orthogonal projection matrices onto the span of {x 1 , .., x i }, that is, Algorithm 1 takes each input vector x i and assigns it {−2q, 2(1 − q)} signs online to maintain low weighted discrepancy with probability 1 − δ.

With probability
Note that η i is defined only when |w i−1 x i | ≤ c. Therefore, η is defined with probability at least 1 − δ.
Remark 1. If w i is a is a mean zero (σ, P i ) sub-Gaussian random vector, then Cov(w i ) ≤ σ 2 P i . Similarly, we have that if w i is a is a mean zero (ν, α, P i ) sub-exponential random vector, then Cov(w i ) ≤ 3 2 ν 2 P i . A proof is given in Lemma 6.
We will use Threorem 1 to derive results on the average treatment effect using the framework developed by Harshaw et al. [2020]. Harshaw et al. [2020] give a linear algberaic expression for the error of HT-estimators in terms of µ. In particular, they show Lemma 1 (Lemma A2 and Corollary A1 in Harshaw et al. [2020]).
and hence, Theorem 1 immediately implies that assignments generated by Algorithm 4 are well balanced. But optimizing just for balance can lead to accidental bias [Efron, 1971]. We have from Lemma 1 that Var(τ ) = 4 n 2 λ max (Cov(z)) µ 2 in the worst case when µ is along the top eigenvector of Cov(z). Therefore, to control accidental bias we need to make sure λ max (Cov(z)) is not high. We achieve this by augmenting the original covariates x i by √ φe i to get parameter which controls the extent of the covariate balance, and e i is a basis vector in dimension n. By a simple calculation, we can see that running Algorithm 4 on x i augmented with √ φe i is equivalent to running it with x i and replacing w i−1 x i by √ 1 − φw i−1 x i everywhere. Therefore, we don't have to explicitly augment the covariates in the algorithm.
We note that with augmented covariates, w n = √ φη √ 1 − φX η is a sub-Gaussian or a sub-exponential random vector as in Theorem 1.

Error bounds
Proposition 5 (Concentration of ATE). Algorithm 4 when run with augmented covariates, generates a random assignment z such that with probability 1 − δ.
Lemma A10 in Harshaw et al. [2020] shows that µ Qµ = min β∈R d 1 φ µ − Xβ 2 + β 2 (1−φ) . Proposition 6. The worst-case mean squared error of the online balancing walk design is upper bounded by Proof of Proposition 6. This follows from Lemma 1 and Proposition 3. We note that Q σ 2 φ I. Proposition 7 (Ridge Connection). The worst-case mean squared error of the online balancing walk design is upper bounded by an implicit ridge regression estimator with regularization proportional to φ. That is,

Algorithm with Restart
We saw earlier that η, z are defined only with probability 1 − δ. This is because with our choice of c, only with probability 1 − δ we have for all i, |w i x i | ≤ c. This means our treatment assignment fails with probability δ. There is a simple way to make sure that the algorithm never fails and have same error bound guarantees with a slightly worse constant. This is achieved by slightly modifying Algorithm 4 so that whenever we have |w i−1 x i | > c for a particular i, we start a new instance of the algorithm for covariates x i+1 , ...x n . This is equivalent to setting w i−1 = 0 and continuing with the algorithm.
Since for any treat assignment procedure E (τ − τ ) 2 just depends on Cov(z), and Algorithm 4 fails with probability ≥ δ, we can show that Therefore, for the modified algorithm, we will have error guarantees as in Propositions 6 and 7 (but worse by at most a factor of 2) and with probability one.

Multiple Treatments
In this section, we consider an online multi-treatment setting, where each vector is assigned to a group in M = {m 1 , m 2 , . . . , m k } immediately on arrival. For each 1 ≤ i ≤ k, group m i is associated with a weight α i . The goal is to minimize the multi-treatment discrepancy: is the sum of all vectors assigned to treatment m. Notice that by setting α 1 = 1 1−q and α 2 = 1 q , we can recover the definition given for the weighted discrepancy between two treatments m 1 and m 2 .
Our algorithm can leverage any oracle (we call it BinaryBalance in Algorithm 2) that minimizes the weighted discrepancy for two treatments. Our results are obtained by using Algorithm 1.
We first build a binary tree where each leaf of the tree corresponds to one of the k treatments in M . Let h be the smallest integer such that 2 h ≥ k. We start with a complete binary tree of height h, and then remove 2 h − k leaves from the tree such that no two siblings are removed. Note that this is possible by the definition of h. We further contract each internal node with only one child to its child. This process does not change the number of leaves in the tree. Let T be the obtained tree. By construction, all internal nodes of T have 2 children and T has exactly k leaves. We then associate each leaf of T with a treatment in M . For each vector assigned to treatment m i , we also say that it is assigned to the leaf corresponding to m i , ∀1 ≤ i ≤ k. For each node v ∈ T , denote by s(v) the sum of all vectors assigned to leaves under v. In addition, let α(v) be the sum of all weights assigned to leaves under v. For each internal node v of T , the weighted discrepancy vector at v is defined as: where v l and v r are the left and right child of v respectively.
For each internal node v in T , we maintain an independent run of a two-treatment algorithm that minimizes w(v) ∞ . At a high level, we minimize the weighted discrepancies at all internal nodes simultaneously. When a new vector x arrives, we first feed it to the algorithm at root r. If the result is +, we continue with the left sub-tree of r. Otherwise, we go to the right sub-tree. We continue in that manner until we reach a leaf l. x will then be assigned to l (and the treatment associated with l).
Algorithm 2 KGroupBalance takes each input vector x i and assigns it to one of the groups online to maintain low discrepancy with probability 1 − δ.
Input: x, k, α. h ← smallest integer such that 2 h ≥ k. T ← complete binary tree with height h. Remove 2 h − k leaves from T such that no two siblings are removed.
Associate each treatment to a leaf of T . Contract each internal node in T with one child to its child. for node v in T do α(v) ← sum of all weights of groups associating with leaves under v. end for for internal node v of T do Instantiate BinaryBalance(v) ← oracle for weighted discrepancy problem at v with weighted signs α(vr) α(v l )+α (vr) and − α(v l ) α(v l )+α(vr) . end for for i from 1 to n do v ← root of T .
for v is an internal node of T do

Experiments
In this section, we provide simulation evidence on the efficacy of our proposed methods. In particular, we use a wide variety of data generating processes, many of which do not assume that units arrive i.i.d. as is often standard in simulation settings for this problem. Subjects are unlikely to arrive truly i.i.d. in the real world. Earlier arrivals will typically be more active than late-arriving units, for example.
All data generating processes used in simulations are shown in Table A1. If not otherwise specified, the sample size is 1000 subjects, the number of groups is two and the marginal probability of treatment is 1 2 . Methods compared in the simulations are simple randomization (Bernoulli coin flips), complete randomization (fixedmargins randomization), the biased coin designs of Efron [1971] and of Smith [1984], QuickBlock of Higgins et al. [2016], and Alweiss et al. [2020]. These are compared to our proposed methods which are generalized versions of the discrepancy minimization procedure of Dwivedi and Mackey [2021]. We provide three versions of our proposed algorithms, called BWD for "Balancing Walk Design". The most basic version (BWDRandom) reverts to simple random assignment for all remaining periods once |w i−1 x i | > c. Our preferred approach restarts the algorithm in this case as described in section 4.3. We examine this with two levels of robustness, φ: 0 (purely balancing) and 0.5 (a uniform mix between randomization and balancing): BWD(0) and BWD(0.5) respectively. For further details, see section A.4.
In these comparisons, BWD need not out-perform all methods in all data-generating processes. QuickBlock, for instance, is a fully off-line method, so comparable performance by an online method is noteworthy. In general, Alweiss and BWD will be most effective when the true relationship between covariates and outcome is linear, since they seek linear balance. All plots incorporate 95% confidence intervals.

Binary treatment
Timing. Our proposed method (BWD) is highly efficient, scaling substantially better than other balancing methods. This analysis of runtime directly compares online methods to a widely used offline balancing method (QuickBlock). It's important to note that the QuickBlock algorithm cannot be used in the online setting, even if its runtime did not make that prohibitive. While QuickBlock is O(n log n), the proposed online balancing methods are all linear-time. Given QuickBlock's runtime, the following simulations only include it for comparisons up to sample sizes of 10 4 .

MSE.
Next, we demonstrate in Figure 2 how imbalance minimization translates to improved estimation of causal effects, measured by the mean squared-error of our estimate of the average treatment effect. We normalize this graph based on n, the rate of convergence of the difference-in-means estimator under simple randomization. The results depend strongly on the true nature of the data-generating process. In short, on non-linear data generating processes, offline blocking performs better than anything else, but in many settings BWD converges to similar error rates as QuickBlock. On linear or near-linear data-generating processes, our proposed algorithms perform very strongly, outperforming QuickBlock even in small sample-sizes. When there is a break from the purely i.i.d. stochastic setting (such as LinearDriftDGP and LinearSeasonDGP), BWD behaves well, as expected.

Multiple Treatments
MISE. BWD gracefully extends to the multiple-treatment setting, which we demonstrate in Figure 3. This chart measures the mean integrated squared-error of our estimates of the ATEs (relative to a single control group). Figure 3 shows the results. BWD consistently outperforms existing online assignment methods by substantial margins.
An array of additional simulation results may be found in Appendix B.

Conclusion
Experiments are a crucial part of how humans learn about the world and make decisions. This paper is aimed at providing a way to more effectively run experiments in the online setting. Practitioners must commonly operate their experiments in this environment, but due to the lack of suitable options for design fall back to simple randomization as the assignment mechanism. In this paper, we have shown how the Balancing Walk Design can be an effective tool in this setting. It is efficient, effective at reducing imbalance (and, therefore, the variance of resulting causal estimates), robust and it is fully suited to the particularities of online treatment assignment. Simulations have shown it to work well across a range of diverse settings. The Balancing Walk Design can improve the practice of large-scale online experimentation. Therefore, it is sufficient to prove Also note that, since v ≤ 1 and C ≥ 1, we can absorb C into v (by taking v := v/ √ C) and therefore without loss of generality assume that C = 1 and v ≤ 1. Define: We want to show Since v 2 ≤ 1, we have α + β ≤ 1. Now considering the difference of the LHS from the RHS after multiplying both sides by β we arrive at Proof.
For all i and u ∈ R d , we have with probability 1 −q.
Definition 3 (Sub-Gaussian). A mean zero random variable X is sub-Gaussian with parameter σ if for all λ ∈ R, A mean zero random vector w is (σ, P ) sub-Gaussian if for all unit vectors u and λ ∈ R, E[exp λw u ] ≤ exp λ 2 σ 2 u P u 2 .
The following concentration bounds for sub-Gaussian and sub-exponential vectors are obtained from standard bounds for scalar sub-Gaussian and sub-exponential random variables Wainwright [2019] by scaling σ, ν and α by appropriate factors.
We get the result with C = 1.5/e. Lemma 6. If w is a is a mean zero (σ, P ) sub-Gaussian random vector, then Cov(w) σ 2 P . If w is a is a mean zero (ν, α, P ) sub-exponential random vector with 2 Proof. For all unit vector u, we have By definition, if w is a is a mean zero (σ, P ) sub-Gaussian random vector, w u is σ 2 u P u sub-Gaussian. Therefore, u Cov(w)u = Var(w u) ≤ σ 2 u P u as desired.
Proof. For the first claim, note that a random variable bounded in [a, b] is sub-Gaussian with σ 2 = (b−a) 2

4
. To prove second claim, we have for 2λ < B. The last step follows from e x ≤ 1 + x + Ax 2 for x < B. Recall thatq < 2q, we have for λ < B/2 as desired.
Proof. We have First, suppose that w i−1 is (σ, P i−1 ) sub-Gaussian, we will prove that w i is (σ, P i ) sub-Gaussian. By Lemma 7, i is 1-sub-Gaussian. Therefore, We now consider the exponent (divided by the common factor λ 2 ). It is sufficient to show Using Proposition 8 (with L ← 1, C ← c/2q) and assuming σ 2 = c/2q ≥ 1, we have Therefore, w i is a (c/2q, P i ) sub-Gaussian random vector. Now suppose that w i−1 is (ν, α, P i−1 ) sub-exponential, we will prove that w i is (ν, α, P i ) sub-exponential. Again, by Lemma 7, i is (4 √ Aq, 2/B) sub-exponential. Therefore, . Note that u P i u is greater than both u P i−1 u and u x i x i u. We have This follows from Proposition 8 by substituting L ← 16Aq and C ← c/2q and noting that ν 2 = 8Ac = Lc.

The result then follows by a union bound over i ∈ [n].
Now suppose c = 9.3 log(2n/δ). Note that A = 0.5803 and B = 0.4310 satisfy e x < 1 + x + Ax 2 for x < B. Let ν 2 = 8Ac and α = 2/B. We have Therefore, by Lemma 4, Again, the result follows by union bounding over i ∈ [n].

A.2 Proof of Theorem 2
Proof of Theorem 2. Let D(δ) be the discrepancy obtained by BinaryBalance as a function of the failure probability δ. We will show that (2 log k)D(δ/k) is the corresponding discrepancy obtained by Algorithm 2. Theorem 2 will then follow from Proposition 4. Notice that with probability δ/k, each run of BinaryBalance at an internal node of T has discrepancy D(δ/k). By union bounding over O(k) internal nodes, we have that with probability 1 − δ, all of the discrepancies are bounded by D(δ/k).
Assume all the discrepancies at the internal nodes in T are bounded, we show how to bound the discrepancy between any two treatments. Let l and l be any two leaves in T . The goal is to show that where v l , v r are the left and right children of v. By definition, Therefore, Hence, both are bounded by D(δ/k).
Now consider v 1 , v 2 and v 3 in T such that v 1 is a child of v 2 and v 2 is a child of v 3 . We have, by triangle inequality, Let l be a leaf in T and let l, v 1 , v 2 . . . r be the path from l to the root r. Repeatedly applying the above relation along the path gives Since there are at most log k nodes in the path from l to r, Finally, for any two leaves l and l , Remark 2. If all weights are uniform, the summation in (3) becomes a geometric series and can be bounded by a constant. Therefore, we can remove the factor log k in Theorem 2.

A.3 Other Proofs
Proof of Proposition 1. Let B = √ φI √ 1 − φX . We have from Threorem 1 that w n = Bη = √ φη √ 1 − φX η is a mean zero (c/2q, P ) sub-Gaussian random vector, where Therefore, by sub-Gaussianity of w n , for any vector u, we have We therefore have the sub-Gaussian claim. The sub-exponential result follows similarly.
Proof of Proposition 3. Suppose c = log(2n/δ)/q. We have from Proposition 1 that z is a ( c/2q, Q) sub-Gaussian vector. This implies that Cov(z) c 2q Q.
From Lemma 4, we have The result then follows by a union bound.
Proof of Proposition 7. This follows from Proposition 3 and the proof of Theorem 3 in Harshaw et al. [2020].

A.4 Robustness
Proposition 3 immediately gives a bound on λ max (Cov(z)) and hence bounds accidental bias. Remark 3 (Accidental Bias). With probability ≥ 1 − δ the maximum eigenvalue of Cov(z) satisfies    Figure A2: Imbalance. BWD is highly effective at eliminating linear imbalance between groups.

B.2 Binary Treatments
Bias. Figure A1 shows that none of the examined methods are biased (but that does not imply that they are robust [Efron, 1971, Harshaw et al., 2020).
Imbalance. To measure linear imbalance, we calculate the L 2 norm of the difference in covariate means. While this is far from the only measure of imbalance, it serves as an effective metric to demonstrate how well various metrics serve to eliminate linear imbalances, a common diagnostic used my experimenters. We further normalize across samplesize by multiplying by the sample size. Since all methods are unbiased, this has the effect of showing parametric convergence rates as a flat line in the graph. Unsurprisingly, methods which directly optimize for linear imbalance perform very well in Figure A2. Our proposed algorithm has better finite sample imbalance minimization than does the algorithm of [Alweiss et al., 2020] due to the finite sample improvements of [Dwivedi and Mackey, 2021].  Figure A3 shows the resulting mean squared-error attained by methods which support marginal probabilities not equal to one-half. All methods perform well, with DM performing effectively on nearly linear processes, and QuickBlock performing slightly better when the true process is highly non-linear.
Marginal probability. The evaluation in Figure A4 examines how closely each method hews to the desired marginal probability of treatment. All methods do a good job of ensuring the appropriate marginal distribution.  Figure A5: Entropy

B.4 Multiple Treatments
Entropy. To ensure that treatment is being assigned with the correct marginal probability, we can measure the normalized entropy of the empirical marginal treatment probabilities. If the values are perfectly uniform, then the value will be exactly one. As the normalized entropy decreases, the marginal probabilities are more uneven, indicating a failure to match the desired marginal distribution. BWD performs very similarly to complete randomization (which almost exactly matches the desired marginal probabilities), slightly out-performing [Smith, 1984]. Note that while Smith [1984] only seeks to optimize the marginal probability of treatment for each unit, BWD additionally balances covariates.