Stratified Sampling Using Cluster Analysis

The Role of Models in Generalization

In developing a method for sample selection, it is helpful to begin by reviewing what would happen in an ideal study aimed at causal generalization. In this ideal study, first a well-defined inference population P of size N would be carefully enumerated and defined. For example, in a large-scale educational experiment, this might be a list of the 65,134 regular middle schools in the United States obtained from the Common Core of Data (National Center for Education Statistics). Second, a sample S of n sites would be randomly selected from this population. For example, in an educational experiment, the sample would typically include between 20 and 60 schools or districts, depending on the study design (Spybrook 2012). Third, within these n schools, units would be randomly assigned to treatment conditions. Depending upon the design of the study, this randomization might happen at the cluster level (e.g., schools) or at a lower level (e.g., classrooms or students). This dual randomization would ensure that both the site and the treatment selection processes were ignorable or noninformative, where an ignorable selection process is one in which unobserved covariates have no effect on the conditional distribution of outcomes given the observed covariates (Rubin 1976; Smith and Sugden 1988). The result would be an “experiment within a survey,” which would clearly enable causal generalizations (Smith and Sugden 1988; Shadish, Cook, and Campbell 2002).

In practice, this dual randomization procedure is generally infeasible (Bloom 2005; Rubin 1974; Shadish, Cook, and Campbell 2002). In fact, Olsen et al. (2013) found that random site selection was implemented in only 7 of the 273 experiments reported in the Digest of Social Experiments (Greenberg and Schroder 2004). Instead, study designers and analysts often choose only one level of randomization, resulting in either a probability survey or an experiment (Fienberg and Tanur 1987; Imai, King, and Stuart 2008). In experiments, while treatment is assigned randomly, the typical practice is to select a convenience sample of n sites, where here by “convenience” we mean without clear reference to a well-defined population (Shadish, Cook, and Campbell 2002).

Given the infeasibility of random site selection, in experiments causal generalizations are typically made through the use of statistical models. To see why models are necessary, note that we can decompose a population average treatment effect (PATE) as follows (Imai, King, and Stuart 2008),

P A T E = ( n / N ) S A T E + ( 1 − n / N ) N A T E ,

where the sample average treatment effect (SATE) and nonsample average treatment effect (NATE) are the average treatment effects for those in the sample and not in the sample, respectively, the sample includes n units and there are N total units in the population. While we can estimate the SATE directly from the n units in the experiment, estimating the average treatment effect for units not in the experiment (NATE) is difficult. If the sites were selected using probability sampling, the NATE could be estimated by taking into account the selection probabilities for the units in the sample. Since the sites in an experiment are typically not selected randomly, the NATE must be predicted based upon a model relating the sample to the population.

In experiments, when generalizations beyond the sample of n sites are of interest, the results are typically analyzed using random effects models (Kirk 1995; Raudenbush 1993); these multilevel models are generalizations of the analysis of variance models typically used to analyze single-site experiments. These super-population models make generalizations through use of random effects (with an assumed distribution), not based upon selection or random assignment probabilities. Specifically, two super-population models are common: the hierarchical model (for cluster-randomized designs) or the random block model (for multisite trials). Whether estimated using traditional hierarchical linear model methods (Raudenbush and Bryk 2002) or the Neyman model (Schochet 2013), the average treatment effect estimated is therefore considered generalizable to schools or districts that are “similar” to those in the study. These “bottom-up” generalizations are difficult since what is meant by similar is typically vaguely defined.

