{'fraction_non_alphanumeric': 0.03592734986079809, 'fraction_numerical': 0.01144950767117014, 'mean_word_length': 4.387222438644332, 'pattern_counts': {'":': 1, '<': 8, '<?xml version=': 0, '>': 14, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 21, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Benefits and costs of matching prior to a Difference in Differences analysis when parallel trends does not hold

The Difference in Difference (DiD) estimator is a popular estimator built on the"parallel trends"assumption. To increase the plausibility of this assumption, a natural idea is to match treated and control units prior to a DiD analysis. In this paper, we characterize the bias of matching prior to a DiD analysis under a linear structural model. Our framework allows for both observed and unobserved confounders that have time varying effects. Given this framework, we find that matching on baseline covariates reduces the bias associated with these covariates, when compared to the original DiD estimator. We further find that additionally matching on the pre-treatment outcomes has both cost and benefit. First, it mitigates the bias associated with unobserved confounders, since matching on pre-treatment outcomes partially balances these unobserved confounders. This reduction is proportional to the reliability of the outcome, a measure of how coupled the outcomes are with these latent covariates. On the other hand, we find that matching on the pre-treatment outcome undermines the second"difference"in a DiD estimate by forcing the treated and control group's pre-treatment outcomes to be equal. This injects bias into the final estimate, analogous to the case when parallel trends holds. We extend our bias results to multivariate confounders with multiple pre-treatment periods and find similar results. Finally, we provide heuristic guidelines to practitioners on whether to match prior to their DiD analysis, along with a method for roughly estimating the reduction in bias. We illustrate our guidelines by reanalyzing a recent empirical study that used matching prior to a DiD analysis to explore the impact of principal turnover on student achievement. We find that the authors' decision to match on the pre-treatment outcomes was crucial in making the estimated treatment effect more credible.


Introduction
The difference in difference (DiD) estimator is a popular method in causal inference for measuring treatment effects in panel data (Abadie, 2005;Ryan, Burgess and Dimick, 2014;Wharam et al., 2007). The DiD estimator uses the change of a control group to estimate how much change a targeted treated group would have experienced over time, absent treatment. In such an approach, researchers assume that the control group's change over time is equivalent to the treatment group's change over time absent any treatment effect. This assumption is often referred to as the "parallel trend assumption" (Abadie, 2005) and is used as the main identifying assumption justifying the DiD estimator. One core concern with this assumption is whether the chosen comparison group is suitable.
Given this concern, one might wish to select a control group where parallel trends seems more plausible. An intuitive tool for this would be to use matching methods. Matching methods are a class of methods in causal inference that allows direct comparison of the treated and control groups by selecting units from a pool of possible control units that resemble the treatment group in terms of baseline, observable characteristics (Stuart, 2010;Rubin, 2006). The hope is, if our observed characteristics are informative enough, that the control units would be, other than the receipt of treatment, just like the treatment units in other aspects. For the DiD case, in particular, we would hope our matched controls would change just like the treated units would have, absent treatment. Given this intuition, many previous empirical works have indeed used matching prior to their DiD analysis (Heckman et al., 1998;A. Smith and E. Todd, 2005;Imai, Kim and Wang, 2021); in these studies, the authors believe that the key parallel trend assumption is more likely to hold conditional on similar observed covariates and lagged outcomes. This assumption is often referred as "conditional parallel trends": the evolution of treated and control units would be the same, absent treatment, for units with the same observable characteristics. Although conditional parallel trends may be more likely to hold in certain scenarios, (Daw and Hatfield, 2018) show that when unconditional parallel trends hold, matching on pre-treatment outcomes actually can bias a perfectly unbiased estimator. We are thus faced with a dilemma: matching intuitively gives us a more plausible control group, but the act of matching injects bias into our estimator. We explore this tradeoff further in this paper.
Engaging with this dilemma is not new. Many works analyze matching prior to a DiD analysis and show that matching may actually hurt or help depending on different scenarios (Daw and Hatfield, 2018;Ding and Li, 2019;Lindner and McConnell, 2019;Zeldow and Hatfield, 2021;Chabé-Ferret, 2015;Chabé-Ferret, 2017). These existing works can be roughly grouped into three categories. The first are those that characterize the bias when either parallel trends or conditional parallel trends perfectly holds (Ding and Li, 2019;Daw and Hatfield, 2018). The second are those that characterize the bias in more general settings via simulations (Lindner and McConnell, 2019;Zeldow and Hatfield, 2021). The third give sufficient conditions for when matching combined with DiD gives perfectly unbiased estimates, as opposed to comparing which of the biased estimates may lead to the minimum bias (Chabé-Ferret, 2017;Chabé-Ferret, 2015). 1 We seek to provide exact mathematical characterizations of the bias under a general framework that allows for both imperfect conditional and unconditional parallel trends driven by unobserved and observed confounders. In other words, all estimators, regardless of matching, are biased in our setting as neither parallel trends nor conditional parallel trends necessarily holds, and we explore how to identify which strategy produce the least (as opposed to zero) bias.
In particular, we characterize the bias of different matching DiD estimators in a linear structural equation model (detailed in Section 3) as a function of the underlying relationship between covariates, pre-treatment outcomes, breakage in parallel trends, and amount of confounding. In Section 4, we present the main results of the bias for the simple unmatched DiD estimator, the DiD estimator after matching on observed covariates, and the DiD estimator after matching additionally on the pre-treatment outcome. We find that matching on observed covariates generally leads to a reduction in bias. This improvement is guaranteed unless unlikely scenarios occur such as pre-existing biases cancelling each other out. This result is also consistent with the "pre-treatment criterion" that suggests practitioners always control for observed covariates when estimating a causal effect (Rosenbaum, 2002;Rosenbaum and Rubin, 1983;Shpitser, VanderWeele and Robins, 2010), although our additional results do undermine this principle.
In contrast to matching on our baseline covariates, matching on the pre-treatment outcome exhibits a trade off. On one hand, matching on the pre-treatment outcome does indirectly partially match on the latent confounders, thus mitigating their bias contribution. The amount in reduction is proportional to a key quantity known as reliability, a measure of how coupled the outcomes are with these latent covariates (Prasad, Vaidya and Vemula, 2016;Trochim, 2006). On the other hand, matching on the pre-treatment outcome undermines the second "difference" in the DiD estimator by forcing the treated and control group's pre-treatment outcomes to be equal. Therefore, if parallel trends originally held, matching on the pre-treatment outcome breaks an initially unbiased estimator. On balance, as we show in Lemma 4.2, if the reliability is higher than the breakage of parallel trends, it is better to match on the pre-treatment outcome.
If we view pre-treatment outcomes as just another baseline covariate, the above results connect to matching on covariates measured with error Webb-Vargas et al. (2017); Lenis, Ebnesajjad and Stuart (2017); Rudolph and Stuart (2018). In particular, when we match on such variables, we only imperfectly match on the underlying, latent, variable directly connected to our post-treatment outcomes, leaving residual imbalance. The DiD estimator can then exacerbate this residual imbalance by, as described in Daw and Hatfield (2018), via a regression to the mean effect. These results show the "pretreatment criterion" does not always apply for difference-in-difference contexts.
For clarity of exposition we initially focus on the simple pre-post (or two time point) DiD. After we present these results, however, we extend to a more general setting with multivariate covariates and multiple pre-treatment outcomes in Section 5. We find that matching on multiple pre-treatment outcomes helps further recover the latent confounder, leading to a greater reliability and reduction in bias. In addition to characterizing the bias of the DiD after matching, we also provide heuristic guidelines in Section 6 to practitioners on whether matching prior to their DiD analysis is sensible, along with a method for roughly estimating the reduction in bias. Finally in Section 7, we apply our results to a recent empirical application (introduced in Section 2) that explores whether principal turnover in schools have any impact on student achievement (Bartanen, Grissom and Rogers, 2019). In this application, Bartanen, Grissom and Rogers combined matching with a DiD analysis to show that principal turnover has a negative effect on student achievement the following year. We show that matching, particularly on the pre-treatment outcomes, helped reduce the bias of the estimate, providing further support for the authors' justification to match, and lending further credibility for their reported treatment effect.

Empirical Application -Impact of Principal Turnover
Principal turnover rates in the United States are high, with some estimates at 18% nationally, which has led to much policy attention and research to understand the effect of principal turnover on student achievement (Goldring and Taie, 2018;Barnes and Benjamin, 2007;Grissom, Kalogrides and Loeb, 2015). Prior studies suggest principals play a crucial role in building a learning climate, supporting teacher improvement, and instructional leadership Boyd et al., 2011;Grissom, Loeb and Master, 2013;Coelli and Green, 2012); we might believe that turnover would undermine such benefits. In this vein, a recent study by (Bartanen, Grissom and Rogers, 2019) indeed found short term negative effects of principal turnover on public schools in Missouri and Tennessee. The Missouri data, initially obtained by the Department of Elementary and Secondary Education, contains 2,400 public schools in 565 districts, and span the years of 1991 to 2016. These data record whether there was a turnover in any given year, and also contain outcomes such as school-average math scores on statewide exams in Grades 3 to 8. The original authors also drew on school demographics from the Common Core of Data (CCD) that contained information about student enrollment size, proportion of Black and Hispanic students, and the proportion of students receiving free lunch (Bartanen, Grissom and Rogers, 2019). The main years of analysis were between 2001-2015 to allow for multiple lagged outcomes. As our running example, we use these Missouri data, generously provided by the initial authors.
In our analysis, the main treatment variable is a binary principal turnover variable, which takes value of 1 if a principal in a school in year t is not the principal in the school in year t + 1 and 0 otherwise. We focus on estimating the immediate impact of principal turnover on student's math scores. To answer this, the initial authors estimated a causal effect using a DiD analysis after matching. Bartanen, Grissom and Rogers matched prior to DiD because "even conditional on school fixed effects and the other controls in the model ... parallel trends do not hold between treatment and comparison schools." The authors argue that since principals could be leaving due to factors that are also negatively impacting the schools, the schools that have principals departing most likely already have a declining student achievement trend in the years leading to the principal turnover.
To address this challenge, the authors construct a matched comparison group, matching on school demographics and pre-treatment outcomes. Then they obtain an estimated treatment effect via a DiD analysis through regression. We replicated their analysis, with results on Figure 1. Details of our replication of this matching step and the model we used for the DiD analysis are in Appendix E. The difference in Year 1 (first post-treatment year) between the treated and control schools indicate a substantial negative impact of principal turnover. Furthermore, matching provided a control group (in green) with a similar downward trend prior to the principal turnover as we see for the treatment group (in red). Consequently, the authors argue that matching on the pre-treatment outcomes was essential in forming a valid control and treatment group comparison.
In this paper, we work to formalize this choice. On one hand, Daw and Hatfield (2018) shows matching can break parallel trends, but on the other, parallel trends is not initially holding. We investigate this bias-bias trade-off in Section 7 with our practical guidelines that we present in Section 6. As part of this we also calculate rough estimates of the reduction in bias obtained from matching due to the covariates, and due to the covariates along with pre-treatment outcomes. We find that although the authors decision on matching on observed covariates was potentially inconsequential, matching on the pre-treatment outcomes was critical to obtain a more credible treatment effect.

Linear Structural Equation Model
We present our main results assuming a linear structural equation model following similar studies related to DiD and matching (see (Zeldow and Hatfield, 2021;Daw and Hatfield, 2018) for examples). Although our model assumes a specific parametric setting, we believe it nevertheless serves as a useful starting point for understanding the main tradeoffs at play. Let Y i,t be the measured outcome of individual i at time t. We first focus on the classic DiD setting with one pre-treatment period at t = 0 and one post-treatment period at t = 1. In this pre-post setting, Y i,1 is our post-treatment outcome and and Y i,0 is the pre-treatment outcome. No one has received treatment at t = 0, and some of our units have received treatment at t = 1. We extend our results to multiple pre-treatment periods in Section 5.
Denote Z i as a binary treatment indicator for individual i, with Z i = 1 indicating membership in the treatment group. We have an observed covariate X i for each individual i and an unobserved latent covariate θ i . We first focus on the case when both X i and θ i are one-dimensional and extend our results to a multivariate setting in Section 5. We assume a model of both the latent and observed covariates with the following conditional multivariate Gaussian distribution, Both θ and X are imbalanced across treatment groups when µ θ,1 = µ θ,0 and µ x,1 = µ x,0 . The above model implies a selection mechanism of where N z (θ, x) is the probability density of (θ, X) in Equation (1) given Z i = z. I.e., our model is equivalent to first generating (θ i , X i ) pairs from some distribution and then assigning treatment to the units with the above probabilities (see Remark 3.2 for alternate assignment mechanisms found in the literature).
We use the potential outcomes framework, and assume the Stable Unit Treatment Value Assumption, i.e., no interference between units and no multiple forms of treatment (see Splawa-Neyman, Dabrowska and Speed (1990) for their historical introduction and Rosenbaum (2009) for an overview). In particular, we generate the potential outcomes Y i,t (z), which denotes the outcome we would see at time t had unit i received treatment z, and the observed outcome Y i,t , with the following linear model, we have assumed homoscedastic errors across time for simplicity), i,t ⊥ ⊥ i,t for all t, t , and β 0,t , β θ,t , β x,t are fixed constants that denote the intercept at time period t and the slopes of θ and X at time t, respectively. We assume the vector of random variables ( are independently and identically distributed according to the data generating process described in Equations 1 and 3 across individuals. While both θ and X have time varying effects on the outcome as represented by the β * ,t , θ and X themselves are time invariant, i.e., they do not change over time (although the pre-treatment outcome, viewed as a covariate, does in the case of multiple pre-treatment periods). See Remark 3.1 for further details and connection to time varying confounders.
Our assumption of a constant treatment effect of τ is without loss of generality; because a DiD analysis is fundamentally about imputing a missing control counterfactual, and because we fully observe all Y i (1) in the treated group, the treatment model is arbitrary. We have assumed no anticipation of treatment effect, with Y i,t (0) = Y i,t (1) for t = 0; this is ensured by the indicator 1(t = 1) in the second line.
We critically allow the slopes of θ, X (β θ,t , β x,t ) to differ across time, which allows the confounders θ and X to have a time varying effect on the response, thus breaking the parallel trends assumption. To capture the breakage in parallel trends succinctly, define the imbalance of θ and the time varying change in the effect of θ as δ θ and ∆ θ , respectively, where, We define δ x , ∆ x similarly for the imbalance and time varying effect of X, respectively. Under this notation, the growth of a unit i from pre to post is a function of the observed and latent covariates, with an expected growth of ∆ θ θ i + ∆ X X i ; if the treatment and control groups have different covariate distributions, they will have systematically different growth rates.
Although we simplify our model by not including interactions between θ and X, we analyze the interaction case in Appendix C; we find that it does not offer more insight than that already provided in the no-interaction case. We end this subsection with two remarks connecting our setup to existing literature.  (3) is nearly equivalent to the one considered in (Zeldow and Hatfield, 2021) except our covariates (θ i , X i ) are time invariant while the aforementioned authors consider time varying (θ it , X it ) covariates. We do not consider the time varying covariate case (other than the pre-treatment outcome) for two reasons: 1) As also noted by Zeldow and Hatfield (2021), the time invariant covariate case is interesting by itself because it allows imperfect parallel trends and leads to surprising and illustrative results when some of these covariates are unobserved as shown in Section 4. 2) The time varying covariate case requires more structural assumptions on how the covariates evolve over time. Zeldow and Hatfield (2021) detail that time varying covariates may evolve either a) completely unrelated to the treatment (like an AR(1) process), b) related to the treatment, or c) with a combination of the treatment and other factors. All three cases would require structurally modeling the time varying confounder, where settings b) and c) poses thorny issues due to post-treatment matching. For these reasons, we choose to present our paper with the more simplified, but still widely used, time invariant covariates that have time varying effects in Equation (3) as a starting point (see discussion in Section 8 for further details in time-varying settings).
Remark 3.2 (Assignment Mechanism). The assignment mechanism we consider in Equation (2) presumes that the units select themselves into treatment or control as a function of (X i , θ i ), their baseline characteristics (both observed and latent), as opposed to the actual pre-period outcome Y i,0 . Other evaluations use alternative selection mechanisms; for example, (Chabé-Ferret, 2017;Chabé-Ferret, 2015) instead consider assignment mechanisms that are purely a function of the preperiod outcome with no covariates. In other words, this alternate mechanism can assign units to treatment that by chance observed an extreme outcome in the pre-period (this is the mechanism that would generate, for example, Ashenfelter's dip (Ashenfelter, 1978)). We leave the connection and extension of our results for this setting to future work.

Parallel Trends
The key assumption behind DiD is the parallel trend assumption. The parallel trends assumption states that the change between the pre-treatment and post-treatment control potential outcomes for the treatment is the same as that of the control. More formally, (4) Under our linear model in Equation (3), Equation (4) is equivalent to, Since the pre-treatment period's potential outcome is unaffected by the treatment status, the potential outcome Y i,0 (0) is observed for both the treatment and control group. The post-treatment potential outcome, Y i,1 (0), is also observed for the control group. However, Y i,1 (0) for the treatment group is unobserved, making it the main quantity to identify in a DiD analysis. Equation (4) allows us to impute this quantity, giving an estimate of the average treatment effect for the treated (ATT): Any departure from Equation (4), i.e., if leads to a biased DiD estimator. It is well known that parallel trends holds despite the presence of observed or unobserved imbalance between the covariates as long as the expected change across time in both the treatment and control is still same. Consequently, Zeldow and Hatfield (2021) claim that covariates that are typically confounders in longitudinal studies are not confounders in a DiD study. Instead, they say covariates "that differ by treatment group and are associated with outcome trends are confounders" in a DiD study (pg. 2). In our model, not only are (θ, X) typically confounders in a longitudinal study but also (time invariant) confounders in a DiD study since (θ, X) break parallel trends due to their time varying effects (β θ,t , β x,t ) on the response. Therefore, we refer to (θ, X) as confounders throughout this paper. For concreteness, consider the principal turnover application. Say we observe the proportion of students on Free and Reduced Price Lunch for each school (this is an indicator of the proportion of the student body that live in lower income families). This may be imbalanced since schools that have more principals departing may be schools with higher or lower proportion of FRPL students, thus δ x = 0. Now say there was increased state-wide government sponsored funding for schools above some threshold of FRPL in the post-treatment years. This could change the relationship between the proportion of FRPL students and overall student achievement in the pre-treatment and post-treatment years, i.e., ∆ x = 0. An example of an unobserved confounder θ could be how well the school is able to adapt to statewide changes in the curriculum. Now, if the state testing were changing, this could change the relationship of our outcome and θ across time. Furthermore, we would expect θ to be correlated, ρ = 0, with many observed covariates X (e.g., some underfunded schools may be less able to adapt to changing curriculum).

DiD and Matching
We now formally introduce the relevant DiD and matching estimators. For the rest of this section, we work with the expected values of each estimator (as opposed to the finite-sample estimators) since we care about the bias. We first introduce the classic unmatched DiD estimatorτ DiD (Abadie, 2005). The expected value ofτ DiD is, where the expectation is over the joint distribution of (X i , θ i , i,0 , i,1 ). We refer toτ DiD as the naïve DiD estimator as it does not use matching at all. The naïve DiD estimator simply takes the difference of the differences between each control and treated group. It is unbiased if the parallel trend assumption in Equation (4) is satisfied. For a matched DiD estimator, we would first find a control unit for each treatment unit that shares the same value of X (or both X, Y 0 , if matching on both). We do this for each treatment unit and fit the resulting DiD estimator to this matched data. For the purpose of the paper, we assume arbitrarily close to perfect matching. In other words, if we are matching on X i , we assume that for every treated individual i such that Z i = 1, we obtain a control individual j such that X i = X j and T j = 0. Although perfect matching rarely occurs in practice, it allows us to isolate bias from matching as an approach vs. implementation. It can also be viewed from an asymptotic argument, if the matching variable has finite support; see, e.g., (Abadie and Imbens, 2006). Additionally, many empirical works that assume conditional parallel trends also implicitly assume they have achieved perfect matching since conditional parallel trend requires parallel trends to hold given both the treated and control group have the exact same values of the matched variable(s) (Imai, Kim and Wang, 2021). We leave fully characterizing the bias under imperfect matches to future work.
Before introducing the matched DiD estimators, we define some additional notation to formally denote the matching step. Suppose we are matching a treatment unit with value X i = x. We denote the expected observed outcome of the matched control unit with an expectation of E(Y i,1 | Z i = 0, X i = x). Now, we will have a control unit for each treated unit, so the average of these control units would be the average of the above over the distribution of X i in the treated group, giving an overall expected average outcome of We have three estimators, depending on how we match: matching on X i , matching on pretreatment outcomes, and matching on both. We found matching only on the pre-treatment outcome and not on an available covariate X predictive of outcome was both practically unsound and also unhelpful for providing insight to the benefits of matching. Consequently, we relegate the analysis for this estimator to Appendix B, although we do analyze the case of matching on a pre-treatment outcome when there is no X.
The expected outcomes of the two remaining matching estimators (assuming perfect matches) of matching on the observed covariates X and matching on both X and the pre-treatment outcome Y 0 are then: Importantly, the expression for E τ X,Y 0 DiD reveals one consequence of matching: we no longer have a second difference correction term because we have matched it away. We will later see that if the matching did not achieve true balance on the latent aspects, this lack of correction term leads to bias.
As a final baseline of comparison, we also consider the simple difference in means estimator , which takes the difference in means for the treated and control using only the post-treatment outcomes. We will refer to this estimator as the simple DiM estimator. We now present the main results of when matching may help or hurt prior to a DiD analysis.

The Biases of the Estimators
We first present the bias results in full generality and then consider special cases to gain intuition on what is driving the bias.
Theorem 4.1 (Bias of DiD and Matching DiD estimators). If (Z i , X i , θ i , Y i,t ) are independently and identically drawn from the data generating process as shown in Equations 1-3, then the bias of our estimators are the following, The proof is provided in Appendix A.1. Before moving to simpler cases, we make a few remarks. First, the bias of the naïve DiD estimator is the degree to which parallel trend is broken, i.e., ∆ PT = ∆ θ δ θ + ∆ x δ x . Second, the simple DiM estimator is biased proportionally to the confounding in θ and X (δ θ , δ x , respectively). Lastly, the r θ|x term is related to the reliability of a measure, which is often used in the statistics literature to denote the consistency of a measure (Prasad, Vaidya and Vemula, 2016;Trochim, 2006). The reliability of a measure Y is the proportion of the total variance of Y explained by a true signal rather than noise. More formally, Definition 4.1 (Reliability). The reliability r θ of a random variable Y with respect to a random variable θ is, In particular, we use r θ|x , the (conditional) reliability of the pre-treatment outcome with respect to θ after controlling for X within the control group, i.e., Y i,0 | (X i , Z i = 0). In our linear framework, r θ|x can be interpreted as the population R-squared statistic if we were able to regress the pretreatment outcome Y i,0 on θ i after accounting for X i within the control group. To formally justify this, elementary probability calculations show that r θ|x = Cor 2 (Y i,0 , θ i | X i , Z i = 0) for our linear model. Therefore, r θ|x can also be interpreted as the square of the correlation between θ and the pre-treatment response after accounting for X within the control group. Throughout this paper, we often refer to r θ|x as simply the reliability term. This reliability cannot generally be directly estimated from the data and may differ greatly depending on the specific application and context. We provide a heuristic way to estimate the reliability under some further assumptions and multiple pre-treatment time points in Section 6.
Because it is not immediately obvious from the expressions in Theorem 4.1 how matching on each variables is helping/hurting the DiD estimator, we will build intuition for these results by going through simple sub-cases first. We first build intuition when there exist no covariates X in Section 4.1. We then bring back covariate X but allow zero correlation with θ (ρ = 0) in Section 4.2 before returning to the general case with correlation in Section 4.3.

No Covariate Case
We start with the simplest case where we do not have any observed covariate X to consider the bias of matching only on Y 0 prior to a DiD analysis. If there is no X (i.e., all β x,t = 0), then the bias expressions simplify.
Corollary 4.1.1 (Bias with no covariates). Assume the same conditions in Theorem 4.1, but with no effect from X, i.e., β x,0 = β x,1 = ρ = 0. The bias expressions in Theorem 4.1 then reduce to the following: The naïve DiD, as explained in Section 3.2, has two interpretable terms contributing to the bias. The first term, δ θ , is the confounding effect of θ, i.e., how imbalanced θ is. The second term, ∆ θ , represents a potentially time varying relationship with the outcome that can break parallel trends. The matched estimator, if the reliability is greater than zero, uniformly dominates (has less bias than) the simple DiM estimator, but the comparison to DiD is less direct.
The key intuition of why matching on pre-treatment outcome before DiD could give more credible causal estimates than simple DiD is that matching should reduce imbalance between the treatment and control groups by reducing the differences driven by δ θ . Unfortunately, θ is not directly observed, but we hope that matching on the pre-treatment outcome will indirectly match on θ. The reliability term, r θ|x , is a natural measure of how good the pre-treatment outcome proxies θ.
To compare matching with DiD to DiD without matching, we rewrite the bias of E τ Y 0 DiD as This expression confirms that matching helps reduce the bias proportionally to how high the reliability r θ|x is (note the bias offsetting term of ∆ θ δ θ r θ|x ). It also reveals an additional bias term due to matching, β θ,0 δ θ (1 − r θ|x ), that grows larger as reliability declines. Indeed, if Y 0 is completely dominated by noise, i.e., r θ|x = 0, then the matching estimator achieves the maximum possible bias of a simple difference in means with no adjustment or DiD approach. On the other hand, if the pre-treatment outcome perfectly proxies θ, i.e., r θ|x = 1, then our bias is exactly zero as expected. Finally, if parallel trends perfectly held (∆ θ = 0), then matching on the pre-treatment outcome would result in a non-zero bias of β θ,0 δ θ (1 − r θ|x ); this is the regression to the mean phenomenon discussed in (Daw and Hatfield, 2018).
On one hand, matching on the pre-treament outcome reduces the bias compared to the naïve DiD by recovering the latent confounder θ proportional to the reliability r θ|x . On the other hand, matching also incurs additional bias by effectively the second "difference" in the DiD estimator in terms of θ. We next provide, for this simple no-covariate case, necessary and sufficient conditions to determine when it is better to match on the pre-treatment outcomes as compared to not matching at all.
Lemma 4.2 (Sufficient and necessary condition to match on pre-treatment outcome). Under the same scenario as Corollary 4.1.1, the absolute bias of matching on the pre-treatment outcomes is smaller or equal to the absolute bias of the naïve DiD, i.e., The proof comes from an algebraic simplification of the main result. Intuitively, we should expect the trade off to depends on the reliability vs. how well parallel trends is satisfied. The ratio of the pre-treatment outcome and post-treatment outcome slopes, s, is a scale invariant measure of how non-parallel the outcome-covariate relationship is before and after treatment. For example, s = 1 implies parallel trends hold and the condition reduces to r θ|x ≥ 1, showing that matching on the pre-treatment outcome is never recommended when parallel trends hold.
Parallel trends can break in two different ways. The first is when the post-treatment slope is bigger than the pre-treatment slope β θ,1 > β θ,0 , and they have the same sign, i.e., 0 < s < 1. In this case, the absolute value in Equation (7) is irrelevant and the condition reduces to r θ|x > s, directly illustrating the trade off between reliability and the breakage in parallel trends. For example, if the breakage in parallel trends is very minimal, i.e., s = 0.99 (close to 1) then the reliability has to be at least 99% for the matching estimator to have lower bias than the naïve DiD. In the second case when s > 1, there is a similar interpretation. For example, if s = 1.01, we need our reliability to be at least 99% to have a lower bias. Lastly, if the pre and post slopes of θ have differing signs, i.e., s < 0, or when s > 2, then Equation (7) always holds because the right hand side is less than zero and reliability is always greater than zero. We summarize the above relationship in Figure 2, which shows the two-dimensional regions of (s, r θ|x ) such that Lemma 4.2 holds. As expected, when parallel trends is almost satisfied (near one), one would generally not want to match (note red) unless reliably is very high. On the other hand, when parallel trends is far from perfect (far from one), we benefit from matching more often, even with modest reliability.

Uncorrelated Case
We next consider the more general situation when the observed covariates do affect the outcome, i.e., β x,0 = 0 and β x,1 = 0, but the covariate is uncorrelated with the latent confounder (ρ = 0). This model is fairly general as we could define the latent confounder as only those unexplained differences beyond what can be explained by X.
Corollary 4.2.1 (Bias with zero correlation). Assume the same conditions in Theorem 4.1 but the observed covariate is uncorrelated with θ, i.e., ρ = 0, then the bias expressions in Theorem 4.1 reduce to the following: Unlike matching on Y 0 , matching on X directly gets rid of the bias terms related to X (∆ x δ x ). Although this is generally an improvement from the naïve DiD, the naïve DiD bias could still be lower than the bias from matching on X if preexisting bias from X and θ were cancelling each other. For example, if the time varying effects of θ and X are the same (∆ x = ∆ θ ) but the imbalance of θ and X were in opposite direction (δ x = −δ θ ), then the bias contribution of X and θ would perfectly cancel without matching, and matching on X would result in non-zero bias. Lastly, the bias of matching on both X and pre-treatment outcome reduces to the same bias expression when matching only on pre-treatment outcome without covariates in Section 4.1. This is not a coincidence since matching on X removes any effect X had on the response and the situation simplifies to the no covariate case.

General Correlated Case
Given the results explained in Section 4.1-4.2, we return to understanding the results presented in Theorem 4.1, where there exists a correlation between θ and X (ρ = 0). Theorem 4.1 has two additional complications as compared to Corollary 4.2.1. The first is regarding the reliability r θ|x , i.e., the proportion of variance of the pre-treatment outcome explained by θ after accounting for X. Because accounting for X via conditioning gives information about θ due to the non-zero correlation, it consequently reduces the total conditional variance of the pre-treatment outcome, V ar(Y i,0 | X i , Z i = 0). The consequence of this "extra" information is captured by the (1−ρ 2 ) term in the expression for r θ|x . This means that as X is increasingly correlated with θ, the conditional reliability decreases, suggesting an increase in the cost of matching additionally on lagged outcomes.
The second complication is a countervailing force where the more correlated X and θ are, the smaller (in general) the post-matched imbalance in θ because as we match on X we also can reduce δ θ towards zero. This is not guaranteed, however, as we can illustrate with the maximal ρ = 1 case. As a reminder, ρ = 1 does not imply that θ and X are perfectly correlated but instead implies that θ and X are perfectly correlated within each of the treatment and control groups. Therefore, even if ρ = 1, matching on the same value of X does not necessarily lead to the same value of θ for both the treatment and control group because the group averages of θ differ. Consequently, the usefulness of the extra information we gain from matching on a correlated X depends not only on the strength of the correlation ρ but also on how much the imbalance of X (δ x ) is similar to the imbalance of θ (δ θ ) (after scale adjustment with the σ θ σx term putting X on the scale of θ). In other words, if the imbalances δ θ and δ X are the same (once δ X is rescaled by σ θ σx ) then the bias from matching on a perfectly correlated X is zero, as stated by Theorem 4.1.
Overall, the correlation term can decrease the conditional reliability, and thus increase the bias as discussed initially, but it also can allow matching on X to indirectly match on θ. We generally find that the more correlated X and θ, the less benefit of matching additionally on pre-treatment outcomes. In the extreme, Theorem 4.1 shows that the bias will always increase if ρδ x is in the opposite direction of δ θ , i.e., sign(ρδ x ) = sign(δ θ ), where sign(.) is the sign function. In such a case, the correlation actually harms our matching DiD estimator and fails to recover θ by pushing it further (in the opposite direction) from the desired value. Additionally, even if the confounding effects were in the same direction as the correlation, if δ x δ θ then X can "over-correct" the confounding effect of θ, leading to possibly greater bias. We formalize the necessary conditions that removes these "edge cases" in Section 6.1, where we give a rule-of-thumb guide for determining whether to additionally match on lagged outcome.  Figure 3: Bias of the naïve DiD estimator (blue), matching on X DiD estimator (green), and matching on both X and the pre-treatment outcome (red) according to the results in Theorem 4.1. We fix the following parameters β θ,1 = β x,1 = 1.5, β θ,0 = β x,0 = δ θ = δ x = σ θ = σ x = 1, r θ|x = 0.5, and ρ = 0. We then vary one of the fixed parameters for each of the three panels in Figure 3. For example, in the left panel we vary δ x and keep all the remaining parameters fixed to their mentioned values.

Takeaways
There is unfortunately no simple condition such as Lemma 4.2 that guarantees that it is always better to match on X or both X and the pre-treatment outcome under a more complicated general setting. We instead visually illustrate the findings above by plotting the bias of the naïve DiD estimator and the two matching estimators across three panels in Figure 3. To produce the plots, we initially fix β θ,1 = β x,1 = 1.5, β θ,0 = β x,0 = δ θ = δ x = σ θ = σ x = 1, r θ|x = 0.5, and ρ = 0. We then vary one of the fixed parameters for each of the three panels in Figure 3. For example, in the left panel we vary δ x . Although the results of the plot are sensitive to the initial parameter values (which were chosen based off our application in Section 7), we show the plots to visually summarize the main findings and show the complex ways the biases can differ. The left panel in Figure 3 shows that matching on X generally reduces the bias relative to doing nothing since the bias of the naïve DiD estimator (in blue) increases as the confounding strength of X increases, while both matching estimators are unaffected by the confounding strength of X (in this uncorrelated case). This confirms the guideline on why it is generally better to match on observed covariates. However, as mentioned in Section 4.1, matching on X is not uniformly better than the naïve DiD: when δ x is negative, the bias of X cancels with the positive bias of θ to allow the naïve DiD (in blue) to have a lower bias than that of the matching estimator that matches on X only (in green).
The second panel shows that as the reliability grows, it is better to additionally match on pre-treatment outcome as illustrated by the red bias curve decreasing to zero.
The third panel shows how matching on X improves as the correlation of X and θ increases. As expected, the bias of both matching estimators is exactly zero when the correlation is one since in this scenario the imbalance of θ, X are the same, allowing matching on X to perfectly recover θ. When the correlation is negative, however, the bias increases for both matching estimators since matching on X actually pushes θ in the opposite direction, as explained in Section 4.3.
Overall, given the above, we find that, regardless of whether parallel trends hold or not, matching on observable covariates X is generally always advised for the following reasons: 1. Matching on X directly gets rid of all biases and confounding that results directly from X (∆ x δ x term).
2. Matching on X may further help recover the unobserved confounder θ when they are correlated. This shows up in the bias by offsetting the effect of δ θ by ρ σ θ σx δ x . That being said, matching on X does not always reduce the bias since pre-existing bias from θ and X may have been cancelling each other out. Furthermore, a non-zero correlation between θ and X does not always help reduce the bias also since the correlation ρ and the confounding effects, δ x , δ θ , may act in opposite directions or matching on X may "over-correct" the bias.
Matching on the pre-treatment outcome Y 0 has a more complex trade-off that depends on the following two quantities: 1. The reliability of the pre-treatment outcome for latent variation beyond what is explained by X. If the conditional reliability is high it will generally favor matching on Y 0 since it helps us recover the latent confounder θ.
2. The amount of breakage in parallel trends of the unobserved confounder θ. Because matching on the pre-treatment outcome erodes the second "difference" in the DiD estimator, matching on the pre-treatment outcome is only favored if parallel trends does not hold.
This means that the decision to additionally match on pre-treatment Y 0 depends on the size of these quantities. Given these tensions, we provide some guidance on how to decide when to match on what variables in Section 6.

Multiple Time Periods and Multivariate Confounders
So far we have focused on the setting when there is only one pre-treatment time point and univariate confounders X, θ. It is common, however, to have outcome measures for more than one pretreatment time point and multiple covariates. For example, in our motivating empirical examples introduced in Section 2, Bartanen, Grissom and Rogers matched on seven observed covariates and six pre-treatment outcomes. In cases such as these, we may want to match on all available pre-treatment outcomes and multiple covariates prior to the DiD analysis. We generalize our linear structural equation model introduced in Section 3.1 to account for T pre-treatment outcomes and multivariate confounders with the following model, We index 0 ≤ t ≤ T − 1 for the T periods prior to post-treatment period t = T . Additionally, β x,t is a p-dimensional column of slopes at time t for X i = (X i,1 , . . . , X i,p ), where X i,j represents individual i's value of covariate j for up to p observed covariates. We define β θ,t , θ i similarly for a qdimensional number of latent confounders. We keep the data generating model for our confounders similar to Equation (1) but allow arbitrary covariances across the multivariate confounders, where Σ XX ∈ R p×p and Σ θθ ∈ R q×q are the covariance matrix for X i and θ i within each treatment and control group respectively, and Σ θX ∈ R q×p is the covariance matrix between θ i and X i . As before each individual's (θ i , X i ) is drawn independently and identically according to Equation (9), where µ θ,z is a q-dimensional column vector for the average of θ in treatment group z and µ x,z is defined similarly. Lastly, we denote Y T i = (Y i,0 , . . . , Y i,T −1 ) to represent the vector of T pretreatment period outcomes.
It is common to use a linear regression framework when there are multiple time periods (Callaway and Sant'Anna, 2021;Wooldridge, 2021;Card and Sullivan, 1987;Tyler, Murnane and Willett, 2000). In particular, one might estimate the treatment effectτ using a two-way fixed effects linear regression, where W t is an indicator that is 1 if in the post-treatment period (t = T ) and zero otherwise, α t , γ i are fixed effects for time and unit, respectively, and i,t is the error. The estimateβ ZW for the interaction term would then be take as the DiD estimate of the ATT. It is well known that in a balanced panel data with multiple pre-treatment outcomes,β ZW is equivalent to the estimate one would obtain using a classical two-period DiD using the average of all the T pre-treatment outcomes, Y i,0:(T −1) = T −1 t=0 Y i,t /T , for the pre-treatment measure (Callaway and Sant'Anna, 2021;Wooldridge, 2021). We will use this balanced property to derive the costs and benefits of matching before DiD using Model 10.
Given our model and the balanced panel setting, the expected values of the generalized DiD estimators we would obtain, if using the above two-way fixed effects model for estimation, is then: where Y T i = (Y i,0 , Y i,1 , . . . , Y i,T −1 ) represents the collection of all pre-treatment outcomes. Under our model, giving the biases of our estimators as stated in the following theorem: are independently and identically drawn from the data generating process in Equations 8-9, then the bias of our estimators are the following: with vector analogs to baseline imbalance and breakage of parallel trends of δ θ := µ θ,1 − µ θ,0 and, with δ x , ∆ x defined similarly, and Σ Y T Y T ∈ R T ×T the covariance matrix of all pre-treatment periods Y 0 , . . . , Y T −1 , with Σ θY T ∈ R q×T , Σ XY T ∈ R p×T defined similarly.
The proof is provided in Appendix A.2. The interpretation for the naïve DiD and the DiD estimator that matches only on X remains similar to that already discussed in Section 4. For example, the naïve DiD is biased again by the amount the multivariate confounders θ, X are breaking the parallel trends. Matching on X similarly reduces bias contribution from X and offsets the confounding effect of the latent confounder by the correlation Σ θX and the shift δ x . However, it is unclear from the expression above how matching on the T additional pre-treatment outcomes reduces/harms the bias. To explore this further, we apply Theorem 5.1 to a simpler case that allows us to simplify the bias expression.

Stable Pretreatment Case
To understand how matching on multiple pre-treatment outcome may impact bias, we first assume a univariate unobserved confounder θ, i.e., we let q = 1. We further assume β x,t and β θ,t are unchanging across the pre-treatment period outcomes, i.e., parallel trends hold in the pre-treatment periods, as given in the following assumption: Assumption 1 (Unconditional parallel trends for pre-treatment outcomes). For all pre-treatment periods, 0 ≤ t, t ≤ T − 1, we assume there is no time varying effects of both confounders θ and X.
This assumption effectively states that the researcher has T independent and identical measurements of pre-treatment outcome Y available to them, i.e., for all the pre-treatment outcomes the parallel trends hold perfectly (unconditional on X or θ).
The univariate θ and pre-treatment stability are working assumptions. First, unless one has strong subject-matter knowledge of regarding missing confounders, one loses little generality by representing all latent confounding with a single confounder θ with arbitrary breakage in parallel trends (∆ θ ) and confounding (δ θ ). Second, the effects of θ and X may generally vary in the pretreatment periods, violating Assumption 1. However, even if that were the case, it does not change how additional pre-treatment period gives more measurements to proxy θ. Therefore, we emphasize that we assume this to give simplified expressions that illustrate this point.
The matching estimator that only matches on X will remain unchanged compared to Theorem 5.1, so we state the main bias results only for the matching estimator that involves matching additionally on the pre-treatment period outcomes.
Theorem 5.2 (Bias: Multiple Time Periods). If (Z i , X i , θ i , Y i,t ) are independently and identically drawn from the data generating process in Equations 8-9 and Assumption 1 holds with q = 1, then the biases of our matching estimators are then The proof is provided in Appendix A.2. The result is similar to Corollary 4.2.1 (and exactly equivalent when T = 1 and θ and X are uncorrelated), except our new "reliability" term, r T θ|x , increases as a function of T . 2 The tildes in the expressions in Theorem 5.2 represent a "residualization" of the outcome and θ by X driven by the correlation of X and θ. Namely, δ θ and σ 2 θ are adjusted by the extra information gained from matching on correlated X (see Section 4.3 for detailed explanations). For example, the more tightly coupled X and θ, the more adjusting on X might give little benefit for additionally matching on pre-treatment outcomes, as represented by a smallerσ 2 θ and generally lower r T θ|x . Note that if p = 1, then Σ θX = ρσ θ σ x ,δ θ = δ θ − ρ σ θ σx δ x , and σ 2 θ = σ 2 θ (1 − ρ 2 ), which recovers Theorem 4.1. As expected, the more pre-treatment period outcomes we match on, the more matching recovers our latent confounder θ. To illustrate how much our reliability increases, suppose r T =1 θ|x = 0.5. Then if T = 2, i.e., we have one more additional pre-treatment period, then r T =2 θ|x ≈ 0.67. If T = 3, then r T =3 θ|x = 0.75, which is a 50% increase of reliability with only two additional time periods. Additionally, the bias for matching on both X and all pre-treatment period outcomes tends to zero as the number of pre-treatment periods goes to infinity. In other words, having greater numbers of pre-treatment observations allows for more and more precise assessment of the latent θ, eventually leading to zero bias. This result is similar to synthetic control results in (Abadie, Diamond and Hainmueller, 2010;Bouttell et al., 2018) and has also been replicated for simulations presented in (Illenberger, Small and Shaw, 2020).

Determining When to Match
Using our model as a working approximation, we now use the theoretical results in Section 4 to provide heuristic guidance on what to match on, along with a means of roughly estimating the reduction (or increase) in bias due to matching. We additionally provide a publicly available script to run our proposed guidelines. 3 In particular, we now give a set of guidelines on whether to match on the observed covariates X and additionally on the pre-treatment outcome(s). We assume a univariate θ, and allow for arbitrary breakage in parallel trends and degree of imbalance.

Guidance for Matching on Covariates
Guideline 1 (Matching on X guideline). Always match on X.
2 Strictly speaking r T θ|x is not the reliability of our outcome as defined in Definition 4.1 However, it can still be roughly interpreted as the ratio of how much total variance from the T pre-period outcome is explained by θ.
One can estimate the reduction in bias (relative to the naïve DiD) from matching on X as: whereˆ δ x is the difference in means of X between the treated and control group andˆ β x,t are the p-dimensional regression coefficients from linear regressions of Y i,t on X i within the control group, one regression for each time point.
This advice is consistent with the current advice on how practitioners should generally account for as many observed covariates as possible (Rosenbaum, 2002;Rosenbaum and Rubin, 1983;Shpitser, VanderWeele and Robins, 2010;Ding and Miratrix, 2015), which is often referred to as the "pre-treatment criterion." ∆ τx is a rough estimate of the degree of bias reduction (relative to the naïve DiD). More formallŷ ∆ τx estimates ∆ τx , where Section 4.2 shows that, when X and θ are uncorrelated, matching directly on the observed covariates X reduces the bias contribution from X, i.e., gets rid of the ∆ x δ x term in the bias. As long as biases are not cancelling, this term would then be the bias reduction. It turns out that, even if Σ θX = 0, naïvely taking the difference of the estimated slopes from a linear regression still accounts for how the correlation affects the bias, because the correlation with θ gets picked up by the estimated slope coefficients themselves being biased. We formalize how this bias estimate works in the following theorem: Theorem 6.1 (Conditions for matching on X reducing bias). Suppose a data matrix of D i = (X i , Y i,0 , . . . , Y i,T , Z i ) for i = 1, 2, . . . , n follows the linear structural equation model in Equation (8)-9 and the following three sign conditions hold: . Then ∆ τx ≥ 0, i.e., matching on X will have less bias than the naïve DiD. Additionally,∆ τx is a consistent estimator for ∆ τx .
The proof is provided in Appendix D.1. In general, while matching on X removes bias from the observed covariates, doing so does not guarantee overall bias reduction due to the "edge cases" mentioned in Section 4.4. The three sign conditions in Theorem 6.1 remove these edge cases. The first sign condition says the pre-existing biases of θ and X are not in opposite directions. The second sign condition ensures the imbalance of δ θ is reduced, not increased, by the extra information we gain about θ by matching on a correlated X. In other words, the second sign condition does not allow the confounding effects of θ and X to go in the opposite direction of the correlation. The third sign condition does not allow the additional reduction in bias gained from matching on a correlated X to over-correct the bias. We believe that these three sign conditions are plausible in reasonable applications, but we leave it to future empirical work to verify this. Importantly, these conditions are not necessary, in that there are many cases where they do not hold, but matching on X is still beneficial.

Guidance for Matching on Pre-treatment Outcome(s)
While the guidance for matching on X is relatively straightforward, the same is not true for matching on the available pre-treatment outcome(s). On one hand, matching on pre-treatment outcomes reduces the effect of the unobserved confounder δ θ by a factor proportional to the reliability. On the other hand, matching on pre-treatment outcomes erodes the second "difference" in the DiD analysis, which can add bias to the overall estimate. Since these quantities are consequences of how latent parameters relate to the outcome, a general guideline will necessarily rely on a set of strong untestable assumptions and heuristics. Here, we use our simplified model from above to provide such a guideline.
Guideline 2 (Matching on X, Y T guideline). Match on all available pre-treatment outcome(s) if the following inequality holds.
where the reliability and regression coefficients can be estimated as described below.
We can estimate the reduction in bias from matching additionally on T pre-treatment outcomes (relative to matching only on X) aŝ where∆ τx,y is an estimate for ∆ τx,y , formally defined as i.e., the additional bias reduction from matching on both observed covariates pre-treatment outcomes relative to just matching on observed covariates. We estimate the components of the above as described below.
Guideline 2 extends Lemma 4.2, which provided the changeover point between matching and not matching in the no-covariate case. In particular, Equation (11) comes from directly utilizing the results in Theorem 5.2 to extend Lemma 4.2 to our multiple time-point, covariate case: Lemma 6.2 (Sufficient and necessary condition to match on T pre-treatment outcomes in addition to X). Suppose a data matrix of D i = (X i , Y i,0 , . . . , Y i,T , Z i ) for i = 1, 2, . . . , n follows the same setting as that listed in Theorem 5.2 and we define the additional reduction in bias from matching on both observed covariates and pre-treatment outcome(s) compared to just matching on observed covariates, ∆ τx,y , with Equation 12. Then we have ∆ τx,y > 0 if where β θ,t is defined in Equation (8) and r T θ|x is defined in Theorem 5.2. In other words, when Equation (13) holds, matching on X, Y T has less bias than only matching on X.
The proof is in Appendix D.2. The intuition of this lemma is that by first taking out the predictive element of X both directly and through its correlation with θ, we can reduce our covariate case to the no-covariate case, and then follow the ideas in Section 4.1, with δ θ replaced withδ θ , the imbalance of the "residualized" latent confounder.
To estimate the quantities needed conduct the guideline check, we first residualize our outcomes for each time point (we require multiple pre-treatment observations, i.e., T > 1): where theˆ β x,t are the regression coefficients for regressing the outcomes at timepoint t onto the covariates, and then estimate the residual variance of the residualized model and finally calculate the empirical estimates of the quantities given in Lemma 6.2's Equation (13): The additional quantities in the formula for bias reduction are estimated aŝ ∆ θ =β θ,T −β θ,T −1 and The above estimation formulae may seem as if we are forgetting to account forσ θ . We show in Appendix D.2 that the final estimate of ∆ τxy does not require estimatingσ θ becauseσ θ cancels out.
Our guideline and estimation technically relies on two simplifying assumptions. The first assumption requires that the parallel trends assumption holds in the pre-treatment period, i.e., β θ,t = β θ,t for all t, t < T , with T ≥ 2, and the same is true for the slopes of X. This pretreatment parallel trend assumption allows us to obtain consistent estimates of β θ,t across time t, and thus a consistent estimate of the reliability r T θ|x . The second assumption is that our data are independently and identical drawn from our linear structural model. This is important because our estimation procedure utilizes the linear structural model to provide consistent estimation.
Although these assumptions are unlikely to hold exactly in practice, they do provide a heuristic for deciding on the benefits of matching. We believe our guidelines to be reasonable even under some deviation from these assumptions, but we leave the robustness of our guidelines, along with the case of T = 1, to future work. In particular, if our assumptions are perfectly valid, then one could make a bias correction by just subtracting out the bias induced by θ since we have consistent estimates for β θ,t and δ θ . However, we choose to present it this way to emphasize that our guidelines are a heuristic that should extend to other settings even when the linearity and stable pre-treatment period assumptions fail. Additionally, we acknowledge that other approaches, e.g., ones that incorporate additional information or data, may lead to superior performance.
We demonstrate estimating these quantities in our applied example in Section 7. We also formally summarize the assumptions that theoretically justify Guideline 2 and show consistency in the above estimation procedure in the following theorem: Theorem 6.3 (Conditions for beneficial matching on X, Y T ). Suppose a data matrix of D i = (X i , Y i,0 , . . . , Y i,T , Z i ) for i = 1, 2, . . . , n follows the same setting as that listed in Theorem 5.2. Then the rule in Guideline 2 assures ∆ τx,y > 0 by Lemma 6.2, i.e., matching on X, Y T has less bias than only matching on X.
If we further suppose T > 1, then∆ τx,y is a consistent estimator for ∆ τx,y ,r T θ|x is a consistent estimator for r T θ|x , andβ θ,T −1 β θ,T is a consistent estimator for 7 Application -The Impacts of Principal Turnover As introduced in Section 2, Bartanen, Grissom and Rogers (2019) are interested in determining whether the impact of principal turnover has any causal impact on student achievement a year after the principal has changed. We follow the guidelines in Section 6 using the same data the authors used with a few adjustments. We first include all seven standardized observed covariates as X (further details in Section 2). 4 Although Figure 1 shows the original matching DiD estimates by Bartanen, Grissom and Rogers, we do not actually need to perform any matching to use our guidelines detailed above. Nevertheless, we also empirically match to show how the treatment effects may vary for illustrative purposes (further detailed below). This application has staggered adoption, in that we have a set of years, and in a given year some schools are treated (lost a principal) and other schools are not. We, therefore, apply our guidelines to each year in turn, and then average the results, weighting by the number of treated units in a given year. Our provided script implements this aggregation automatically.
To quantify the reduction in bias from matching on X only, for example, we estimate the bias contribution by following Guideline 1. More specifically, we estimate ∆ x by obtaining the post-treatment and pre-treatment slopes of X through seven separate linear regression of Y i,t for t = 0, 1, . . . , 7 on X i in the control group, respectively. We also estimate δ x by taking the simple difference in means of X from the treatment and control group across the seven covariates. We repeat this for years 2005 through 2016, and then weight across all years, weighting by number of treatment units in each year.
Final estimates are shown in the top left entry of Table 1. The estimates correspond (left to right) to school enrollment size, proportion of students with free lunch, proportion of Hispanic students, proportion of Black students, proportion of new-to-school teachers, principal experience, and average number of principal transition in the past five years.
The authors reported an estimated treatment effect of aboutτ = −0.035 (Bartanen, Grissom and Rogers, 2019). The first column of Table 1 shows that matching on observed covariates X reduces bias by about∆ τx = 0.012. The estimated breakage in parallel trends for X,∆ x (estimated by taking the difference between the linear regression slope of X on the post-treatment outcome and the average slope of X across all the pre-treatment outcomes), is generally low but not zero.  Table 1: Results from applying guidelines in Section 6 on principal turnover empirical application. Column 1 and 2 contains results for guidelines on matching on X and Y T respectively. The first row contains key parameters that summarize how the matching may help or not. The second row is a binary answer on whether the practitioner should match for the respective variables prior to their DiD analysis. The last row contains an estimate of the reduction in bias from the resulting matching (if second row is yes).
Consequently, the naïve DiD without matching on X would have likely given a similar result to matching on X.
To illustrate this point, we compared naïve DiD and simple matching on these data. The weighted-average point estimate of the naïve DiD, obtained from doing a simple DiD on every year and then taking the weighted average of the resulting impact estimates, is −0.023. By contrast, the point estimate obtained from DiD after doing one-to-one matching on X 5 is −0.034. Although we should not expect the difference to be exactly∆ τx because of imperfect matching and probable deviations from the linear model, in this case the difference is very close to our estimated bias difference.
To determine whether to additionally match on the six available pre-treatment outcomes (T = 6), we estimate the conditional reliability and breakage in parallel trends of our latent covariate (see Guideline 2). We being by obtaining "residualized" responsesỸ i,t for all time periods through the residuals of seven separate linear regression of Y i,t on X i within the control group at all available pre-treatment periods t = 0, 1, 2, . . . , 6. Using the residualized outcomes as the new outcome, we then obtain an estimate of the variance of the noise in the outcome,σ 2 E , by leveraging the information from all six pre-treatment outcomes and the fact that the error of the variance is assumed homoscedastic throughout the pre-treatment periods. We next use the equations in Section 6.2 to obtain estimates of the pre-and post-treatment θ coefficients by exploiting the fact that, e.g., the variance of the residualized pre-treatment outcome is a simple function of the pre-treatment slope and σ 2 E . This givesβ θ,T −1 = 0.48 andβ θ,T = 0.72. We also estimate the imbalance of θ, gettingδ θ = −0.22, since the difference in means of the treated and control group's residualized pre-treatment outcome is proportional to δ θ under our linear model. Finally, we estimate the reliabilityr T θ|x = 0.86 by using a simple plug-in estimator, using the estimate ofβ θ,T −1 andσ 2 E . The second column of Table 1 summarizes the key estimated parameters. Our provided script again automates all of these estimation steps. Table 1 shows that matching on the six available pre-treatment outcome leads to a significant estimated reduction in bias of around 0.033, which is almost the same in magnitude of the estimated treatment effect ofτ = −0.035. This is largely because the estimated reliability (given six pretreatment outcomes) is quite high and the breakage in parallel trends,β θ,T −1 β θ,T = 0.66, is also high (the post-treatment slope of θ is approximately 1.5 times larger than the pre-treatment slope).
For illustrative purposes, we calculated the weighted-mean of the DiD point estimate for each year after matching on both X and the six additional pre-treatment outcomes using the same matching procedure as used for matching on X alone. We found the matching DiD point estimate to be −0.04, different than the naïve DiD point estimate of −0.023. These results are consistent with (Bartanen, Grissom and Rogers, 2019) as the authors emphasize the importance of matching on pre-treatment outcomes as schools who had principals departing were observed to have different pre-treatment trends than other schools. Our guidelines show that the authors decision to match on pre-treatment outcomes was likely essential for obtaining a credible treatment effect.

Concluding Remarks
We explore the bias of matching prior to a DiD analysis when parallel trends does not initially hold by mathematically characterizing the bias under a linear structural equation. We verify that matching on observed covariates likely leads to a reduction in bias. Further, we find that matching on pre-treatment outcomes exhibits a bias-bias trade off between recovering the latent confounder and undermining the parallel trend.
We use our results to create guidelines for determining whether matching is recommended, and further provide an estimate of the resulting reduction in bias. We apply this strategy to a recent application of a DiD analysis involving matching to evaluate the impact of principal turnover on student achievement (Bartanen, Grissom and Rogers, 2019). We find evidence that the authors' decision to match on all available pre-treatment outcomes did in fact reduce bias and create a more credible causal estimate.
Our work, however, is not comprehensive. First, our results are specific to the linear structural equation model. Although this is a useful starting point that allows for flexible exploration, it still contains strong parametric limitations. We did, however, consider interaction effects in Appendix C; these findings did not substantively differ from those found in Theorem 4.1. We also only consider an assignment mechanism where units self-select into treatment or control based on observed and unobserved confounders as opposed to selection based on the actual pre-period outcome as detailed in Remark 3.1. Additionally, our setup does not explicitly allow for time varying covariates, instead focusing on time varying relationships between covariate and outcome (Remark 3.1).
The tests for whether to match, discussed in Section 6, are based on our linear model and additional stability assumptions. We acknowledge that there are potentially more interesting directions and robust ways to determine whether matching is suitable for a specific application; this could be fruitful area for future work. Assessing how rigorous our proposed guidelines are in the face of model misspecification is also an interesting question.
Lastly, we only consider the theoretical bias under a perfect matching scheme. In most applications, perfect matches do not exist and in some cases having even an approximately good match is difficult. An interesting future direction would be to quantify the bias under imperfect matching, perhaps even accounting for systematic imperfect matching based on covariates.

A Proofs
A.1 Proof of Theorem 4.1 The proofs for the bias of the naïve DiD and the difference in means estimator are straight forward, thus we omit them. We also leave out all time intercept term β 0,t when evaluating the relevant expectations because they cancel out. Before proving the bias of the various matching estimators, we will use a well known and useful fact of the multivariate normal distribution.

A.1.2 Matching on covariates and past outcome
Proof.
The key is that we have the following multivariate Gaussian, Next, we apply Lemma A.1 to obtain where we define the following to avoid long algebraic expressions 6 : The rest follows from simple algebraic manipulation.

A.2 Proof of Section 5 (Multiple covariates and time points)
We first prove Theorem 5.1. We omit the proof of the bias of the naïve DiD estimator and matching on only X as the results come directly from elementary calculation or applying Lemma A.1 to find . We now prove the bias of matching on both observed covariate and the T pre-treatment outcomes.
Proof. To compute this bias, we must find . Similar to the above proofs, the key is to realize that: We again apply Lemma A.1 to obtain the desired result.
We now prove Theorem 5.2. We first state two useful Lemmas about matrix inversions: Lemma A.2 (Matrix Inversion). Suppose we have a matrix Q (T +1)×(T +1) such that: where x, y ∈ R such that (x − y) = 0 and (x + T y) = 0 then we have that, . . .
The proofs of these Lemmas are omitted because these are well established properties of matrix inverses.
Proof. We first assume that Σ θX = 0, i.e., θ and X are uncorrelated. We prove the case when Σ θX = 0 later. We remind readers that we are under Assumption 1 that states that all slopes of θ and X are the same in the pre-treatment periods. Consequently, without loss of generality we write all the pre-treatment slopes of θ and X as β θ,0 and β x,0 , respectively. We start the proof by writing the bias as, Focusing on the inner expectation E(θ i | Z i = 0, X i = x, Y T i = y t ), we notice the following multivariate Gaussian distribution using the fact that θ and X are uncorrelated, and A is the same matrix structure as the matrix Q in Lemma A.
, we need to invert matrix D. This can be done by first applying Lemma A.3 to matrix D giving us, The key is that Thus, M is again the same matrix structure as matrix Q in Lemma A.2 with x = (β θ,0 ) 2 σ 2 θ +σ 2 E and y = (β θ,0 ) 2 σ 2 θ . Consequently, we can apply Lemma A.2 to get M −1 , which is again the same form as the matrix Q in Lemma A.2 with x = Putting everything together we have that: where a is some matrix that will be irrelevant to further calculations. Returning back to calculating We also have from simple matrix multiplication that Returning, to the more general case for arbitrarily correlated θ and X, i.e., Σ θX = 0. We first defineθ that captures the information left from θ after residualizing how much the correlated X can explain θ in the following way,θ This new latent confounder is purposefully constructed to be uncorrelated with X. In other words, we have the following facts aboutθ, We can then rewrite our new latent confounder by, Then, using Equation (14), we can write the bias (in this more general setting) as, It remains to find E(θ i | Z i = 0, X i = x, Y T i = y t ). However, the key is that Cov(θ i , X i | Z i = 0) = Σ θX − Σ θX = 0 (by construction). Thus, we are back in the uncorrelated case. More formally, we have the following multivariate normal distribution, and C is defined similarly as above and A is again the same matrix structure as the matrix Q in Lemma A.2 with x = (β θ,0 ) 2σ2 θ + β x,0 Σ XX β x,0 + σ 2 E and y = x − σ 2 E . Finally, following the same steps as the above proof for the uncorrelated case, we realize that Putting this back into Equation (15) gives the desired result since β θ,0 = β θ,T −1 under Assumption 1.

B Matching only on Pre-treatment Outcome
In this section we explore the matching estimator that matches only on the pre-treatment outcome Y 0 . We show that it does not offer more insights than that already offered by Theorem 4.1. We first introduce the estimator that matches only on pre-treatment outcome and state our results under the same univariate confounder setting in Section 4.
We now provide the bias result for this matching estimator that matches only on the pretreatment outcome under the same linear model described in Section 3.
Theorem B.1 (Bias of DiD and Matching only on Pre-treatment Outcome). If (Z i , X i , θ i , Y i,t ) are independently and identically drawn from the data generating process as shown in Equations 1-3, then the bias of matching only on the pre-treatment outcome is the following, In the case where we do not match on the observed covariate X, both X and θ are treated as "unobserved" since we do not directly match on X. Therefore, the pre-treatment outcome aims to proxy both of these confounders. Consequently, matching on Y 0 partially recovers both θ and X proportional to their respective contribution to the reliability of Y 0 as represented by the r θ , r x , respectively. r θ , r x terms approximately captures the proportion of variance of Y 0 each variable captures respectively. If θ accounts for more variance of Y 0 , then matching on Y 0 reduces the confounding effect of θ through a larger term of r θ . Finally, even if both X and θ do not have time varying effects on the outcome (∆ θ = ∆ x = 0), then matching on the pre-treatment outcome does not result in zero bias since it again erodes the second difference in the DiD estimator as shown by the bias expression above.
The main difference in this matching estimator compared to the matching estimator that matches on both X and Y 0 is the third and fourth terms in the above bias expression. Unlike the case when we are able to match on X and perfectly get rid of the effect of X, matching only on Y 0 allows the confounders' imbalance (δ x , δ θ ) to continue making an impact on the match, leading to the third and fourth terms. To illustrate this point, we consider the following hypothetical scenario. Suppose that θ accounted for all the total variance of Y 0 , i.e., r θ = 1, r x = σ E = 0. Therefore, we should hope to perfectly recover θ when matching for Y 0 . Further suppose that µ x,1 = 5 and µ x,0 = 3 while β x,1 = β θ,1 = 1. Consequently, Y i,0 | (Z i = 0) = θ i + 3 for the control while Y i,0 | (Z i = 1) = θ i + 5 for the treatment (ignoring the intercept). If we observe a treatment unit with Y i,0 | (Z i = 1) = 6, thus θ i = 1, and find a match j from the control unit Y j,0 | (Z j = 0) = 6, then θ j = 3 = θ i . We see that the difference between θ i and θ j is proportional (in this case exactly equal since all slopes are one) to δ x . Even if θ accounted for the total variance of the pre-treatment outcome, matching on the pre-treatment outcome does not recover θ due to the imbalance of X. Therefore, the third term in the bias expression captures this shift of θ due to the imbalance of X. The fourth term is the symmetrical shift for X due to the imbalance of θ.
In summary, the main insight from matching on just the pre-treatment outcome is similar to that presented in Section 4. On one hand matching on the pre-treatment helps correct for the bias contribution of both θ and X through the reliability terms of the respective variable. On the other hand, if parallel trends originally held matching on the pre-treatment outcome biases the DiD estimator. From the point of the practitioners, it is insensible to not match on the observed covariate X while matching on the pre-treatment outcome as the pre-treatment outcome now has to proxy both X and θ. Putting all this together, we chose to not present this matching estimator as it is both undesirable and does not lead to further insight than that already provided in Section 4.

C Extending to Include Interactions
In this section, we analyze the case when we add an interaction between θ and X. The idea of an interacted model is allow for the relationship between the latent confounder and outcome to change for different levels of the observed covariate X. The interacted model is Y i,t (0) = β 0,t + β θ,t θ i + β x,t X i + β θx,t X i θ i + i,t , with Y i,t (1), Y i,t defined as they were in Equation (3).
In this section, we show the bias results of the naïve DiD and the two main matching estimators when ρ = 0 for this model. We focus on the uncorrelated case to clearly show how the interaction term does not lead to any new insights than those provided in Section 4. Adding a correlation between θ and X only contributes by further offsetting the bias of θ with the correlation, as described in Section 4.3.

C.1 Interpretation of Bias Results
We first focus on the bias of matching on X only. Section 4 shows matching on the observed covariates helped mitigate the bias contribution of X but does not fully reduce the bias to zero due to the latent confounder. The interaction case is similar, except matching on X partially removes the bias contribution of the interaction by allowing the bias contribution of the interaction term to come through the imbalance of θ (as opposed to both the imbalance of θ and X via δ θx ). Contrary to what one may have expected, matching on X when X and θ are interacted does not improve the plausibility of the conditional parallel trend (conditional on X) compared to the no-interaction case. Similar to the results in Section 4, the conditional parallel trends (conditional on X) only holds when the unconditional parallel trends originally held due to the latent confounder θ. When matching additionally on the pre-treatment outcome, we again have a similar story to the no-interaction results presented in Section 4. First,r θ|x acts as the new "reliability" when we have interactions. In this case, matching on the pre-treatment outcome reduces the bias by both main effect and interaction effect of θ through this new reliability term. The only difference is the "leftover" terms (β θ,1 + β θx,1 µ x,1 )L in the bias expression. These terms appear because matching on X does not cleanly get rid of all contribution of X due to the interaction similar to the results presented in Appendix B. Therefore, we see δ x and δ θx continuously making an impact on the bias by acting as a hindrance to how well we recover θ. This "leftover" effect of X is captured by the L term and will be exactly zero if the interaction effect did not exist, i.e., if β θx,0 = 0. In summary, because the main ideas presented in Section 4 still hold in the interaction case, we choose to leave this interaction case in the Appendix.

D Proofs of Guidelines
D.1 Proof of Theorem 6.1 Proof. We prove the first part of Theorem 6.1, i.e., that the three sign conditions guarantee that ∆ τx > 0. Although Section 6.1 states the results assuming θ is univariate, the result holds generally for multivariate θ, thus we prove it in this more general case. In other words, we want to show that given our three sign conditions we have the following, We have two cases. The first case is if ∆ θ δ θ > 0. We can then drop the absolute value for the left side because of the first sign condition. ∆ θ δ θ > 0 can hold in two ways. The first is when ∆ θ > 0, δ θ > 0. In this case, we can drop the absolute value for the right side because of the third sign condition, leading to the following inequality By the second sign condition, the left hand side is positive while the right hand side is negative, thus we are done. ∆ θ δ θ > 0 holds also when ∆ θ < 0, δ θ < 0. In this case, we can similarly drop the absolute value for the right side because of the third sign condition, leading to the same steps as above, thus proving the first case. The second case is if ∆ θ δ θ < 0. Using the first sign condition, we can drop the left hand side's absolute value after multiplying by negative one. Following the same steps as above, ∆ θ δ θ < 0 Consequently, instead of the old latent confounder θ i , we have a new latent confounderθ i = θ i − Σ θX Σ −1 XX X i . We next show the consistency ofσ 2 E . First, with a slight abuse of notation we define the population version ofỸ i,t asỸ