Before we present our main arguments, it is useful to introduce the notation we will use going forward. Following the notation in [22], we assume that the experiment of interest is $2^{1}$, meaning there is one experimental factor with two levels (i.e., a treatment and control level). Let $W_{i}=1$ if unit $i$ is treated and $W_{i}=0$ if unit $i$ is control. Say $i\in\{1,...,n\}$. The vector ${W}_{n\times 1}$ denotes treatment assignment. Let $\mathbf{X}_{n\times k}$ denote the data matrix containing $k$ pre-treatment covariates relevant to the outcome of interest. We are interested in the outcome matrix $\mathbf{Y}_{n\times 2}$, where ${Y}(0)_{n\times 1}$ denotes the complete potential outcomes under control and ${Y}(1)_{n\times 1}$ denotes the complete potential outcomes under treatment. If $\left(Y_{\textrm{Obs}}({W})\right)_{n\times 1}$ denotes the vector of observed outcome values, $Y_{\textrm{Obs},i}=Y_{i}(1)W_{i}+Y_{i}(0)(1-W_{i})$. If we assume the sharp null hypothesis (zero treatment effect for every unit) and if we leave “${Y}_{\textrm{Obs}}$ fixed and simulating many acceptable randomization assignments, ${W}$, we can empirically create the distribution of any estimator,” $g(\mathbf{X},{W},{Y}_{\textrm{Obs}}({W}))$ conditional on the null hypothesis [22]. We make the additional assumption that the experimenter is interested in using Fisherian inference and that the treatment effects are constant and additive. Next, we define randomization-based fiducial intervals because these intervals are the core of both [22]’s insight and our contribution to the literature. Researchers can also exploit the duality between intervals and tests to form a fiducial interval from a randomization test.

Is it possible to take rerandomization “too far” in the sense that exact tests are no longer informative? We argue that it is indeed possible, in the sense of accepting so few randomizations that it is no longer possible to perform a meaningful randomization test. This important method of assessing significance in exact tests becomes non-informative.

Consider the following example: there are 8 experimental units, and we would like to test the effect of one factor at two levels (one treatment group and one control group). There is one background covariate, $h$, which denotes the hour a test is conducted. For simplicity, assume that there are 8 observed values for $h$, ranging from $1$ to $8$. We can quantify the randomization balance with the square root of the quadratic loss:

where smaller values of $m$ are more desirable because they indicate better balance on $h$. If the experiment is completely randomized and no rerandomization done, the smallest possible $p$-value is $1\div{8\choose 4}=0.015$ and the average $m$ value is $1.4$. If the experiment is pair-matched and complete randomization occurs within pairs for whom $h\in\{1,2\},...,h\in\{7,8\}$, the average $m$ improves to $0.375$. However, the minimum possible $p$-value in randomization inference is now $1\div{2\choose 1}^{4}=0.063$. The minimum possible $p$-value does not vary linearly with the maximum acceptable $m$ in completely randomized experiments (see Table 2).

This example illustrates that low acceptance thresholds can invalidate the estimation of uncertainty in the sense that the minimum possible $p$-value will eventually increase above $0.05$. In addition, the minimum possible $p$-value varies nonlinearly with the level of strictness the researcher adopts in accepting randomizations. These considerations suggest that it may be possible to systematically establish an acceptance criterion that will minimize the minimum possible $p$-value while maximizing the balance improvement from rerandomization.

This insight formalizes the intuition, described in [10], that including every possible combination of treatment allocations will inevitably involve including a large number of irrelevant treatment allocations that do not help the researcher understand a possible treatment effect. “Randomization, after all, is random,” Deaton writes, “and searching for solutions at random is inefficient because it considers so many irrelevant possibilities.” One intuitive way to understand our approach to optimal rerandomization is that it has the advantage of, as Deaton writes, “provid[ing] the basis for making probabilistic statements about whether or not the difference [between treatment and control groups] arose by chance,” without the disadvantage of including many irrelevant treatment allocations, including ones which the researcher would specifically like to avoid.

An important theoretical principle here is that

This expression is true because the exact $p$-value is defined as the fraction of hypothetical randomizations showing results as or more extreme than the observed value. If there are $r$ acceptable randomizations, then at least $1$ of the $r$ randomizations is as or more extreme than the observed value, so the minimum fraction is $1/r$. When the set of acceptable randomizations gets smaller, the minimum $p$-value invariably gets larger.

Setting $n_{\textrm{Cand}}={10\choose 5}$ in Figure 2, we see the clear non-linear relationship between $p_{a}$ and the minimum $p$-value. This figure is consistent with the later Monte Carlo results of Figure 3.