Recent work in the analysis of experiments has focused on how to improve generalizations within this super-population framework. The methods begin by carefully defining an inference population of interest, using a population frame like the Common Core of Data or a state longitudinal data system (e.g., Stuart et al. 2011; Tipton 2013). An important feature of these population frames is that they both enumerate all units in the population and include a rich set of covariate information on these units. Next, the N schools in the population and the n schools in the sample used in the experiment are compared on a set of p covariates using a propensity score (Rosenbaum and Rubin 1983). These covariates are selected so as to meet a sample selection ignorability condition; here, ignorability is met if this set includes all covariates associated with both the site selection process and the treatment effect variability (Stuart et al. 2011; Tipton 2013). Propensity score methods are then used to reweight the sample and the population so that the two groups are balanced on this set of p covariates, where balanced means that the means (and higher order moments) of these covariates are similar for the two groups. These reweighting estimators include subclassification or poststratification estimators used in combination with the multilevel models given earlier (e.g., Hedges and O’Muircheartaigh 2011). The result is an estimate of the average treatment effect given for a well-defined population based upon a super-population model.

Retrospective generalizations, like those given earlier, are helpful in that they shift generalizations from vaguely defined to well-defined populations. However, as Tipton (2013) shows, sometimes the effectiveness of these methods is limited, particularly when there are coverage errors. Coverage errors arise when there exist segments of the population for whom there are no similar units in the experiment; this, we argue, is the lasting effect of the bottom-up generalization approach. In this article, we focus instead on developing a method for site selection that makes these generalizations “top-down.”

Instead of carefully defining the inference population after the study is completed, our goal here is to begin the study with a well-defined inference population and then design a site-selection process that makes model-based generalizations possible.

Balanced Sampling

Just as in the ideal study, the goal here is to begin by first carefully defining and enumerating an inference population P of size N and then to select the n units in the sample S strategically, using a more formal sample selection plan. The goal is for the sample selection process to be noninformative or ignorable (Rubin 1976) by which we mean that the resulting sample and population are compositionally similar on the set of covariates that explain treatment effect variability (Stuart et al. 2011; Tipton 2013).

The approach we develop here builds on theory and results found in model-based sampling in the sample survey literature (Valliant, Dorfman, and Royall 2000). Model-based sampling is a purposive sampling alternative to design-based random sampling methods; while random selection can be used within the model-based framework, it is not required. The model-based sampling approach, while not as commonly used in survey sampling, has much in common with the model-based super-population approach commonly used to analyze experiments (Fienberg and Tanur 1987; Rao 2005). As such, this method will allow inferences from a sample to super-population using the same random effects approach currently used in the analysis of large-scale experiments; the key difference here is that now the definition of the super-population will be clearly and carefully defined.

In order to develop this strategy, we first define the causal effects of interest in the potential outcomes framework. For each unit in the population P, let W = 1 if a unit is assigned to the treatment condition. Then assume that for each unit, there exists two potential outcomes, Y(1) = Y(W = 1) and Y(0) = Y(W = 0), where Y(1) is the unit’s potential outcome under treatment and Y(0) is the unit’s potential outcome under some specified alternative condition. Now for each unit in the population P, let the potential treatment effect Δ = Y(1) − Y(0). Note that as a result of the Fundamental Problem of Causal Inference, both potential outcomes—and by extension the potential treatment effect—can never be observed for a particular unit (Holland 1986). However, the goal of an experiment is generally to estimate the PATE τ^P = E[Δ], where the expectation is across all units in the inference population P. Finally, let Z = 1 if a unit in the population P is selected into the sample S. Since only units in the sample (Z = 1) are in the experiment, we focus here on causal impact estimators that are based only on the n units in S.

An important question is under what conditions an estimator T that is unbiased for the SATE, τ ^S = E(Δ|Z=1), is unbiased for the PATE, τ ^P = E(Δ). Stuart et al. (2011) and Tipton (2013) show that bias arises in relation to covariates that explain variation in potential treatment effects. Since we do not know, a priori, what these covariates are (since we have yet to conduct the experiment), this requires us to propose a model. We might posit a simple model

Δ = β 0 + β 1 X ,

1

where X is a single covariate with known values for all units in the population. For example, based on theoretical or previous empirical findings, we may believe that the effect of a school-based reading intervention (Δ) linearly increases or decreases in relation to last year’s school average reading test scores (X).

