Analysis of Regression Discontinuity Designs with a Binary Moderating Variable

Abstract

Regression discontinuity (RD) designs are increasingly used for causal evaluations. However, the literature contains little guidance for conducting a moderation analysis within an RDD context. The current article focuses on moderation with a single binary variable. A simulation study compares: (1) different bandwidth selectors and (2) local polynomial regressions with interactions to local regressions on subsets of the data defined by values of the moderating variable. We find that existing bandwidth selectors optimized for main effects will choose bandwidths that are too small for moderation analysis. Additionally, choosing an optimal bandwidth for the subset regression approach may not be feasible for small to moderate sample sizes unless the moderator prevalence is near 0.5 and correlation with the assignment variable is small. We conclude that when sample sizes are small a global regression approach is likely to be preferred to utilizing bandwidth selectors optimized for main effects.

Keywords

regression discontinuity designs causal inference moderation analysis effect heterogeneity

Introduction

Randomized controlled trials (RCTs) are often viewed as the preferred design for rigorously evaluating the efficacy and effectiveness of interventions. For instance, the Common Guidelines for Education Research and Development jointly issued by the U.S. Department of Education and the National Science Foundation says that efficacy and effectiveness research should, “Generally and when feasible … use designs in which the treatment and comparison groups are randomly assigned” (U.S. Department of Education & National Science Foundation, 2013; p. 21). The preference for RCTs exists because a well-executed RCT ensures there are no systematic differences between groups at baseline so that observed effects can be causally attributed to the intervention with a high degree of confidence (Bloom, 2005; Morgan & Rubin, 2012; Shadish et al., 2002). However, RCTs are not always a feasible design to implement as they tend to be more costly, take more time to manage and conduct, and sometimes are not politically palatable to leaders and policymakers. For instance, when educational programs are being evaluated, RCTs may not be viable because educational leaders do not wish to deny students, teachers, or schools a potentially beneficial intervention, even if only for a short time as in a delayed-intervention design (Shadish et al., 2002).

It may be particularly difficult to justify a randomized evaluation design when working with a local education agency (LEA) or state education agency (SEA) to evaluate a policy or program targeting schools and students most in need (Center on School Turnaround, 2017; Hitt et al., 2018). These entities are unlikely to want to deny services to any of the targeted schools or students, and they are often under pressure to make decisions and put strategies in place sooner rather than later, which may not leave time for a prospective research design such as an RCT.

When all the neediest schools or students within a given LEA or SEA receive an intervention, it can be hard to make causal inferences about the intervention based on observational data. These schools and/or students are often located in disadvantaged communities that suffer from a lack of resources and supplies, insufficient facilities, and low morale. Additionally, they can struggle to recruit highly qualified teachers (Podolsky et al., 2016). In short, students and schools in need possess a unique profile all their own. As such, it is impossible to form a group of matched comparison units with similar profiles. One alternative is to identify students or schools in a neighboring school or SEA/LEA, but obtaining data for these purposes and controlling for between-school or SEA/LEA differences adds further complexity. Evaluators might consider leveraging longitudinal data and conducting an interrupted time series analysis (e.g., Bloom, 2003). However, without a comparable within-district comparison group measured over the same time frame, the validity threat posed by history effects is often difficult to dismiss.

In such situations, a regression discontinuity (RD) design (Cook et al., 2008; Hahn et al., 1999; Imbens & Lemieux, 2008; Shadish et al., 2002) can be a particularly attractive alternative. If need at the school or student level is measured by a continuous indicator and services are provided for those below a certain cutoff value relative to that indicator, an RD study with strong causal validity is possible (Aiken et al., 1998; Berk et al., 2010; Shadish et al., 2011). Regression discontinuity evaluations are possible whenever units are selected for treatment exposure based on whether their value on a continuous measurement falls above (or below) a threshold referred to as a cut-point. Units with values on one side of the cut-point receive the intervention while units on the other side of the cut-point comprise a comparison group. We refer to the continuous measurement that determines treatment selection in this article as the assignment variable.

When the selection process just described occurs, RD designs allow for the identification of the local average treatment effect (LATE¹) at the cut-point. However, as a practical matter it can be difficult to obtain accurate estimates of the identified LATE parameter. It will usually be the case that the assignment variable is correlated with the outcome of interest. This implies that accurate estimates of the LATE will need to condition on the assignment variable. However, by construction, there is no overlap between the treated group and the comparison group with respect to this variable. The lowest values of the assignment variable occur in the intervention group, and the highest values occur in the comparison group (or vice versa if the intervention is delivered to those above rather than below the cut-point). This lack of overlap implies that estimation of the LATE depends on accurate estimation of a regression function at a boundary point (i.e., at the cut-point where the regression relationship is assumed to change discontinuously). Obtaining accurate estimates of regression functions at boundary points is notoriously difficult. Thus, as a practical matter, the estimated LATE can vary widely based on different assumptions about the functional relationship between the assignment variable and the outcome variable (Lee & Lemieux, 2010). The use of polynomial functions is standard, but the appropriate degree of the polynomial to be used is not clear (Lee & Lemieux, 2010). Subject to other assumptions that are standard in all evaluations with causal objectives, when the functional form that relates the assignment variable and the outcome variable is correctly specified, the obtained estimate of the LATE is unbiased. Obtaining an unbiased estimate of the average causal effect is, of course, the sine qua non of causal evaluations.

One way to limit sensitivity to functional form assumptions is to limit the data used to estimate the LATE to a small region around the cut-point. The size of this region is known as the “bandwidth” in the RD literature. A smaller bandwidth will generally reduce bias in the estimate obtained but will result in an estimator with a larger standard error. Therefore, choosing a bandwidth for estimation of the LATE involves a trade-off between bias and variance.

An emerging consensus in the RD literature is that estimates should be obtained using data within a bandwidth chosen by using one of the “bandwidth selectors” that is available in the literature (Calonico et al., 2014; Imbens & Kalyanaraman, 2012). Put differently, global regressions that utilize all available data should be avoided. For instance, the U.S. Institute of Education Sciences’ What Works Clearinghouse review standards will only give the highest rating to RD studies that compute estimates within a “justified bandwidth” (WWC, 2022, p. 74). Additionally, a recent review paper discussing RD studies states that “the recommendation is to always use a data-driven procedure that is optimal for a given criterion” (Cattaneo & Titiunik, 2022, p. 834). Finally, Cattaneo and Vazquez-Bare (2016; p. 134) say of RD designs that “the most important task in practice is to select the appropriate neighborhood near the cutoff, that is, to correctly determine which observations near the cutoff will be used.” Estimating treatment effects in an RD design within an empirically justified bandwidth is desirable because these bandwidths are developed to minimize the mean squared error (MSE) of the LATE when the functional form of the equation relating the assignment variable and the outcome is mis-specified. As such, using an empirically justified bandwidth reduces sensitivity to functional form assumptions by ensuring that the estimate is robust to mis-specified models.

Evaluators are not only interested in determining “what works” but are also focused on understanding “for whom and under what circumstances” programs work (Dong, Kelcey & Spybrook, 2023). This has resulted in a growing literature emphasizing the importance of understanding heterogeneity in treatment effects (e.g., Athey & Imbens, 2016; Bryan et al., 2021). Interest in this nuance is fueled by a need to identify interventions that might be more impactful for some subgroups of individuals and less so for others, ultimately allowing policymakers to efficiently allocate resources to accomplish the greatest good (Haegerich & Massetti, 2013; Supplee et al., 2013). Perhaps the simplest approach to exploring treatment effect heterogeneity involves estimating statistical models to determine the extent to which the average treatment effect varies as a function of a single covariate measured prior to any units experiencing the intervention. This process is often called “moderator analysis” (Baron & Kenny, 1986). Note that additional assumptions are required for drawing causal inferences about moderators that have not been randomized (Bansak, 2021; Bansak & Nowacki, 2022).

Given the increasing interest in exploring treatment effect heterogeneity, it is surprising that the literature appears to contain no advice on the topic of selecting a bandwidth for conducting a simple moderation analysis within the context of an RD study. The current article aims to fill this gap. For simplicity, we focus attention on moderation with respect to a single, binary, variable. Even this seemingly simple case involves multiple complexities that the evaluator must address to obtain a good estimate of a moderation effect.

Both authors were part of a National Center for Education Research-funded evaluation exploring the impact of Accelerating Literacy for Adolescents (ALFA) Lab on ninth-graders’ reading achievement. That project served as the genesis for the two research questions that guide the current work. The first research question is, “Should we conduct the moderator analysis: (1) by estimating separate LATEs for the two subsets of the data defined by the different values of the moderating variable or (2) by estimating a single regression model containing interactions with the moderating variable?” We refer to the first option as the subgroup regressions approach and the second option as the interacted regression approach. The second research question is, “Should we use an existing bandwidth selector if the quantity of interest is not the LATE itself but rather the extent to which a binary variable serves as a moderator of the LATE?.” Furthermore, “if we should use a bandwidth selector, which one should we use?.” Neither of these questions has been addressed in the existing literature.