The preceding expression shows how there is a point after which, as experimenters reduce further $p_{a}$, fiducial intervals get wider as well, since it becomes more and more difficult to reject $H_{0}:\tau=\tau_{0}$ for any $\tau_{0}$. In the extreme, if the minimum $p$-value is above $\alpha$, randomization-based fiducial interval will become non-informative at $(-\infty,\infty)$.

For the sake of illustration, we briefly examine an example from simulations that we will later explore more carefully (see §5.1 for design information). For each candidate randomization, we calculate $M$, where

where $p_{w}$ denotes the fixed proportion of treated units and $\textrm{Cov}(\mathbf{x})$ denotes the sample covariance matrix of $\mathbf{X}$. We then accept the top $p_{a}$-th fraction of randomizations (in terms of balance) across all possible randomizations.²²2Note: We could also accept $M$ with probability $p_{a}$ based on the inverse CDF of $M$. Due to the Multivariate Normality of X, we see that $M\sim\chi_{k}^{2}$; results are similar across approach.

We see in Figure 3 how, in small samples, different rerandomization acceptance thresholds induce different expected exact $p$-values, where the expectation averages across the data-generating process. More stringent acceptance probabilities reduce imbalance and therefore finite-sample bias and sampling variability of estimated treatment effects (see [19], also [6]). However, at a certain point, the size of the acceptable rerandomization set no longer supports significant results, as eventually there is only one acceptable randomization (leading to a minimum and expected $p$-value of 1).

As a consequence of this discussion, as the acceptance threshold gets more stringent, despite improved sampling variability of the treatment effect estimate, exact fiducial intervals will eventually become non-informative as the interval width goes to $\infty$ because the size of the rerandomization set can only yield a $p$-value above 0.05. The exact threshold when non-informativity is reached will vary with the experimental design, but the general pattern is theoretically inevitable: if we shrink the size of the set in which a random choice gets made, we also shrink the reference set for assessing uncertainty. This simulation shows that the interval widths do not necessarily get monotonically smaller until the point of non-informativity but can start to increase well before that point (the exact change point will require future theoretical elucidation). In this context, how can we characterize “optimal” ways to rerandomize? What sources of uncertainty affect this optimization?

There is an interesting connection to the discussion of randomness in rerandomization approaches here to the literature on matching contains a non-informativity of a similar nature. In that context, researchers must decide how much dissimilarity to allow between matches. This dissimilarity is often set with a caliper. [9], [1], and others have made recommendations about the optimal caliper size in the context of propensity score matching. These recommendations seek to reduce the imbalance as much as possible without ending up with no actual matches: if the caliper were exactly 0 and if the sample space were continuous, the analysis would be non-informative as no units could be matched. In both matching and rerandomization, then, we would like to create a situation where treated and control units are as similar as possible to one another. However, whereas achieving no acceptable matches affects point estimates, achieving a single acceptable randomization affects the evaluation of null hypotheses using exact inference.

This section illustrates an important tradeoff in rerandomization involving small samples: more stringent acceptance thresholds improve balance and reduce sampling variability in treatment effect estimates, but too stringent thresholds lead to non-informative hypothesis tests. Balancing these considerations, the choice of the rerandomization acceptance probability $p_{a}$ (or equivalently, the threshold $a$) is critical. In the next section, we will explore several approaches for setting the acceptance threshold.

As Morgan and Rubin write, “for small samples, care should be taken to ensure the number of acceptable randomizations does not become too small, for example, less than $1000$” (p. 7, [22]). This view suggests that $p_{a}$ should be set as low as possible so long as the number of acceptable randomizations is greater than a fixed threshold that is determined a priori. In this sense, one decision rule might be to set $p_{a}$ such that

The value of $p_{a}$ which yields $\beta$ can easily be found by inverting the formula for the minimum $p$-value. On the one hand, it is difficult to know whether a given threshold is reasonable given the number of units available, the number of covariates observed pre-treatment, and the computational resources at researchers’ disposal. On the other hand, a threshold is easy to interpret and does not involve additional optimization steps.

Our simulations illustrated how the fiducial interval widths do not get monotonically smaller with $p_{a}$ until the point of non-informativity; rather, they can start to increase before that point in small samples. Thus, in selecting $p_{a}$, there is a tradeoff between the variance gains from rerandomization and the decreasing size of the acceptable randomization set. We could therefore propose an alternative procedure to determine a $p_{a}$ value that explicitly trades off the costs and benefits of rerandomization.

Following [22], let $e_{p_{a}}(p_{a})=e_{a}(a):=\textrm{minimum $p$-value given $p_{a}$}$. Let $v_{p_{a}}(p_{a})=v_{a}(a):=v_{a}$, where $v_{a}$ loosely denotes the remaining variance between treatment and control covariate means on a percentage basis. This notation allows us to emphasize the duality between setting $a$ and $p_{a}$, where $p_{a}$ is defined as the acceptance probability and $a$ is defined as the acceptance threshold for $M$.