If in fact the potential treatment effects do vary in relation to X, then by selecting the sample so that the average value of X in the sample and population is the same (i.e., E(X|Z = 1) = E(X)), an estimator T of the SATE τ ^S would be model unbiased for the PATE τ ^P . We call this a balanced sample since the sample and population are balanced on the covariate X. Importantly, this idea of balance is similar to the goal of balance found in retrospective methods; the key difference is that retrospective balance is achieved through reweighting, while here balance is achieved through the strategic selection of the sample.

The idea of a balanced sample is in many ways similar to the naive sense of a “representative sample” (Royall and Herson 1973), where here by representative we mean that the sample is a like a “miniature” of the population (Kruskal and Mosteller 1980). Importantly, while random sampling will result in a balanced sample on average, it is not the only or best method for achieving such balance, particularly when the sample is small. Developing a method that performs well in small samples is of particular concern here since in most cluster-randomized or multisite studies in education, the number of higher level units (e.g. schools or districts) tends to be between 20 and 60 (Spybrook 2012).

Bias-Robust Balanced Sampling

Based on the model 1.1, we might propose to develop a sample selection plan that results in a sample that is balanced for the covariate X (i.e., E(X|Z = 1) = E(X)). This necessarily leads us to ask: what happens if our model is wrong? For example, suppose that potential treatment effects actually vary in relation to the model

Δ = β 0 + β 1 X + β 2 X 2 + β 3 W .

2

For example, it may be that the potential treatment effects vary nonlinearly in relation to school average pretest scores (X) and vary also in relation to the proportion of the school that is minority (W). Now the bias of an estimator T for τ ^P can be written as

E T − τ P = β 2 E X 2 | Z = 1 − E X 2 + β 3 E W | Z = 1 − E W .

Clearly, T is only unbiased for τ ^P if E(X ²|Z = 1) = E(X ²) and E(W|Z = 1) = E(W). This means that a sample that is balanced on X is only adequate if our former model (1) holds but not under this new model (2).

In model-based sampling, the goal is to select the sample so that the method is bias robust, where bias robust is shorthand for “bias-robust-against-model-failure” (Valliant et al. 2000). In this framework, in causal generalization this requires us to first propose multiple models relating the potential treatment effects to possible moderators, and second, to develop a sample selection plan that guards against selection of the wrong model. Here the key tool guarding against model failure is the selection of a balanced sample. Subsequently, we define this more formally.

Definition: balanced sample (of order R)

Let X = {X₁ , …, X _p } be a set of covariates that might explain variation in the potential treatment effects Δ. A sample S is said to be a balanced sample (of order R) if for each covariate X_h (for h = 1 … p) and for each power r = 1 … R, E(X_h ^r |Z = 1) = E(X_h ^r ), where E(X_h ^r |Z = 1) is the average value of X_h ^r in the sample S and E(X_h ^r ) is the average value in the population P

When a balanced sample (of order R) is selected on p covariates, it is easy to show that an unbiased estimator T of the SATE τ^S is unbiased for the PATE τ^P if the true model is of the form

Δ = δ 0 β 0 + ∑ h = 1 p ∑ r = 1 R δ h r β h r X h r ,

3

where for h = 1 … p, β _hr is the regression coefficient associated with the covariate X_h ^r and δ _hr = 1 if the coefficient β _hr is included in the model and zero otherwise. For example, if a sample is balanced (of Order 2), then T is unbiased when the true relationship is linear

Δ = β 0 + β 11 X 1 + β 21 X 2 + … + β p 1 X p ,

4

quadratic,

Δ = β 0 + β 11 X 1 + β 21 X 2 + … + β p 1 X p + β 12 X 1 2 + β 22 X 2 2 + … + β p X p 2 ,

5

when it includes only a subset of the covariates,

Δ = β 0 + β 11 X 1 + β 21 X 2 ,

6