The rest of the article proceeds as follows. The next introductory subsection explains how moderation analysis in RD studies differs from moderation analysis in RCTs. The following introductory subsection presents a brief review of the existing literature on treatment effect heterogeneity in RD studies. The next section describes the methodology of a simulation study that was designed to help answer the two questions posed in the previous paragraph. The subsequent subsection presents the results of the simulation study, followed by an empirical illustration of the analytic approaches described in the methodology section applied to data from the ALFA Lab evaluation mentioned above. The final section concludes with a discussion.

Differences Between Moderation Analysis in an RCT Versus an RD Study

To fix ideas, assume that there is an educational evaluation that aims to understand the causal effect of an intervention on student achievement. The binary moderator of interest is eligibility for free or reduced lunch (FRL). Let Y_k be the value of the outcome variable for student k, T_k be an indicator for intervention group status (T_k = 1 if individual k is in the intervention group and 0 otherwise) and M_k be an indicator for moderator status (M_k = 1 if individual k is FRL eligible and 0 otherwise). Let n_ij be the number of students in the sample for whom T_k = i and M_k = j and σ_ij be the standard deviation of the outcome variable among students for whom T_k = i and M_k = j.

In an RCT, the standard approach to moderator analysis is to fit the interacted regression model given in equation (1):

E (Y_{k}) = β_{0} + β_{1} T_{k} + β_{2} M_{k} + β_{3} T_{k} M_{k}

(1)

The estimate of

β_{3},

{\hat{β}}_{3}

is the estimated moderating effect of

M_{k}

and the test statistic,

t_{i n t} = \frac{{\hat{β}}_{3}}{S E ({\hat{β}}_{3})}

(2)

can be compared to the relevant t-distribution to conduct a statistical test.

Alternatively, a data analyst could obtain the subset of the data for which M_k =0 and fit the model in equation (3):

E (Y_{k}) = α_{0} + α_{1} T_{k}

(3)

and then obtain the subset of the data for which M_k = 1 and fit the model in equation (4):

E (Y_{k}) = γ_{0} + γ_{1} T_{k} .

(4)

Now

\hat{γ_{1}} - \hat{α_{1}}

estimates the moderating effect of

M_{k}

. Furthermore, if we assume

σ_{11} = σ_{12} = σ_{21} = σ_{22}

and

n_{11} = n_{12} = n_{21} = n_{22}

then the statistic

t_{s u b} = \frac{\hat{γ_{1}} - \hat{α_{1}}}{\sqrt{V A R ({\hat{γ}}_{1}) + V A R ({\hat{α}}_{1})}}

(5)

is equivalent to the t-statistic in equation (2). Even if the equal variance and/or equal sample size assumptions do not hold, t_sub and t_int will not be too far apart unless one subset of the data has a particularly high variance combined with a particularly low sample size (Bloom, 2005). In other words, there is no need to decide between the two approaches in an RCT because they will provide equivalent answers.

Additionally, because there is no need to model the relationship between the outcome variable and an assignment variable in an RCT, there is no need to choose a bandwidth to reduce sensitivity to functional form assumptions. The data analyst might as well use all the data. In sum, neither of the two questions that motivate the current study are germane to an RCT.

However, the need to choose an optimal bandwidth to reduce sensitivity to functional form assumptions makes things very different for RD designs. The first step in bandwidth selection is to determine the order of the polynomial function that relates the assignment and outcome variables. We will assume that, after choosing an appropriate bandwidth, the data analyst will fit a local linear model consistent with recent advice in the literature (Cattaneo & Titiunik, 2022; Gelman & Imbens, 2019). For instance, Cattaneo and Titiunik (2022) state, “The order of the polynomial should always be low, to avoid overfitting and erratic behavior near the cutoff point. The default recommendation is p = 1 (local linear regression).”

The next step is to choose a kernel function that will be used to weight observations within the chosen bandwidth. Common choices are the triangular kernel, which places more weight on observations near the cut-point, and the uniform kernel, which weights all observations equally. Our study uses the triangular kernel as it is asymptotically optimal for point estimation when used in combination with one of the MSE optimal bandwidth selectors described in the next paragraph (Cattaneo & Titiunik, 2022).

Once the functional form and the kernel have been chosen, it is possible to choose the bandwidth. The earlier writings on bandwidth selection (Calonico et al., 2014; Imbens & Kalyanaraman, 2012) derived bandwidth selectors by approximating the unmodeled curvature in the polynomial function (to obtain a bias estimate) and then minimizing the asymptotic MSE of the LATE estimator. More recent writings have derived bandwidths that optimize criteria other than MSE, such as the coverage error of confidence intervals (Calonico et al., 2020). Our study focuses on point estimation and so considers only the Imbens and Kalyanaram (IK) and Calonico et al. (CCT) approaches. We also consider the older “rule of thumb” (ROT) bandwidth from the kernel density estimation literature (Silverman, 1986).

It is also necessary to decide whether the bandwidth selector will constrain the bandwidth to be the same on both sides of the cut-point or if separate bandwidths will be estimated on either side of the cut-point. Imbens and Kalyanaraman (2012) recommend a common bandwidth, but there is no consensus in the literature that a common bandwidth is preferred. Our study allows the bandwidth to vary on either side of the cut-point, as that is the default setting in the rdrobust R package used for estimation (Calonico et al., 2022).

To the best of our knowledge, all existing work on bandwidth selection has focused on choosing a bandwidth that is optimal (relative to some criterion) for estimation of the LATE. The bandwidth selector proposed in Calonico et al. (2014) has been updated to allow for the inclusion of covariates (Calonico et al., 2019). However, the target of estimation is still assumed to be the LATE, as opposed to a moderator of the LATE. As a result, given the current state of the methodological literature, it is necessary to proceed in an ad-hoc manner when choosing a bandwidth for moderator analysis.

One reasonable approach is to use an existing bandwidth selector for the LATE and to conduct moderator analysis by fitting the necessary statistical model within the bandwidth selected to estimate that LATE. A possible downside of this approach is that moderation effects are typically estimated with less precision than main effects (Aguinas, 1995), so the chosen bandwidths may be smaller than what is desirable to estimate a moderator effect. Exploring this possibility is the main motivation of the current article.

In the interacted model approach, once a bandwidth is chosen, the data analyst fits the model described in equation (6) using data within the bandwidth. This model, of course, includes a term for the assignment variable, which we denote X_k. We assume that this variable has been centered around the cut-point so that X_k > 0 implies $T_{k} = 1$ and X_k < 0 implies $T_{k} = 0$ (or vice versa if lower scores on the assignment condition result in $T_{k} = 1$ ):

E (Y_{k}) = β_{0}^{'} + β_{1}^{'} T_{k} + β_{2}^{'} M_{k} + β_{3}^{'} T_{k} M_{k} + β_{4}^{'} X_{k} + β_{5}^{'} T_{k} X_{k} + β_{6}^{'} M_{k} X_{k} + β_{7}^{'} T_{k} M_{k} X_{k}

(6)

The estimate of