We want both the minimum $p$-value given $p_{a}$ and $v_{a}$ to be small: when the minimum $p$-value is small, we can form tight fiducial intervals; when $v_{a}$ is small, we achieve a large reduction in variance for each covariate. Furthermore, both $e_{p_{a}}(p_{a})$ and $v_{p_{a}}(p_{a})$ are bounded between $0$ and $1$. Thus, we want to set $p_{a}$ (or, equivalently $a$) such that

where $\lambda\in[0,1]$ determines the tradeoff between the variance reduction and the minimum possible $p$-value. It is possible to derive an efficient solution to this optimization problem using gradient descent.

The optimization framework described here has several attractive features. First, given a minimum $p$-value, we can also back out the $\lambda$ value that this implies. Thus, this framing therefore allows for the evaluation of threshold choices. Second, when $w=0.5$ and there is only one covariate, the minimum $p$-value will always be less than 0.05. Third, conditional on number of covariates, $p_{a}^{*}$ describes the same relative place on the objective function no matter the sample size, so in this sense, it is comparable across experiments. Finally, $p_{a}^{*}$ has an intuitive interpretation in that it captures the point at which we have gained most of the variance-reducing benefits of rerandomization but have not incurred significant inferential costs.

The preceding discussion of heuristic tradeoffs describes dynamics regarding the choice of randomization threshold when no prior information is available. However, when prior design information on units is accessible, we can incorporate that information explicitly to help improve exact tests in a design-based manner.

Considerations of optimal precision have a long history in the analysis of treatment effects. For example, [23] and [24] discuss minimizing confidence interval width via stratum-size informed stratified sampling [18]. In general, in the Neymanian framework, optimal experimental design is that which generates unbiased and minimal variance estimates, usually under the condition of no prior assumptions on the data-generating process.

While discussions of power in the rerandomization framework are developed in [5] and the weighting of covariates in the rerandomization [20], there is less guidance in the literature about how to consider optimality in randomization choice from a Neymanian-informed Fisherian perspective—a task we turn to in this section.

Because rerandomization-based inference is simulation-based, the a priori calculation of the optimal acceptable balance threshold is difficult (although as mentioned asymptotic results are available [19]). We here show how we can use prior information on background covariates and (b) prior information about the plausible range of treatment effects to select the covariate balance threshold so as to explicitly minimize the acceptable randomization threshold minimizing the expected $p$-value. The expectation is taken over the aforementioned sources of prior information, with the intuition being that when we have an estimate for how prior information is relevant for predicting the outcome, we can use this estimate for determining how to select an optimal point on the tradeoff between better balance and sufficient randomness to perform exact rerandomization inference.

Formally, we select the rerandomization threshold using knowledge of the design—in particular, the covariates and prior assumptions on plausible distributions over treatment effect and relationships between the baseline covariates and outcome:

with $\boldsymbol{\tau}$ being the vector of individual treatment effects, $\beta$ defining the parameters of the potential outcome model, denoted using $\widetilde{\mathbf{Y}}$, and $\mathbf{D}(p_{a})$ defining the distribution over the treatment assignment vector as a function of the acceptance probability, $p_{a}$. The approach here involves simulating potential outcomes under prior knowledge, calculating the expected $p$-value under various acceptance thresholds, and selecting the threshold minimizing this quantity (thereby maximizing power to reject the exact null).

By deciding on the acceptance threshold before analyzing the data, we use the principles of design-based inference pioneered by Neyman and contemporaries. And, by incorporating rerandomization, we improve treatment effect estimates while maintaining the ability to calculate valid exact intervals in a way that incorporates information about the structure of randomness in the randomization inference.

In summary, selecting the best rerandomization threshold in a given experiment is difficult to do in small experiments given without guidance from asymptotic results. However, we can use simulation-based methods to form an approximation for how $p$-values will on average operate under the prior design knowledge. These approximations are then useful in pre-analysis design decisions, especially as investigators can select the acceptance threshold minimizing a priori expected $p$-values.

We explore the impact of varying levels of randomness in rerandomization on inference through simulation. We simulate potential outcomes under a linear model:

where covariates, $\mathbf{X}_{i}$ are drawn from a Multivariate Gaussian with diagonal covariance and where $\epsilon_{i}^{Y(t)}$ are independent Gaussian error terms. We study an environment where errors are $N(0,0.1)$, allowing us to assess robustness to noise in the potential outcomes. We also study results across low and high treatment effect strength conditions (where $\tau\in\{0.1,1\}$). Finally, we vary the number of observations, $n\in\{6,12,18\}$.