or any other combination of the models. Since we cannot and do not know the true model, in the language of model-based sampling, the goal is to select a sample using a bias-robust strategy. In causal generalization, this is a strategy that leads to an unbiased estimate of τ ^P = E(Δ) under a variety of possible models.

A sample selection strategy can be bias robust in two ways. First, the sample becomes more bias robust as the dimension p of X increases. This is because in practice when balance is achieved on a wide variety of covariates it is often approximately achieved on other covariates, including those that may have been omitted from the model (Stuart 2010; Smith and Sugden 1988; Royall and Herson 1973; Brewer 1999; Rubin and Thomas 1996). Second, for a fixed set of covariates X, a sample becomes more bias robust as the order of R increases. For example, a sample which is balanced on not just the first moments but also the second moments of X—a balanced sample of Order 2—will lead to an unbiased estimate whether the true model is linear or quadratic (Valliant, Dorfman, and Royall 2000).

The idea that the sample should be selected using a bias-robust strategy with the goal of achieving balance on both the means and the higher moments of multiple covariates is familiar—it is exactly the post hoc approach used in the propensity score literature (Rosenbaum and Rubin 1983). Propensity scores are commonly used to achieve balance under a variety of models when random processes are not possible (e.g., quasi-experiments; post hoc generalizations from experiments) or when they fail (e.g., attrition in experiments; nonresponse in surveys). While the method we develop here does not use a propensity score approach, the goals of the procedure are the same—to replace a random process with a model and to design the study in such a way that the results do not depend heavily on any one model.

Once the sample of n units is selected, the robustness of the sample to model misspecification can be evaluated by comparing the R moments of the p covariates under study in the realized sample and population. Ideally, these differences will be small, enabling the use of the simple multilevel random effects estimator commonly used in experimental analysis. If differences remain, post hoc methods like regression adjustment or those introduced by Hedges and O’Muircheartaigh (2011) and Tipton (2013) can be used to decrease bias. Importantly, since the sample was selected with an eye toward balance, however, the achieved sample and population are likely to be more similar than if a bias-robust selection method had not been used. The fact that the sample and population are more similar will mean that there will be fewer coverage errors and that, if adjustment is needed, the cost in terms of variance inflation will be smaller (Tipton 2013).

Defining a Stratified Estimator

Recall that in the balanced-sampling framework, the goal is simply to select the sample so that it is like a “miniature” of a well-defined population. To do so, we begin by positing a model explaining variation in potential treatment effects (e.g., Model 1.1), and then propose alternative models (a bias-robust strategy). To this end, we need a method for sample selection that allows for balance on orders greater than one (i.e., not only on first but also on second or higher moments) and enables X to include a large and varied set of covariates. While many possible methods exist, in model-based sampling, the simplest method that achieves this goal is that of stratified sampling with proportional allocation (Valliant, Dorfman, and Royall 2000).

Stratified sampling is already widely used in both survey sampling and in large-scale experiments. In probability sampling, strata are used to reduce the variance of an estimate; since the focus is on variance reduction (not bias), here it is common for the strata to be created on only one or two covariates (Lohr 1999). The use of one or two covariates for strata creation is also sometimes used for nonrandom site selection in experiments. For example, experiments sometimes attempt to include sites in rural and urban locations (i.e., urbanicity) or in regions throughout the country (e.g., northeast, southeast, midwest, west). Here while strata are created in order to improve generalizations, the method and framework for making these generalizations are typically informal. In contrast, in this article we propose to create strata with the goal of reducing bias through the inclusion of many covariates—hopefully all covariates that explain treatment effect heterogeneity—and in order to balance the sample and population on both means and higher order moments. In contrast to prevailing practice, our goal is to explicitly state, develop, and evaluate a bias-robust strategy for sample selection aimed at generalization.