β_{3}^{'},

{\hat{β}}_{3}^{'}

is the estimated moderating effect of

M_{k}

on the LATE. We refer to this moderating effect as the differential LATE (or DLATE) in the following. The test statistic

t_{i n t}^{'} = \frac{{\hat{β}}_{3}^{'}}{S E ({\hat{β}}_{3}^{'})}

(7)

can be compared to the relevant t-distribution to conduct a statistical test.

The analyst might be tempted to omit the interaction terms X_kM_k and X_kM_kT_k to simplify the model. However, omitting these terms will in general introduce bias when there is a correlation between M_k and the assignment variable. We note as an aside that to interpret the estimate of $β_{1}^{'}$ as the LATE for the $M_{k} = 0$ group it is necessary to center the moderator variable M_k around its mean (Calonico et al., 2019). We caution the reader that the moderating effect of $M_{k}$ on the LATE should also be interpreted locally. That is, the $T_{k} X_{k}$ , X_kM_k, and X_kM_kT_k terms are included to allow maximal flexibility when modeling the relationship between the assignment variable and the outcome variable. However, these terms should not be given a causal interpretation. For instance, we cannot interpret $β_{7}^{'}$ as the extent to which the moderating influence of $M_{k}$ varies as a function of the assignment variable.

The subgroup regression approach involves first obtaining the subset of the data for which M_k =0 and fitting the model in equation (8):

E (Y_{k}) = α_{0}^{'} + α_{1}^{'} T_{k} + α_{2}^{'} X_{k} + α_{3}^{'} T_{k} X_{k}

(8)

and then obtaining the subset of the data for which M_k= 1 and fitting the model in equation (9):

E (Y_{k}) = γ_{0}^{'} + γ_{1}^{'} T_{k} + γ_{2}^{'} X_{k} + γ_{3}^{'} T_{k} X_{k} .

(9)

Again

{\hat{γ}}_{1}^{'} - {\hat{α}}_{1}^{'}

estimates moderating effect of

M_{k}

and

t_{s u b}^{'} = \frac{{\hat{γ}}_{1}^{'} - {\hat{α}}_{1}^{'}}{\sqrt{V A R ({\hat{γ}}_{1}^{'}) + V A R ({\hat{α}}_{1}^{'})}}

(10)

can be used for inference. However, due to the different bandwidths that will be used and the correlation of the moderator with the treatment indicator variable,

t_{s u b}^{'}

and

t_{i n t}^{'}

will not be equivalent.

Greater flexibility in bandwidth choice is provided if the analyst adopts the subset regression approach (i.e., separate bandwidths can be chosen for different values of the moderating variable). While this allows for the choice of optimal bandwidths for the separate estimation of the LATE in the $M_{k} = 1$ and the $M_{k} = 0$ groups, this does not necessarily lead to an optimal estimate of the differential LATE (i.e., the moderating effect of $M_{k}$ ). Additionally, there may be downsides to this approach. A crucial difference between RCT and RD designs is that the moderator variable will not be correlated with treatment assignment in an RCT design, but it could potentially be highly correlated with treatment assignment in an RD design. In such situations, there may be very few cases with $M_{k} = 1$ (or $M_{k} = 0$ ) either above or below the cut-point. This will make it very difficult for a bandwidth selector to accurately estimate the curvature necessary to approximate the bias and choose a bandwidth. The problem would seem to be particularly acute when the bandwidth is not constrained to be the same on both sides of the cut-point. Indeed, our simulations suggest that this is the case.

Existing Literature on Treatment Effect Heterogeneity in RD Studies

This section reviews the limited existing literature on the topic of understanding treatment effect heterogeneity in RD designs. Cattaneo et al. (2022) discuss the use of covariates in RD studies generally, including a section on the use of covariates for understanding treatment effect heterogeneity. The paper mentions that either the subgroup regression or interacted regression approach can be chosen, but it does not discuss factors that would lead data analysts to choose one or the other approach. Reguly (2021) presents a machine learning approach similar to regression trees that can be used to explore heterogeneity with respect to a large number of possible covariates. Thus, the Reguly paper contemplates a situation very different from the current effort, where the goal is to explore moderation with respect to a single potential moderating variable or, at most, a handful of potential moderating variables.

Hsu and Shen (2019) propose a family of tests for three sorts of hypotheses relating to treatment effect heterogeneity: (1) that the intervention is beneficial for at least some subpopulations, (2) that the intervention has any impact for at least some subpopulations, and (3) that the intervention effect is heterogeneous across the defined subpopulations. The tests define the null and alternative hypotheses using moment inequalities, which can be estimated with local linear polynomials. However, the paper assumes a particular bandwidth has already been chosen and does not provide guidance regarding how to choose this bandwidth. Additionally, the Hsu and Shen (2019) approach is rarely used in applied work. Even papers that cite Hsu and Shen (2019) seem to choose instead to adopt either the subgroup regression or the interacted regression approach (Wasserman, 2023).

Bansak and Nowacki (2022) distinguish between studies of effect heterogeneity where the goal is simply to understand how the LATE systematically differs across levels of a third variable and where the goal is to understand if this variable causes the observed differences in the LATE. They term the first type of heterogeneity analysis Heterogeneity-in-Discontinuities (HiD) and the second type Moderation-in-Discontinuities (MiD). Importantly, the current article does not make this distinction. That is, it is agnostic as to whether the objective of the analysis is to estimate HiD or MiD effects and so the term moderator is used to apply to any analysis exploring variation in LATE effects with respect to a single binary variable. In particular, our work does not assume that interest is in causal moderation in the sense of Dong et al. (2023). Causal moderation requires additional assumptions that we do not adopt in this article. The Bansak and Nowacki (2022) paper discusses the use of either the subgroup regression or the interacted regression approach to estimate HiD effects. However, it does not discuss bandwidth choice or provide guidance regarding how best to choose between these two approaches.

Our informal impression suggests that most papers adopt the subgroup regression approach to produce separate estimates of the LATE for different subgroups (e.g., Card & Giuliano, 2016; Hansen, 2015). However, these papers do not conduct a formal test of the equivalence of the two estimates (i.e., they do not conduct moderation analysis). Nonetheless, some papers have adopted the interacted model approach (e.g., Barrow et al., 2020). Hsu and Shen (2019) confirm our impression, reporting results from a review of 15 papers in economics journals during the 2015 and 2016 publication years that explored treatment effect heterogeneity and used an RD design. Thirteen of the 15 papers used a subgroup regression approach whereas two of the 15 papers used an interacted regression approach. Regardless of whether the subgroup or interacted regression approach is chosen most existing papers appear to explore treatment effect heterogeneity by fitting models within arbitrarily chosen bandwidths. The current article will explore whether this practice could be improved by using previously developed bandwidth selectors.

Methods

We conducted a simulation study to gain insight into the two motivating questions of the article. Those questions are: (1) “Is it preferable to adopt a subset regression approach or an interacted regression approach for estimating the moderating effects of a binary variable?” and (2) “How should we determine the best bandwidth to use for each approach?” We first describe the data generation process for the different conditions that were explored. Next, we describe the different analytic approaches and the evaluation metrics used to understand the performance of these methods.

Data Generation Process

To simulate the assignment variable (X), we generated a random variate from a truncated normal distribution with mean 0, standard deviation 1, and truncation points −3 and 3. This allowed the assignment variable to function much like a normal random variable without needing to worry about extreme values. The cut-point was set at 0 for all simulation runs with intervention assigned to units with X values below 0.

Many applications of RD utilize large administrative datasets with sample sizes in the tens of thousands (e.g., Barrow et al., 2020; Hansen, 2015). However, in applied education research and evaluation work, primary data collection activities rarely yield datasets of this size (e.g., Baker et al., 2015; Louie et al., 2016). We chose to explore three different sample sizes with values of (n = 200, 400, and 600) based on each author's experience engaging in primary education evaluation data collection. We refrained from exploring large sample sizes as bandwidth selection algorithms have asymptotically optimal properties that can be expected to work well with large sample sizes but may not work as well when only smaller sample sizes are available.

We simulated potential outcomes (Rubin, 1974) under both the intervention (Y₁) and control (Y₀) conditions for each unit in the sample. We then created the observed outcome variable (Y) by retaining the intervention potential outcome (Y₁) for units with an X value below 0 and the control potential outcome (Y₀) for units with an X value above 0. This approach not only accurately simulates the assignment process in an RD study but also makes transparent the assumptions that are made when the functional relationship between the assignment variable and the outcome variable is different in the intervention and control conditions.

The equations for generating the potential outcomes are as follows:

Y_{0 k} = A_{0} X_{k} + B_{0} X_{k}^{2} + C_{0} X_{k}^{3} + ε_{0 k}

(11)

for the control potential outcome and

Y_{1 k} = A_{1} X_{k} + B_{1} X_{k}^{2} + C_{1} X_{k}^{3} + τ + τ_{m} M_{k} + A_{2} M_{k} X_{k} + B_{2} M_{k} X_{k}^{2} + C_{2} M_{k} X_{k}^{3} + ε_{1 k}

(12)

for the intervention potential outcome. The

ε_{0 k}

and

ε_{1 k}

parameters are independent truncated normal random variates with mean 0, standard deviation 1, and truncation points −3 and 3. The values of

A_{0}

B_{0},

C_{0}, A_{1}

B_{1}

C_{1}

A_{2}

B_{2}

, and

C_{2}

vary across simulation conditions as described in the subsequent subsection titled Functional form parameters.

M_{k}

is generated differently in different simulation conditions as described in the subsequent subsection titled Moderator generation. The values of τ and τ_m are set to 0.2 and 0.1, respectively. The value of 0.2 roughly corresponds to a “large” effect size for educational interventions using the benchmarks in Kraft (2020). In any event, the size of the overall LATE is not expected to influence the relative performance of the different analytic approaches that are explored.

Little is known about the extent to which LATE values systematically vary as a function of measured, binary moderator values. It seems likely that DLATEs are generally smaller than LATEs. A value of 0.1 was chosen in acknowledgement of this intuition while still ensuring a large enough moderating effect to emphasize possible differences between the analytic approaches considered.

Functional Form Parameters

The values of $A_{0}$ , $B_{0},$ $C_{0}, A_{1}$ , $B_{1}$ , $C_{1}$ , $A_{2}$ , $B_{2}$ , and $C_{2}$ govern the functional relationship between the assignment variable, the outcome variable, and the moderating variable. When $A_{2} = B_{2} = C_{2} = 0$ , we say that this functional relationship is moderator independent (MI). When, additionally, $A_{0} = A_{1}$ , $B_{0} = B_{1}$ , and $C_{0} = C_{1}$ we say that the functional relationship is treatment independent (TI). When, $B_{0} = B_{1} = 0$ and $C_{0} = C_{1} = 0$ , we say the functional relationship is simple (i.e., linear). When any of $B_{0}$ , $B_{1}, C_{0}, C_{1} \neq 0$ , we say the functional relationship is complex (i.e., nonlinear, specifically cubic). Investigating functional relationships of various levels of complexity is important because the primary rationale for restricting estimation to a particular bandwidth is to reduce sensitivity of estimates to misspecification of the functional form governing the relationship between the assignment variable and the outcome variable. Accordingly, we constructed simulation conditions with different levels of complexity as described in Table 1. Figure 1 plots the implied functional relationship between the assignment variable and the outcome variable for the simple, MI, TI (solid, black line) and complex, MI, TI (dotted line) and simple MI, not-TI (dashed line) cases (the dash-dot line represents the relationship in the treatment group with the solid line representing the relationship in the control group). The complex, MI, not-TI relationship looks like the dashed line for the treatment group and the solid line for the control group. The simple, not-MI, not-TI relationship looks like the simple MI, not-TI relationship but with a steeper line in the M = 1, T = 1 group. The complex, not-MI, not TI relationship looks like the complex, MI, not-TI relationship but with a shallower curve in the M = 1, T = 1 group.

Table 1.

Simulation Settings Governing Relationship Between Assignment Variable, Outcome, and Moderator.

Condition	Parameter values
Simple, MI, TI	$A_{0} = A_{1} = 1$ ; $A_{2} = B_{0} = B_{1} = B_{2} = C_{0} = C_{1} = C_{2} = 0$
Complex, MI, TI	$A_{0} = A_{1} = 1$ ; $B_{0} = B_{1} = 0.08$ ; $C_{0} = C_{1} = - 0.06$ ; $A_{2} = B_{2} = C_{2} = 0$
Simple, MI, not-TI	$A_{0} = 1; A_{1} = 2$ ; $A_{2} = B_{0} = B_{1} = B_{2} = C_{0} = C_{1} = C_{2} = 0$
Complex, MI, not-TI	$A_{0} = 1; A_{1} = 2$ ; $B_{0} = 0; B_{1} = 0.08$ ; $C_{0} = 0; C_{1} = - 0.06$ ; $A_{2} = B_{2} = C_{2} = 0$
Simple, not-MI, not-TI	$A_{0} = 1; A_{1} = 2; A_{2} = 1$ $B_{0} = B_{1} = B_{2} = C_{0} = C_{1} = C_{2} = 0$
Complex, not-MI, not-TI	$A_{0} = 1; A_{1} = 2$ ; $A_{2} = 1; B_{0} = 0; B_{1} = 0.08$ ; $B_{2} = 0.04; C_{0} = 0; C_{1} = - 0.06$ ; $C_{2} = - 0.03$

Note. MI = moderator independent; TI = treatment independent.

Moderator Generation

An important difference between moderation analysis in RD as compared with RCT studies is the possibility of an imbalance in the distribution of the moderator between the intervention and the control groups in an RD study due to the correlation between the moderator variable and the assignment variable (randomization ensures that such an imbalance cannot occur in an RCT). Additionally, the overall prevalence of the moderating variable may be a factor in determining the appropriate bandwidth, particularly when the subgroup regression approach is used since in that case a very low (or high) prevalence will result in a small sample size in one of the subgroups. Accordingly, we generated moderator values with three levels of overall prevalence corresponding to Pr( $M_{k} = 1$ ) equal to 0.1, 0.3, or 0.5 and also with three levels of correlation between the moderator value and the running variable, namely Cor $(M_{k}, X_{k})$ equal to 0, −0.2, and −0.5. These values were chosen to represent the range of values that might be observed in practice. For instance, an educational evaluation might find that the proportion of English Language Learner students in the study is around 0.1, the proportion of students eligible for FRL is around 0.3 and the proportion of females in the study is around 0.5. Similarly, certain possible moderating variables may be completely unrelated to the assignment variable (represented by a correlation of 0) whereas a moderator variable like FRL might have a strong negative correlation with an assignment variable measuring achievement, possibly as large as −0.5. The desired prevalence and correlation values were obtained by randomly generating Bernoulli random variables with $l o g i t (P r (M_{k} = 1)) = η_{0} + η_{1} X_{k}$ . The values of $η_{0}$ and $η_{1}$ were adjusted to result in the desired combination of prevalence and correlation values.

Overall, we consider three different sample sizes, six different functional form types, and nine combinations of moderator prevalence and correlation, resulting in 3*9*6 = 162 total simulation conditions. For each condition, we simulate 1,000 datasets and apply the analytic procedures described in the next subsection to each dataset.

Estimation Procedures and Evaluation Metrics

In keeping with the goals of the article, for each simulated dataset, we estimate the DLATE using several analytic approaches. The first dimension we consider is regression type with levels subgroup regression and interacted regression. For each regression type, we consider three different bandwidth selectors.

The first selector is Silverman's (1986) “ROT” bandwidth from the kernel density estimation literature. We refer to this as the ROT selector. This bandwidth selector was not developed specifically for application to RD studies. However, it was in wide use for RD treatment effect estimation before publication of the Imbens and Kalyanaraman (2012) paper (Lee & Lemieux, 2010). Additionally, because the MSE-optimal selectors are estimating bias and variance for the LATE, not the DLATE, the more general ROT approach may end up having better properties. Finally, because the ROT selector does not require an estimate of the asymptotic bias or variance, it will be possible to obtain values of the ROT bandwidth even when the data do not support estimation of one of the MSE-optimal bandwidths. Because it is common in the literature to explore sensitivity to different bandwidth choices, we also utilize 50% of the estimated optimal bandwidth and 150% of the estimated optimal bandwidth (see Wasserman, 2023, for an example of this approach).

The second selector is an MSE-optimal selector with a “first-generation” (Cattaneo & Titiunik, 2022) plug-in rule as proposed by Imbens and Kalyanaraman (2012) and implemented in the R package rdrobust (Calonico et al., 2020). We refer to this as the IK selector. It is considered a first-generation approach because unknown quantities in the bandwidth formula are replaced with inconsistent estimators (Cattaneo & Titiunik, 2022). In our experience, the IK selector tends to choose slightly larger bandwidths than the CCT selector described in the next paragraph. This may be a desirable property in the context of moderation analysis where the DLATE will usually have a larger standard error than the LATE. We again explore sensitivity to different bandwidths by utilizing 50% of the estimated optimal bandwidth and 150% of the estimated IK bandwidth.

The third selector is an MSE-optimal selector using a “second-generation” (Cattaneo & Titiunik, 2022) plug-in rule described in Calonico et al. (2014) and implemented in the R package rdrobust (Calonico et al., 2020). It is considered a second-generation approach because unknown quantities in the bandwidth formula are replaced by consistent estimators. While the second-generation approach can be expected to perform better (at least asymptotically) for estimation of the LATE its optimality with respect to the DLATE is unknown. Indeed, because it tends to choose smaller bandwidths than the ROT and IK approaches it may not work well for DLATE estimation. We again explore sensitivity to different bandwidths by utilizing 50% of the estimated optimal bandwidth and 150% of the estimated optimal bandwidth.

Finally, we estimate the DLATE using a global regression that does not use any bandwidth selection at all. Given the imprecision with which moderator effects are estimated global regression approaches will be relatively advantaged. Thus, for each regression type we estimate the DLATE within 3*3 + 1 = 10 different bandwidths.

We evaluate the different approaches by computing their bias and root MSE (RMSE) relative to the true value of DLATE (0.1). However, as described in the subsequent section, for certain simulation runs it was impossible for a given bandwidth selector to estimate the necessary bandwidth. Even if a bandwidth could be identified, sometimes the resulting regression model had a singular fit and the DLATE could not be estimated. Therefore, let n_val represents the number of simulation runs where an estimate of DLATE could be obtained (with n_val $\leq$ 1,000). We estimate the bias for a particular combination of data generating process and estimation procedure by equation (14):

B i a s = \frac{\sum (E s t . - D L A T E)}{n v a l}

(14)

Furthermore, the RMSE for a particular combination of data generating process and estimation procedure is given by equation (15):

R M S E = \frac{\sum {(E s t . - D L A T E)}^{2}}{n v a l} .

(15)

Results

Tables 2–4 present the main results of the simulation study. Because the subgroup and interacted regression approaches give equivalent estimates when a global regression is used (Bansak & Nowacki, 2022) we present results only for interacted models for the global regression case. The tables are a cross-classification of analysis type (rows) and data generating process (columns). Graphical representations of these results can be found in the supplemental material. Additionally, Supplemental Material Tables SM7 and SM8 present p-values and effect sizes associated with the 10 largest effects from three-way analysis of variance (ANOVA) models estimating the relationship between the simulation design factors and bias and RMSE, respectively. Note that the analytic approach was a factor in 18 of the 20 largest effects listed across the two tables and in eight of 10 and all 10 of the largest effects from the bias and RMSE ANOVA results, respectively. Table 2 presents average bias and RMSE values as a function of analysis type and functional form (averaging across all other simulation settings). Table 3 presents average bias and RMSE values as a function of analysis type, moderator prevalence, and moderator correlation with the assignment variable (averaging across all other simulation settings). Table 4 presents average bias and RMSE values as a function of analysis type and sample size (averaging across all other simulation settings). Because RMSE and bias values were unduly influenced by outlying results, Tables 2–4 are based only on the 99% of simulations with the smallest differences between the estimated and true DLATE values. For completeness, Tables SM4–SM6 in the supplemental material present results without trimming.

Table 2.

Bias and RMSE Varying Functional Relationship Between Assignment Variable, Moderator, and Outcome Variable; Trimmed.

Approach	% of BW	RMSE						Bias
Approach	% of BW	Comp, MI, TI	Comp, MI, not-TI	Comp, not-MI, not-TI	Simp, MI, TI	Simp, MI, not-TI	Simp, not-MI, not-TI	Comp, MI, TI	Comp, MI, not-TI	Comp, not-MI, not-TI	Simp, MI, TI	Simp, MI, not-TI	Simp, not-MI, not-TI
CCT int	0.5	5.62	5.90	8.95	6.20	5.61	6.56	−0.22	−0.22	−0.23	−0.16	−0.20	−0.17
	1	3.49	3.35	3.44	3.56	3.43	3.37	−0.14	−0.15	−0.15	−0.19	−0.17	−0.18
	1.5	1.29	1.30	1.35	1.30	1.27	1.31	−0.07	−0.07	−0.08	−0.07	−0.07	−0.07
CCT sub	0.5	2.56	2.55	2.59	2.53	2.50	2.50	−0.38	−0.38	−0.37	−0.33	−0.36	−0.33
	1	1.16	1.16	1.17	1.19	1.18	1.21	−0.02	−0.02	−0.03	−0.08	−0.08	−0.08
	1.5	1.18	1.06	1.02	1.05	1.06	1.06	−0.04	−0.05	−0.06	−0.03	−0.03	−0.04
IK int	0.5	3.80	2.64	2.86	3.84	2.87	3.79	−0.19	−0.13	−0.16	−0.17	−0.15	−0.14
	1	1.71	1.74	1.75	1.75	1.78	1.77	−0.08	−0.09	−0.10	−0.12	−0.11	−0.12
	1.5	1.02	1.03	1.06	1.04	1.00	1.05	0.01	0.01	0.00	0.01	0.01	0.01
IK sub	0.5	1.75	1.66	1.90	1.60	1.83	1.92	−0.11	−0.05	−0.13	−0.01	−0.10	−0.10
	1	0.81	0.80	0.78	0.83	0.85	0.79	−0.01	−0.02	0.00	−0.04	−0.06	−0.02
	1.5	3.91	3.91	1.12	3.90	2.73	1.10	−0.02	−0.02	−0.05	−0.02	−0.01	−0.03
ROT int	0.5	1.17	1.65	1.90	0.98	1.62	3.63	−0.04	−0.02	−0.04	−0.06	−0.03	−0.14
	1	0.67	0.69	0.67	0.70	0.72	0.68	0.03	0.02	0.00	0.01	0.02	0.01
	1.5	0.64	0.67	0.73	0.64	0.67	0.72	0.02	0.01	0.00	0.02	0.03	0.03
ROT sub	0.5	3.19	1.87	2.10	3.61	2.94	2.00	−0.14	−0.04	−0.06	−0.14	−0.10	−0.06
	1	0.82	0.84	0.83	0.83	0.85	0.84	0.01	0.01	0.00	0.00	0.00	0.00
	1.5	0.79	0.82	0.89	0.77	0.79	0.89	0.01	0.01	0.00	0.02	0.01	0.01
Global int	n/a	0.64	0.64	0.65	0.66	0.65	0.59	−0.01	−0.04	−0.14	0.01	0.02	0.01

Note. Results are based on the 99% of simulations with the smallest differences between the estimated and true differential local average treatment effect (DLATE) values. “Int” refers to an interacted regression fit within a single bandwidth. “Sub” refers to separate regressions fit within separate bandwidths for subsets of the data defined by the values of the moderator variable. BW = bandwidth; MI = moderator independent; TI = treatment independent; CCT = Calonico et al. (2020) approach; IK = Imbens and Kalyanaraman (2012) approach; ROT = Silverman's (1986) rule-of-thumb bandwidth.

Table 3.

Bias and RMSE Varying Moderator Prevalence and Correlation with Assignment Variable; Trimmed.

Approach	% of BW	RMSE						Bias
		Moderator prevalence			Moderator correlation with assignment variable			Moderator prevalence			Moderator correlation with assignment variable
		0.5	0.3	0.1	0	−0.2	−0.5	0.5	0.3	0.1	0	−0.2	−0.5
CCT int	0.5	1.02	4.33	20.26	3.01	3.31	15.50	0.01	−0.19	−0.61	−0.15	−0.26	−0.19
	1	0.50	0.53	11.48	1.81	5.63	2.70	0.01	−0.03	−0.60	−0.15	−0.33	0.04
	1.5	0.35	0.39	3.78	0.59	1.08	2.57	0.02	0.01	−0.30	−0.02	−0.21	0.04
CCT sub	0.5	1.42	2.03	7.10	2.97	2.39	2.13	−0.03	−0.72	−0.38	−0.19	−0.14	−0.90
	1	0.63	0.94	2.87	1.22	1.22	1.06	0.02	−0.06	−0.19	−0.03	−0.11	0.00
	1.5	0.44	0.53	3.40	1.23	1.23	0.60	0.02	0.01	−0.27	−0.08	−0.04	0.01
IK int	0.5	0.84	1.11	10.70	2.71	2.62	5.05	0.03	−0.01	−0.69	−0.13	−0.28	−0.02
	1	0.42	0.47	5.27	0.82	2.20	2.42	0.03	−0.01	−0.41	−0.06	−0.27	0.07
	1.5	0.30	0.33	2.94	0.41	0.69	2.34	0.03	0.02	−0.03	0.00	−0.02	0.06
IK sub	0.5	0.83	1.14	4.54	1.84	1.93	1.47	0.03	−0.06	−0.33	−0.07	−0.17	0.01
	1	0.47	0.51	1.80	0.79	0.86	0.77	0.03	−0.01	−0.14	−0.06	−0.02	0.01
	1.5	0.32	0.36	10.04	0.60	0.47	9.20	0.03	0.02	−0.17	−0.02	0.01	−0.07
ROT int	0.5	0.37	0.41	5.69	0.72	2.01	3.08	0.03	−0.01	−0.23	−0.03	−0.18	0.07
	1	0.20	0.25	1.91	0.29	0.34	1.69	0.02	−0.01	0.03	0.00	−0.01	0.06
	1.5	0.20	0.24	1.88	0.29	0.34	1.65	0.02	0.01	0.03	0.01	0.00	0.06
ROT sub	0.5	0.37	0.42	8.88	1.52	3.16	3.41	0.03	−0.01	−0.37	−0.09	−0.20	0.06
	1	0.19	0.25	2.46	0.31	0.67	1.76	0.02	0.00	0.00	−0.01	−0.04	0.07
	1.5	0.19	0.24	2.43	0.31	0.67	1.73	0.02	0.01	0.00	0.00	−0.03	0.08
Global int	n/a	0.16	0.20	1.84	0.22	0.28	1.67	−0.02	−0.05	−0.01	−0.01	−0.04	−0.03

Note. Results are based on the 99% of simulations with the smallest differences between the estimated and true DLATE values. “Int” refers to an interacted regression fit within a single bandwidth. “Sub” refers to separate regressions fit within separate bandwidths for subsets of the data defined by the values of the moderator variable. BW = bandwidth; CCT = Calonico et al. (2020) approach; IK = Imbens and Kalyanaraman (2012) approach; ROT = Silverman's (1986) rule-of-thumb bandwidth.

Table 4.

Bias and RMSE Varying Sample Size; Trimmed.

Approach	% of BW	RMSE			Bias
		Sample size			Sample size
		200	400	600	200	400	600
CCT int	0.5	8.34	2.29	8.94	−0.33	−0.14	−0.14
	1	7.84	1.21	1.47	−0.52	0.05	−0.04
	1.5	1.57	1.56	0.80	−0.24	0.01	0.00
CCT sub	0.5	3.48	2.83	1.58	−1.08	−0.11	−0.05
	1	1.77	1.20	0.69	−0.02	−0.06	−0.07
	1.5	1.80	1.11	0.43	−0.08	−0.02	−0.03
IK int	0.5	4.61	1.97	3.41	−0.18	−0.21	−0.08
	1	2.93	1.27	1.10	−0.32	0.04	−0.04
	1.5	0.92	1.43	0.75	0.00	0.01	0.00
IK sub	0.5	3.03	1.39	1.04	−0.24	0.04	−0.07
	1	1.27	0.64	0.55	−0.04	0.02	−0.06
	1.5	0.87	0.42	6.97	0.00	0.00	−0.08
ROT int	0.5	2.61	1.57	1.32	−0.15	0.01	−0.02
	1	0.50	0.92	0.65	0.04	0.01	−0.02
	1.5	0.50	0.91	0.62	0.04	0.01	0.00
ROT sub	0.5	4.80	1.17	2.02	−0.27	0.03	−0.04
	1	0.92	0.95	0.65	0.02	0.02	−0.02
	1.5	0.92	0.94	0.63	0.02	0.01	0.00
Global int	n/a	0.41	0.89	0.61	−0.01	0.00	−0.07

Note. Results are based on the 99% of simulations with the smallest differences between the estimated and true differential local average treatment effect (DLATE) values. “Int” refers to an interacted regression fit within a single bandwidth. “Sub” refers to separate regressions fit within separate bandwidths for subsets of the data defined by the values of the moderator variable. BW = bandwidth; CCT = Calonico et al. (2020) approach; IK = Imbens and Kalyanaraman (2012) approach; ROT = Silverman's (1986) rule-of-thumb bandwidth.

Before discussing the results in Tables 2–4, it is necessary to reflect on the implications of the fact that for certain combinations of estimation procedure and simulation run either: (a) it was impossible to obtain an estimate of the bandwidth (denoted BNE for “bandwidth not estimable”) or (b) the model fit was singular (denoted SMF for “singular model fit”). Supplemental Tables SM1–SM3 present the rate of occurrence of each of these problems as a function of analysis type and functional form (SM1); analysis type, moderator prevalence, and moderator correlation with assignment variable (SM2); and analysis type and sample size (SM3).

Table SM1 suggests that the rates of BNE and SMF are not impacted by the functional form (i.e., for a given row, entries are constant across columns). Table SM1 also shows that BNE problems occur only when the subset regression approach is employed. This occurs because there may be a very small sample size in one of the subsets when the moderator prevalence is low and the moderator is correlated with the assignment variable. Notably (and unsurprisingly), the ROT estimation approach does not have any BNE problems because there is no need to estimate a statistical model to compute the bandwidth estimate. On the other hand, SMF problems occur more frequently when interacted models are used. There are likely two reasons for this. The first is that datasets where the bandwidth is hard to estimate are also likely to be datasets where there is a strong correlation between the moderator variable and the assignment variable, increasing the likelihood of an SMF. Since models are never estimated on datasets with BNE problems, there are fewer problematic datasets contributing to the SMF rates when subset regression is used. The second reason for fewer SMF problems in subset regression models is due to the restriction of range of the assignment variable, thus decreasing its correlation with other variables in the model and decreasing the likelihood of an SMF.

Perhaps the most important takeaway from Table SM1 is that if we combine the BNE and SMF rates, subset regression using a ROT bandwidth performs best, with a 4% rate of problems with model estimation. All other analysis procedures have combined rates of 6% or more.

Table SM2 provides further understanding of situations when BNE and SMF problems occur. Almost all problems occur when the moderator prevalence is only 0.1. Additionally, if there is no correlation between the moderator and the assignment variable then problems occur only when a very small bandwidth is used (50% of optimal). Slightly more problems pop up with the correlation between the moderator and assignment variable is −0.2. However, most problems occur only when the correlation between the moderator and assignment variable is −0.5.

Table SM3 illustrates that, while there are certainly more problems at smaller sample sizes, even with a sample size of 600 the combined BNE and SMF error rate is still non-negligible. It ranges from 0.03 for the ROT subset approach to 0.09 for the very small bandwidths associated using 50% of the CCT or IK optimal bandwidths.

Given the different error rates for different estimation approaches, the results in Tables 2–4 must be interpreted with caution. It is possible that approaches that seem to have the best bias and RMSE properties have been unduly advantaged because datasets that are hard to estimate have been eliminated before bias and RMSE are computed. Luckily, the average error rates across estimation approaches do not differ much. However, the differences that do exist may still be large enough to influence the results. With this caveat in mind, we turn attention to the bias and RMSE results in Tables 2–4.

Table 2 suggests that, from an RMSE perspective, global regression is the preferred approach. Global regression also performs best when considering bias, except when there is a complex functional form that depends on treatment assignment. In this case, using the ROT bandwidth results in lower average bias. Interestingly, when bandwidth selectors are used, bias seems to be lowest when 150% of the chosen bandwidth is utilized. However, from an RMSE perspective, 150% of the chosen bandwidth is best for both types of CCT selectors and for interacted models using the other selectors. When using subgroup regression, 100% of the chosen bandwidth is preferred for the IK and ROT selectors. Furthermore, the subgroup regression approach generally appears preferable to the interacted approach when a bandwidth selector is used. However, the general takeaway is that global regression is almost always preferred.

Global regression also emerges as the preferred approach when looking at Table 3. It has the lowest RMSE for all values of moderator prevalence and moderator correlation with the assignment variable, except for a correlation of −0.5. In this case, a subgroup regression approach using either 150% of the CCT bandwidth or 100% of the IK bandwidth appears to be preferred. However, these approaches also had high problem rates at a correlation of −0.5, so that fact may be skewing results. Global regression also has uniformly low average bias values, generally lower than any other estimation approach. However, for some simulation settings, the ROT approach appears preferable (either using interacted models or subset regression).

Global regression continues to perform well when examining Table 4. It is the best approach from both a bias and RMSE perspective at a sample size of 200. However, at sample sizes of 400 and 600, 100% of the IK bandwidth selector using subgroup regression has lower RMSE values (and at a sample size of 400, 150% of the IK bandwidth has even lower RMSE). Any of the approaches based on the ROT bandwidth appear preferable from a bias perspective when the sample size is 600.

Empirical Illustration

Having reviewed the results from the simulation study, we now apply these analytic approaches to data from a recently completed evaluation of an educational intervention. The intervention involved providing additional services to support struggling 9^th-grade readers. The data available from the evaluation consisted of 837 students from diverse high schools in four states. Students falling below a threshold on a pretreatment reading measure were assigned to receive the intervention (n = 287) with students scoring below the threshold receiving business as usual literacy instruction (n = 550). The study made use of the STAR Reading assessment (Renaissance Learning, 2009) as the primary reading achievement outcome of interest. For the purposes of this illustration, we focus on the moderating effects of the variable, minority, which equals 1 if the student in question is a member of a minority racial/ethnic group (n = 611; 73%) and equals 0 otherwise (n = 226, 27%).

Table 5 describes the estimated values of the DLATE. The seven analytic approaches used were the CCT bandwidth selector (for both an interacted model and subgroup regressions), the IK bandwidth selector (for both an interacted model and subgroup regressions), the ROT bandwidth selector (for both an interacted model and subgroup regressions) and a global regression with an interacted model. The results illustrate just how consequential the choice of bandwidth can be, with estimates ranging from as low as 6.9 to as high as 59.1.

Table 5.

Seven Different Approaches to Estimating the Moderating Effect of Minority status.

Approach	n Left	n Right	BW Left	BW Right	Estimate	Lower CI	Upper CI
CCT Int	155	147	0.74	0.44	42.6	−214.0	299.2
CT Sub	(30, 59)	(79, 101)	(0.57, 0.34)	(0.78, 0.43)	13.3	−132.9	159.5
IK Int	236	195	1.33	0.56	44.8	−160.9	250.4
IK Sub	(39, 71)	(120, 130)	(0.75, 0.75)	(0.68, 0.51)	46.1	−86.9	179.1
RoT Int	271	497	1.97	1.91	12.9	−153.1	178.9
RoT Sub	(67, 211)	(140, 371)	(1.77, 2.23)	(1.77, 2.23)	6.9	−100.8	114.7
Global	287	550	−4.19	3.05	59.1	−85.3	203.4

Note. “Int” refers to an interacted regression fit within a single bandwidth. “Sub” refers to separate regressions fit within separate bandwidths for subsets of the data defined by the values of the moderator variable. BW = bandwidth; CCT = Calonico et al. (2020) approach; IK = Imbens and Kalyanaraman (2012) approach; ROT = Silverman's (1986) rule-of-thumb bandwidth. When subgroup regression is used the (x,y) notation denotes that the sample size (bandwidth) for the nonminority group is x and the sample size (bandwidth) for the minority group is y.

An additional reason for the wide range of estimates is the imprecision with which the DLATE is estimated (the study was substantially underpowered due to attrition). Confidence interval lower limits range from −214 to −85 and upper limits range from 114 to 299. Our simulation study suggests that global regression is generally the optimal approach. That would imply putting most of our trust in the 59.1 estimate that comes from the global regression approach. However, the available sample size was slightly larger than the maximum sample size we used in our simulations, and Table 5 suggests that it might be optimal to use the IK bandwidth selector with subset regressions in this case. Interestingly, the global and IK subset estimates are the two largest among the seven procedures explored. This suggests that the larger estimates may be the most trustworthy in this particular application.

Discussion

As noted above, current best practice in RD studies suggest that a data-driven bandwidth selector should be used to pick the subset of data for estimating the LATE. Using data only within a certain bandwidth of the cut point improves the robustness of the resulting estimate to model misspecification. However, existing bandwidth selectors have been derived to minimize an approximation to the MSE of the estimate of the LATE. There do not exist any bandwidth selectors that are optimized to minimize the MSE of the moderating influence of a third variable on the LATE.

Given this gap in the methodological literature an evaluator seeking to conduct a moderator analysis has two choices: (a) conduct the moderator analysis using all of the data (i.e., perform a global regression) or (b) choose a bandwidth and analytic approach in an ad-hoc manner. The current article reports on a simulation study that was conducted to help evaluators make this choice when the moderator in question is binary. The ad-hoc approaches considered in the simulation study involve using bandwidth selectors that have been suggested for estimation of the LATE. Once the bandwidth is identified, the moderator effect is estimated within the chosen bandwidth.

One important result of the study is information pertaining to when it will even be possible to choose an MSE-optimal bandwidth and/or estimate a nonsingular model to assess moderation. When the moderator prevalence is far from 0.5 or when the correlation of the moderator with the assignment variable is high (in absolute value) evaluators may have problems implementing certain approaches unless the sample size is large. Our simulation study suggests that in these cases using a ROT bandwidth with separate regressions for each value of the moderating variable is the approach most likely to avoid these problems. However, as explained below, for the sample sizes and functional forms considered in this simulation study, it is best to not estimate a bandwidth at all and instead to use a global regression.

Bandwidth selectors seek to optimize a trade-off between bias and variance that exists when choosing a subset of the data to use to estimate an RD treatment effect. A smaller bandwidth will tend to result in less bias when the model is mis-specified. However, choosing a smaller bandwidth will also inflate the variance. It is well known that sampling variance of estimates of moderator effects tends to be higher than the variance of estimates of main effects (Aguinas, 1995; Spybrook et al., 2016). Thus, we might suspect that existing bandwidth selectors that are optimized to minimize the MSE for a main effect will tend to choose bandwidths that are too small when they are utilized to estimate a moderating effect.

The current simulation study confirms this suspicion. Results not presented here due to space considerations demonstrate that the CCT approach tends to choose smaller bandwidths than the IK approach, which tends to choose smaller bandwidths than the ROT approach. Of course, a global regression represents the maximum possible bandwidth by using all of the data. The estimated mean squared of these approaches follows closely the inverse ordering. Estimates using the CCT approach have larger MSE than estimates using the IK approach. Estimates using the ROT approach have smaller MSE than estimates using the IK approach, but larger MSE than estimates using a global regression approach. In other words, the larger the bandwidth the smaller the MSE.

As a result, the current study suggests that when the sample size is small (n = 600 or smaller), existing bandwidth selectors should not be used when conducting a moderation analysis. Instead, a global regression approach should be used. This advice holds regardless of whether the relationship between the assignment variable and the outcome is linear or cubic and regardless of whether that relationship depends on the value of the treatment and/or moderator variable. It also holds across different levels of moderator prevalence and different levels of correlation between the moderator and the assignment variable. However, an important caveat is that this finding only holds if the higher rate of singular model fit for the global regression approach has not impacted the RMSE and bias results found in our study. It is hard to know whether or not this might be the case. It is also important to note that singular model fits were only a problem for certain simulation settings. So, for instance, if the moderator prevalence is between 0.3 and 0.7, or if the moderator variable is uncorrelated with the assignment variable, our simulations unambiguously suggest that global regression should be preferred (for the particular functional form and sample sizes that were explored).

Often evaluators will want to estimate models with continuous moderators and/or multiple moderating variables. While our simulation study cannot directly speak to these issues, it is reasonable to assume that existing bandwidth selectors will also choose bandwidths that are too small for these situations. Indeed, if there are multiple correlated moderators, statistical power will be very low for each one and the problem will be even more pronounced.

An important limitation to our findings is that the cubic function we used to represent a complex functional form is relatively linear within a wide range on either side of the cut point. So, the level of model misspecification was not that extreme. It is possible that the bandwidth selectors would perform better relative to global regression if the curvature of the regression function in the vicinity of the cut point was more extreme.

In the future, it may be possible to find an MSE-optimal bandwidth selector specifically attuned to moderation analysis. Until then, we caution evaluators against using bandwidth selectors meant for main effects when conducting a moderation analysis. At least for the small sample sizes and cubic functional form considered in this article, it is likely preferable to use all the data rather than a data-driven bandwidth.

Supplemental Material

sj-docx-1-aje-10.1177_10982140241260384 - Supplemental material for Analysis of Regression Discontinuity Designs with a Binary Moderating Variable

Supplemental material, sj-docx-1-aje-10.1177_10982140241260384 for Analysis of Regression Discontinuity Designs with a Binary Moderating Variable by Jason A. Schoeneberger and Christopher Rhoads in American Journal of Evaluation

Footnotes

Declaration of Conflicting Interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by the Institute of Education Sciences, (grant number R305A180154).

ORCID iDs

Jason A. Schoeneberger

Christopher Rhoads

Supplemental Material

Supplemental material for this article is available online.

Notes

References

Aguinas

(1995). Statistical power problems with moderated multiple regression in management research. Journal of Management, 21(6), 1141–1158. https://doi.org/10.1177/014920639502100607

Aiken

L. S.

West

S. G.

Schwalm

D. E.

Carroll

Hsuing

(1998). Comparison of a randomized and two quasi-experimental designs in a single outcome evaluation: Efficacy of a university-level remedial writing program. Evaluation Review, 22(2), 207–244. https://doi.org/10.1177/0193841X9802200203

Athey

Imbens

(2016). Recursive partitioning for heterogeneous causal effects. Proceedings of the National Academy of Sciences, 113(27), 7353–7360. https://doi.org/10.1073/pnas.1510489113

Baker

Smolkowski

Chaparro

Smith

J. L. M.

Fine

(2015). Using regression discontinuity to test the impact of a tier 2 Reading intervention in first grade. Journal of Research on Educational Effectiveness, 8(2), 218–244. https://doi.org/10.1080/19345747.2014.909548

Bansak

(2021). Estimating causal moderation effects with randomized treatments and non-randomized moderators. Journal of the Royal Statistical Society Series A: Statistics in Society, 184(1), 65–86. https://doi.org/10.1111/rssa.12614

Bansak

Nowacki

(2022). Effect heterogeneity and causal attribution in regression discontinuity designs [Paper]. https://doi.org/10.31235/osf.io/vj34m

Baron

R. M.

Kenny

D. A.

(1986). The moderator-mediator variable distinction in social psychological research: Conceptual, strategic, and statistical considerations. Journal of Personality and Social Psychology, 5(6), 1173–1182. https://doi.org/10.1037/0022-3514.51.6.1173

Barrow

Sartain

de la Torre

(2020). Increasing access to selective high schools through place-based affirmative action: Unintended consequences. American Economic Journal: Applied Economics, 12(4), 135–163. https://doi.org/10.1257/app.20170599

Berk

R. A.

Barnes

Ahlman

Kurtz

(2010). When second best is good enough: A comparison between a true experiment and a regression discontinuity quasi-experiment. Journal of Experimental Criminology, 6, 191–208. https://doi.org/10.1007/s11292-010-9095-3

10.

Bloom

H. S.

(2003). Using “short” interrupted time-series analysis to measure the impacts of whole-school reforms: With applications to a study of accelerated schools. Evaluation Review, 27(1), 3–49. https://doi.org/10.1177/0193841X02239017

11.

Bloom

H. S.

(2005). Randomizing groups to evaluate place-based programs. In Bloom

H. S.

(Ed.), Learning more from social experiments (pp. 115–172). New York, NY: Russell Sage Foundation.

12.

Bryan

C. J.

Tipton

Yeager

D. S.

(2021). Behavioural science is unlikely to change the world without a heterogeneity revolution. Nature Human Behaviour, 5, 980–989. https://doi.org/10.1038/s41562-021-01143-3

13.

Calonico

Cattaneo

M. D.

Farrell

M. H.

(2020). Optimal bandwidth choice for robust bias corrected inference in regression discontinuity designs. Econometrics Journal, 23(2), 192–210. https://doi.org/10.1093/ectj/utz022

14.

Calonico

Cattaneo

Farrell

Titiunik

(2019). Regression discontinuity designs using covariates. The Review of Economics and Statistics, 101(3), 442–451. https://doi.org/10.1162/rest_a_00760

15.

Calonico

Cattaneo

M. D.

Farrell

M. H.

Titiunik

(2022). Rdrobust: Robust data-driven statistical inference in regression-discontinuity designs (R Package version 0.99.8) [Paper]. https://cran.r-project.org/web/packages/rdrobust/rdrobust.pdf

16.

Calonico

Cattaneo

Titiunik

(2014). Robust nonparametric confidence intervals for regression-discontinuity designs. Econometrica, 82(6), 2295–2326. https://doi.org/10.1038/s41562-021-01143-3

17.

Card

Giuliano

(2016). Can tracking raise the test scores of high-ability minority students? American Economic Review, 106(10), 2783–2816. https://doi.org/10.1257/aer.20150484

18.

Cattaneo

Keele

Titiunik

(2022). Covariate adjustment in regression discontinuity designs. In Zubizarreta

J. R.

Stuart

E. A.

Small

D. S.

Rosenbaum

P. R.

(Eds.), Handbook of matching and weighting adjustments for causal inference (pp. 153–168). Boca Raton, FL: CRC Press.

19.

Cattaneo

Titiunik

(2022). Regression discontinuity designs. Annual Review of Economics, 14, 821–851. https://doi.org/10.1146/annurev-economics-051520-021409

20.

Cattaneo

Vazquez-Bare

(2016). The choice of neighborhood in regression discontinuity designs. Observational Studies, 2(2), 134–146. https://doi.org/10.1353/obs.2017.0002

21.

Center on School Turnaround. (2017). Four domains for rapid school improvement: A systems framework (Center on School Turnaround Four Domains Series). San Francisco, CA: WestEd.

22.

Cook

T. D.

Shadish

W. R.

Wong

V. C.

(2008). Three conditions under which experiments and observational studies produce comparable causal estimates: New findings from within-study comparisons. Journal of Policy Analysis and Management, 27(4), 724–750. https://doi.org/10.1002/pam.20375

23.

Dong

Kelcey

Spybrook

(2023). Experimental design and power for moderation in multisite cluster randomized trials. Journal of Experimental Education, 1–16. https://doi.org/10.1080/00220973.2023.2226934

24.

Gelman

Imbens

(2019). Why high-order polynomials should not be used in regression discontinuity designs. Journal of Business & Economic Statistics, 37(3), 447–456. https://doi.org/10.1080/07350015.2017.1366909

25.

Haegerich

T. M.

Massetti

G. M.

(2013). Commentary on subgroup analysis in intervention research: Opportunities for the public health approach to violence prevention. Prevention Science, 14(2), 193–198. https://doi.org/10.1007/s11121-012-0331-0

26.

Hahn

Todd

van de Klaauw

(1999). Evaluating the effect of an antidiscrimination law using a regression-discontinuity design [Working Paper 7131]. Cambridge, MA: National Bureau of Economic Research. https://doi.org/10.3386/w7131

27.

Hansen

(2015). Punishment and deterrence: Evidence from drunk driving. The American Economic Review, 105(4), 1581–1617. https://doi.org/10.1257/aer.20130189

28.

Hitt

D. H.

Robinson

Player

(2018). District readiness to support school turnaround: A guide for state education agencies and districts (2nd ed., Center on School Turnaround Four Domains Series). San Francisco, CA: WestEd.

29.

Hsu

Y.-C.

Shen

(2019). Testing treatment effect heterogeneity in regression discontinuity designs. Journal of Econometrics, 208(2), 468–486. https://doi.org/10.1016/j.jeconom.2018.10.004

30.

Imbens

Kalyanaraman

(2012). Optimal bandwidth choice for the regression discontinuity estimator. The Review of Economic Studies, 79(3), 933–959. https://doi.org/10.1093/restud/rdr043

31.

Imbens

G. W.

Lemieux

(Eds.). (2008). The regression discontinuity design: Theory and applications. Special issue. Journal of Econometrics, 142(2), 611–614. https://doi.org/10.1016/j.jeconom.2007.05.008

32.

Kraft

(2020). Interpreting effect sizes of education interventions. Educational Researcher, 49(4), 241–253. https://doi.org/10.3102/0013189X20912798

33.

Lee

D. S.

Lemieux

(2010). Regression discontinuity designs in economics. Journal of Economic Literature, 48(2), 281–355. https://doi.org/10.1257/jel.48.2.281

34.

Louie

Rhoads

Mark

(2016). Challenges when using the regression discontinuity design in educational evaluations: Lessons from the transition to algebra study. American Journal of Evaluation, 37(3), 381–407. https://doi.org/10.1177/1098214015621787

35.

Morgan

K. L.

Rubin

D. B.

(2012). Rerandomization to improve covariate balance in experiments. The Annals of Statistics, 40(2), 1263–1282. https://doi.org/10.48550/arXiv.1207.5625

36.

Podolsky

Kini

Bishop

Darling-Hammond

(2016). Solving the teacher shortage: How to attract and retain excellent educators. Palo Alto, CA: Learning Policy Institute. https://doi.org/10.54300/262.960

37.

Reguly

(2021). Heterogeneous treatment effects in regression discontinuity designs. ArXiv preprint: arXiv:2106.11640v3. https://doi.org/10.48550/arXiv.2106.11640

38.

Renaissance Learning. (2009). Star Reading technical manual .

39.

Rubin

D. B.

(1974). Estimating causal effects of treatments in randomized and nonrandomized studies. Journal of Educational Psychology, 66(5), 688–701. https://doi.org/10.1037/h0037350

40.

Shadish

W. R.

Cook

T. D.

Campbell

D. T.

(2002). Experimental and non-experimental designs for generalized causal inference. Boston, MA: Houghton, Mifflin and Company.

41.

Shadish

W. R.

Galindo

Wong

V. C.

Steiner

P. M.

Cook

T. D.

(2011). A randomized experiment comparing random to cutoff-based assignment. Psychological Methods, 16(2), 179–191. https://doi.org/10.1037/a0023345

42.

Silverman

B. W.

(1986). Density estimation for statistics and data analysis (Monographs on Statistics and Applied Probability Series). New York, NY: Springer US.

43.

Spybrook

Kelcey

Dong

(2016). Power for detecting treatment by moderator effects in two and three-level cluster randomized trials. Journal of Educational and Behavioral Statistics, 41(6), 605–627. https://doi.org/10.3102/1076998616655442

44.

Supplee

L. H.

Kelly

B. C.

MacKinnon

D. M.

Barofsky

M. Y.

(2013). Introduction to the special issue: Subgroup analysis in prevention and intervention research. Prevention Science, 14(2), 107–110. https://doi.org/10.1007/s11121-012-0335-9

45.

U.S. Department of Education & the National Science Foundation. (2013). Common guidelines for education research and development (pp. 13–126). Alexandria, VA: NSF.

46.

Wasserman

(2023). Gender differences in politician persistence. The Review of Economics and Statistics, 105(2), 275–291. https://doi.org/10.1162/rest_a_01099

47.

What Works Clearinghouse (WWC). (2022). What Works Clearinghouse procedures and standards handbook, Version 5.0. U.S. Department of Education, Institute of Education Sciences, National Center for Education Evaluation and Regional Assistance.

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.00 MB

0.24 MB