We then generate a sequence of acceptable randomizations, ranging from 10 to the set of all possible completely randomized treatment vectors, to evaluate performance at different rerandomization acceptance thresholds. We run the pre-analysis design procedure described in §4.3, which returns an estimated expected $p$-value using the prior design information. We compute the estimated expected $p$-value from each Monte Carlo iteration. We then compare the range of estimates against the true expected $p$-value minimizer, analyzing both bias and relative RMSE, where for interpretability the RMSE estimate has been normalized by the expected RMSE under uniform selection of the acceptable randomization threshold, $p_{a}$. We assess performance averaging over randomness in the covariates and outcome.

Figure 5 displays the main results. In the left panel, we see that the Bayesian selection procedure for the rerandomization threshold described in §4.3 is somewhat downwardly biased under the non-informative priors specified here, both in the small (“S”) and large (“L”) treatment effect cases. In other words, it yields an acceptance probability that is somewhat lower than what we find to be optimal via Monte Carlo. However, in the right panel, we see that the procedure generates lower relative RMSE compared to uniform selection of $p_{a}$ across the interval $[0,1]$. This finding, again robust to treatment effect size, indicates that there is a systematic structure in the randomness surrounding the rerandomization problem that can be leveraged in pre-design analysis decisions.

Overall, we find a large decrease in relative root mean squared error (RMSE) when using the pre-analysis rerandomization threshold selection procedure described §4.3. Rerandomization allows robust causal inference in cases with high noise levels and smaller sample sizes.

Randomization intervals incorporate randomness only through variation in the treatment vector, and do not explicitly make a reference to repeated experiments (as in Neyman’s classical confidence intervals). In this case, we can produce an interval by finding the set of all null hypotheses that the observed data would fail to reject. If we recall the constant and additive treatment effect assumption (${Y}(1)={Y}(0)+\tau$), then the hypotheses can be written as

An $\alpha$-level fiducial interval for $\tau$ consists of the set of $\tau_{0}$ such that the observed test statistic would not lead to a rejection of the null hypothesis at significance-level $\alpha$. Here, when $\tau_{0}\neq 0$, the randomization test is conducted by constructing

Then, keeping $Y_{i}(0)^{*}$ fixed, we permute the treatment assignment vector, and calculate ${Y}_{\textrm{Obs}}^{*}$ by adding $\tau_{0}$ to the treatment group outcomes under the permutation. This procedure generates a distribution of $\hat{\tau}$ under the null $H_{0}:\tau=\tau_{0}$. We can use this distribution to calculate a $p$-value for $\hat{\tau}_{\textrm{Obs}}$ under the null. We can form an $100(1-\alpha)\%$ interval for $\tau$ by finding the values of $\tau_{0}$ which generate a $p$-value greater than or equal to $\alpha$. [13] discuss an efficient algorithm for obtained randomization-based fiducial intervals, which searches for the interval endpoints using a procedure based on the Robbins-Monro search process.

In this section, we define $V:=W^{C}$ for notational convenience. We assume that the covariates are Multivariate Gaussian and the imbalance metric is the Mahalanobis distance between treated and control means. Then, we can appeal to the finite-population Central Limit Theorem to assume $\mathbf{\overline{X}}_{T}-\mathbf{\overline{X}}_{C}$ is Multivariate Normal and therefore that $M(V)\sim\chi_{k}^{2}$ [11].

We argued in the above that Kasy’s optimal assignment procedure is a special case of rerandomization in which the only acceptable randomization is the one that provides optimal covariate balance. Let the $\mathbf{x}$ matrix include all observed covariates that the experimenter seeks to balance. Let $\mathbf{W}$ denote the $n$-dimensional treatment assignment vector indicating the treatment group for each unit. In the notation of Morgan and Rubin, the rerandomization criterion can be written as

With the appropriate $\varphi(\cdot,\cdot)$, Kasy’s non-deterministic assignment procedure can be considered a special case of rerandomization. That is, we will obtain Kasy’s optimal assignment if we define $\varphi(\cdot,\cdot)$ as follows, letting $R(\cdot)$ denote either Bayesian or minimax risk, letting $\hat{\beta}$ denote the choice of estimator, and letting $U$ denote a generic randomization procedure independent of $\mathbf{x}$ and $\mathbf{Y}$,

With $n$ units, we will obtain a $\mathbf{W}$ such that $\varphi(\mathbf{x},\mathbf{W})=1$ with probability $(1/2)^{n}$. The expected wait time until obtaining this $\mathbf{W}$ is given by a Geometric distribution, with mean $2^{n}$. This calculation implies that, as soon as $n\geq 20$, the expected wait-time will exceed $1$ million rerandomization draws if rerandomization is used [28].

This section has put our discussion in closer dialogue Kasy’s work in order to further illustrate how Kasy’s deterministic procedure for optimal treatment assignment is a special case of rerandomization developed in [22].