If a simple Model (1) is of interest, in order to create a stratified sample with proportional allocation, three steps would be involved. First, the values of the covariate X would be divided into strata. If X is categorical with k categories, then k strata could be naturally created. If X is continuous, many possible strategies for strata creation could be used; one of the simplest is to define the j = 1 … k strata so that each contains an equal portion (i.e., w_pj = N/k = N_j/N) of the population units (Cochran 1968). Second, under proportional allocation, stratum j would be allocated n_j = w_pj × n units in the sample. This means that in stratum j, n_j/N_j units would need to be selected into the experiment. Third, in stratum j, the sample would need to be selected so that E(X|Z = 1, j = j) = E(X|j = j), which is to say that the stratum-specific mean of X is the same in the population and the sample.

The main benefit of using a stratified sampling approach with proportional allocation is that the resulting sample is self-weighting. This means that the sample and population are balanced on the covariate X, since

E X | Z = 1 = Σ n j / n E X | Z = 1 , j = j = Σ w p j E X | j = j = E X .

7

The fact that the sample is self-weighting means that no additional adjustments are needed and that the usual multilevel random effects model can be used for estimating the average treatment effect and making generalizations. Since the standard estimator can be used, this means that the sample selection process does not impact the power analysis (used to determine the sample size n). As a result, the issues of statistical power and generalization can be separated, which is logistically helpful in designing the experiment. It is also helpful that stratified sampling with proportional allocation is conceptually easy to understand and explain, making it appealing in the policy context where many of the results of large-scale experiments are used and interpreted.

When X is a single covariate, k strata can be easily created leading to a balanced sample. However, when the dimension of X is large and when X includes both categorical and continuous covariates, strata creation becomes more difficult. For example, if each covariate takes only two values, this would lead to 2 ^p unique strata, which can easily be larger than the total sample size n. For this reason, a method for stratification is needed that allows for the selection of balanced and bias-robust samples with fewer strata.

In the survey sampling literature, the problem of stratified sampling under a multivariate and continuous X has been addressed through the use of cluster analysis. These methods were pioneered by Day and Heeler (1971) and Golder and Yeomans (1973) who used cluster analysis to classify units in a pricing experiment and city wards in a study of seat belt use, respectively. Since 1980 cluster analysis methods have been used to create strata for the Current Population Survey (CPS) and other demographic surveys administered by the U.S. Census Bureau (Mansure and Reist 2010; Murphy 2008) as well as in surveys by census agencies around the world (e.g. Northern Ireland Statistics and Research Agency [2001] 2002). Additionally, in education research clustering methods have been used to classify schools (e.g. ETS 2008; Sleegers, Bergen, and Giesbers 1994; Lipson et al. 2004), and in the analysis of experimental data clustering methods have been developed to combine similar treatment groups (e.g. Cox and Spjotvoll 1982; Scott and Knott 1974; Tukey 1949) and to analyze subgroup average treatment effects in experiments (e.g., Peck 2005; Peck et al. 2012). Here, we propose to use cluster analysis to define strata for a stratified sampling design aimed at improving generalizations from experiments.

Cluster Analysis Method

The goal of cluster analysis is to divide units in the inference population into strata so that units in the same stratum are more similar than units in different strata. We focus here on the k-means partitioning method, though other methods are available (Everitt et al. 2011). The basic idea of k-means clustering is to create k strata and then to assign each unit to one stratum so that a measure of similarity is maximized. There are two main steps involved in a cluster analysis, and we briefly review these here.

Choosing a Distance Metric

In this article, we use cluster analysis to group units into strata that are as close to homogeneous as possible. This means that a measure of distance or similarity is required to define “close.” There are two common distance measures that we argue are useful for our purposes, though others are available (Everitt et al. 2011). The decision to use one of these metrics over the other will largely depend on the type of covariates included in X and on information or assumptions regarding the importance of each covariate.

When all of the covariates in X are continuous, the distance metric most commonly used is the weighted Euclidean distance,

d i i ′ e = ∑ h = 1 p w h X i h − X i ′ h 2 ,

8

where each covariate X_h has weight w_h , and X_ih and X_i’h are the values of the hth covariate for units i and i’. One option for weights is to set w_h = 1 for all covariates h, in which case d^e _ii’ is the Euclidean distance; this gives the most weight to covariates with the largest variances. Alternatively, if there is no information regarding which covariate is a more important or better predictor of treatment effect heterogeneity, then the obvious solution is to use inverse-variance weights, where w_h = 1/V(X_h ). In this framework, the weighted covariates have a common variance of one, and each covariate contributes equally to the distance metric (though other weighting methods are possible).

Alternatively, when X includes both continuous covariates and categorical or dummy variables, any Euclidean-based measure does not perform as well. There are various solutions to this problem. One solution that is commonly used and easily implemented is to use a general similarity (distance) measure proposed by Gower (1971),

d i i ′ g = ∑ h = 1 p w i i ′ h d i i ′ h ∑ h = 1 p w i i ′ h d i i ′ h ∑ h = 1 p w i i ′ h ∑ h = 1 p w i i ′ h ,

9

where d_ii’h is the similarity between units i and i’ on the covariate X_h . Note that this measure of distance, d_ii’h , can differ for different variable types. For dummy and categorical variables, d_ii’h = 1 if the two units i and i’ have the same value and 0 otherwise. For continuous variables, it is standard to use

d i i ′ h = 1 − X i h − X i ′ h R h ,

10

where |.| indicates absolute value and R_h is the range of observations for the covariate X_h . Using distance measures defined this way ensures that d_ii’h ∊ [0,1] for covariates of all types. Additionally, this general distance measure allows for missing data. For example, the weights can be determined so that w_ii’h = 0 if the outcome of the X_h is missing for either or both of units i and i’. Again, other weighting schemes can also be used, particularly if information on the importance of particular covariates in explaining potential treatment effect heterogeneity is available.

Sample Allocation, Selection, and Evaluation

Once the number of strata k has been selected, a strategy for sample selection within each cluster must be developed. In this section, we detail the steps involved in this process.

Calculate and Rank Within Stratum Distances

Since the overall goal is to select a balanced sample (of some order R), a method for sample selection within each stratum is needed. The goal is to simply select a balanced sample (of Order 1) within each stratum. This means selecting a sample such that in each stratum j = 1 … k, E(X_h |Z = 1, j = j) = E(X_h | j = j) for each covariate X_h in X and where Z = 1 indicates units in the sample. One method for meeting this goal is to select the n_j units randomly. This would ensure that on average (over repeated samples) units in stratum j would be balanced on both the observed covariates (in X) and the unobserved covariates (not in X). While this certainly increases the bias robustness of the method by guarding against having the wrong model, in any particular sample this method may not result in balance, particularly when the sample size n_j is small. This is particularly important since the multilevel random effects model used for generalization makes inferences in reference to the realized sample.

An alternative method for selection is as follows. First, in each stratum, calculate E(X_h | j = j) for each covariate X_h . Then for each of the i = 1 … N_j units in stratum j, a measure of distance is calculated. One strategy is to use the weighted Euclidean distance

d i j = ∑ h = 1 p w h X i j h − E ( X h | j = j ) 2 ,

11

where w_h is the weight given to covariate X_h , E(X_h |j = j) is the average value of the hth covariate in stratum j, and X_ijh is the value of the hth covariate for unit i in stratum j. Thus, each unit i in stratum j has one total combined distance measure d_ij . Again, just as discussed previously, different weights could be used, particularly if it is hypothesized that one covariate matters more in terms of explaining potential treatment effect heterogeneity than another covariate. Note that to the degree that the strata are homogenous (having small within-stratum variances), we would expect these distances d_ij not to vary by much.

Based on these within-stratum distance measures, each of the N_j population units within a particular stratum j can be ranked from smallest to largest. This ranked list can then be used for selecting the n_j units into the experiment. For example, a recruiter might start at the top of the list with the unit ranked “1,” then if the unit does not agree to be in the experiment, move on to the unit ranked “2,” and so on until n_j units agree to be in the experiment.

Subgroup Analyses

When treatment effects do in fact vary, a single PATE is clearly inadequate for summarizing the effectiveness of a program or an intervention. One solution is to additionally report subgroup average treatment effects. A benefit of our site-selection approach is that a separate average treatment effect can be calculated for each of the k strata. This strategy is similar to the method for subgroup creation proposed by Peck (2005) in general, and, when the covariates contained in X also predict treatment compliance (or other post-random-assignment groupings), to the methods proposed by Peck (2003), (2013) and Bell and Peck (2013).

Eligibility Problems

It is possible that some of the N units in the population P are not eligible for inclusion in the experiment. For example, certain schools may already use the program under study, or when resources are limited, travel outside a particular area may be infeasible. The goal of bias-robust sample selection is still useful here, particularly since it requires study planners to determine whether the reasons for ineligibility explain variations in treatment effects (warranting inclusion in X). For example, if only large schools are eligible to be in an experiment (because of the design of the study), then generalizations to all schools in the population are only warranted if treatment effects do not vary in relation to school size. Otherwise, the population needs to be redefined to focus only on large schools. Second, stratified sample selection using cluster analysis can also be used here. Again, the strata would be created based on the N units in the population P, though only eligible schools or sites would be included in the ranked list for sample selection. As a result of ineligibility, however, the sampling fractions n _j/N _j may differ by stratum. Tipton et al. (2014) includes a full discussion of sample selection under eligibility constraints, including an alternative bias-robust approach using propensity score methods.

Nonresponse Analysis

One of the biggest practical concerns in sample selection in experiments is the fact that many sites will not agree to be in the experiment. A benefit of this approach is that by articulating at the outset a set of covariates (X) that possibly explain treatment effect variability (thereby leading to bias), information on sites that refuse to be in the study can be tracked and later compared to those that accept and to the inference population of interest using available nonresponse analysis methods (e.g., Little and Rubin 1987). Additionally, this method requires designers to plan for refusals and even allows for the collection of additional information on those who refuse or accept.

Relationship to Post Hoc Methods

The goals and framework for balanced sampling have much in common theoretically with the goals of propensity score matching for post hoc adjustments for generalization. Even when using the stratified sampling method introduced here, differences may remain between the achieved sample and the population, particularly when nonresponse is high. However, as illustrated in the example, these remaining imbalances can be more easily adjusted if generalization is planned for; this is both because any remaining imbalances are typically smaller (Table 2) and since coverage errors are greatly reduced (Figure 2), both of which make reweighting procedures more effective.

Random Selection and Design-Based Inference

Given the infeasibility of random sampling in experiments, this article has provided a method for site selection that is strategic, model based, and nonrandom. However, the method we develop for site selection—stratified sampling with proportional allocation based on cluster analysis—can also be used in a design-based framework with probability sampling.

Data Frame Concerns

A potential weakness of the stratified selection method developed here is that it requires an available sampling frame. This results in two limitations to our approach. First, like the post hoc methods available for generalization, this sampling frame needs to include a rich set of covariates on all units in the population, where here the covariates that matter are those that explain variation in treatment effects. Since most censuses of schools and districts focus on demographics, using any model-based approach could result in omitted variable bias. An important feature of this method is that it requires a thoughtful discussion of the benefits and limitations of the achieved sample for generalizations and to whom, where, and under what conditions or assumptions these generalizations are most warranted.

Second, the fact that the sampling frame must enumerate all N units in the population means that the method we have developed here is most useful when the sampling units are aggregates, since data on aggregates (e.g., schools) are typically publically available, while data on individuals are not. However, while the strata creation method developed here may not be practical when the unit of analysis is the individual, we argue that the bias-robust framework is still useful. Stratified sampling is only one method for creating a balanced sample in the bias-robust framework. Future research should investigate the practicality of methods that do not require such detailed population frames—for example, quota sampling and respondent-driven sampling (e.g., Smith 1983; Watters and Biernacki 1989).