EconBase
← Back to paper

Non-linear Triple Changes Estimator for Targeted Policies

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

44,126 characters · 12 sections · 24 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Non-linear Triple Changes Estimator for Targeted Policies

abstractThe renowned difference-in-differences (DiD) estimator relies on the assumption of `parallel trends,' which does not hold in many practical applications. To address this issue, the econometrics literature has turned to the triple difference estimator. Both DiD and triple difference are limited to assessing average effects exclusively. An alternative avenue is offered by the changes-in-changes (CiC) estimator, which provides an estimate of the entire counterfactual distribution at the cost of relying on (stronger) distributional assumptions. In this work, we extend the triple difference estimator to accommodate the CiC framework, presenting the `triple changes estimator' and its identification assumptions, thereby expanding the scope of the CiC paradigm. Subsequently, we empirically evaluate the proposed framework and apply it to a study examining the impact of Medicaid expansion on children's preventive care.

Introduction

In the domains of econometrics and quantitative social sciences, the difference-in-differences (DiD) estimator and its extensions have emerged as indispensable tools for estimating causal effects in observational studies. With its roots traced back to the work of snow1855mode and further popularized by seminal works such as ashenfelter1984using and card1993minimum, DiD has provided a reliable method for evaluating the effect of policy interventions over time. In particular, DiD stands out in studies where inference based on controlling for confounders or using instrumental variables is deemed unsuitable, and where pre-treatment information is available. Researchers have actively pursued various extensions of the DiD estimator, aimed at enhancing its applicability and robustness athey2006identification, sofer2016negative, callaway2021difference, roth2023parallel.

The DiD estimator hinges on the assumption of parallel trends, which states that the (average) outcome of both the control and treatment groups share exactly the same evolution trend in the absence of the treatment. Mathematically speaking, this assumption translates to

equation[equation omitted — 132 chars of source]

where $Y^0(t)$ denotes the potential outcome associated with the absence of the treatment at time $t$, and $D$ denotes the assigned treatment. In particular, $D=1$ and $D=0$ represent the treatment and control groups, respectively. The validity of this assumption can be challenged, and parallel trends might get violated due to unobserved time-varying factors or dynamic changes in the study context. This violation may result in biased estimates. As a response to potential deviations from parallel trends, researchers have turned to more flexible models, such as the triple difference framework gruber1994incidence. The triple difference estimator can be formulated as the difference between two DiD estimators. Intuitively, the difference of two DiD estimators is unbiased, provided both estimators have the same bias. Indeed, sole purpose of of subtracting the second DiD estimator is to debias the first one. Despite the prevalent use of the triple difference estimator, especially over the last two decades raifman2018association, sakurai2020relationship, han2016effect, chen2020triple, tai2001racial, only recently a formal presentation of the framework and its identification assumptions was provided olden2022triple.

\import{./figures/}{fig0.tex}

The triple difference framework guarantees only the identification of the average treatment effect (on the treated). Another issue that arises is that the triple difference estimator is biased if the outcomes in the pre- and post-treatment or among control and treatment groups are measured on a different scale. The changes-in-changes (CiC) framework athey2006identification, on the other hand, is scale-invariant and yields the identification of the counterfactual outcome probability distribution. To do so, CiC framework requires additional assumptions beyond Eq. (ref). It postulates that there exists a unique monotone mapping $T_d$ which maps the probability measure over $Y^0(t_0)$ to that over $Y^0(t_1)$ within group $D=d$ for $d\in\{0,1\}$ (see Figure (ref).) Moreover, these two mappings are assumed to be identical. That is, $T_1(y) = T_0(y)$ for every $y$, or since the mappings are bijective,

equation[equation omitted — 66 chars of source]

where $\textrm{Id}$ represents the identity map. Eq. (ref) states that there is no drift in the evolution trend across groups. Note that Eq. (ref), was not explicitly stated in athey2006identification as an assumption but can be derived as a consequence of the following four assumptions:

assumption[Model assumption] The potential outcomes $Y^0(t)$ can be modelled using a production function $h$ of a latent variable $U$, which models the individual characteristics: \begin{equation*} \forall t\in\{0,1\}:\quad Y^0(t) = h(U; t). \end{equation*}

In particular, (ref) posits that $h$ does not depend on the group assignment ($D$).

assumption[Strict monotonicity] The function $h(\,\cdot\,; t)$ is strictly increasing in $U$ for every $t\in\{t_0,t_1\}$.
assumption[Time invariance] Within every subgroup, the distribution of the latent variable $U$ does not change over time.
assumption[Latent support overlap] The support of the latent variable $U$ in the treatment group is a subset of its support in the control group.

In the case of a one-dimensional outcome, mappings $T_d$ can be expressed as

equation[equation omitted — 99 chars of source]

where $F_{Y^0(t)\,\vert\, D=d}(\cdot)$ represents the cumulative density function of $Y^0(t)$ in group $D=d$. Under regularity conditions, $T_d$ is the unique monotone map that pushes forward the probability measure over $Y^0(t_0)$ to that over $Y^0(t_1)$ in group $D=d$ villani2009optimal, santambrogio2015optimal. Combining Equations (ref) and (ref), it is straightforward to identify the counterfactual distribution of $Y^0(t_1)$ in the treatment group, $F_{Y^0(t_1)\,\vert\, D=1}(y)$. Specifically, \[

split[split omitted — 164 chars of source]

\] which matches Eq. (9) in the original work of athey2006identification.

The no-drift assumption specified in Eq. (ref) can be challenged in practice, especially in scenarios where the treatment or exposure is directed toward a specific sub-population. This situation arises in studies on the impact of targeted interventions, such as a criminal justice initiative ($D$), on recidivism rates, where this intervention is exclusively administered to individuals with specific criminal histories or risk profiles. Naturally, it is expected that the time evolution of counterfactual recidivism rates will exhibit significant disparities between the control and treated groups. Similar challenges arise when the eligibility criteria is narrow within the context of social programs such as welfare or housing assistance, which target specific demographic groups.

Recognizing the strengths and weaknesses of both the triple difference and CiC frameworks, we propose a novel estimator that combines the best of both worlds. Our proposed `triple changes' estimator aims to overcome the limitations of DiD and the triple difference estimator by leveraging the scale invariance and strong identification results of CiC, while introducing a more flexible mapping assumption that allows us to relax (ref). We briefly present the contributions of this work.

itemize[leftmargin=*] • We formally present and analyze the triple changes estimator as an extension to CiC, and provide the necessary assumptions for its point identification first in the scalar case. We then discuss how to generalize our results to high-dimensional outcomes by harnessing theory of optimal transport. • We provide several partial identification results under relaxed versions of our proposed point identifiability assumptions. Further, we show the validity of analogous results for the classic CiC framework as a special case of our derivations. • We introduce a finite-sample estimator for the average treatment effect on the treated within our framework and analyze its asymptotic behaviour. • We conduct an empirical evaluation of our estimator on both synthetic and real datasets.

This paper is organized as follows. Section (ref) reviews the necessary background and the setup of the study. In Section (ref), the identification of our estimand of interest is studied for a scalar outcome. Section (ref) provides an estimator for the latter and studies its asymptotic properties. In Section (ref), we extend our work to high-dimensional outcomes using theory of optimal transport. Numerical evaluations are presented in Section (ref).

Causal model

We consider a study where we have access to data from two sources, e.g., two states of the united states, or two cities, or any two separate populations. These two sources of data will be denoted by $S=s_0$ and $S=s_1$ throughout. We assume that a treatment (e.g., a health-care policy) is administered in one state, without loss of generality in $S=s_1$, and not in the other. In both states, the individuals are partitioned into two cohorts, namely, $D=d_1$ and $D=d_0$, signifying the individuals that are eligible and not eligible for receiving the treatment, respectively\footnote{In a classic controlled trial, these would correspond to the treatment and control groups, respectively.}. The outcome is measured in two time points, namely $t_0<t_1$, where the eligible individuals in state $s_1$ receive the treatment at an infinitesimal amount of time after $t_0$\footnote{Note that we do not limit our setting to panel data. In particular, the individuals for which the outcome is measured may differ across time points.}. We denote by $Y^{D=d_0}(t)$ and $Y^{D=d_1}(t)$ the potential outcome variables associated with the outcome at time $t$ in the absence, and in the presence of treatment, respectively. To improve readability, we will often use the short-hands $Y^0(t)$ and $Y^1(t)$ for $Y^{D=d_0}(t)$ and $Y^{D=d_1}(t)$, respectively. We denote the observed outcome at time $t$ by $Y(t)$. Throughout, we make the following standard consistency assumption rubin1980randomization.

assumption[Consistency] At each time $t$, the realized outcome $Y(t)$ is determined as \[Y(t)=\sum_d\mathbbm{1}\{D=d\}\cdot Y^{D=d}(t),\] where $\mathbbm{1}\{\cdot\}$ denotes the indicator function.

The estimand of interest is the effect of treatment on the treated, i.e., the group corresponding to $S=s_1,D=d_1$. As such, we target learning the probability measure over the counterfactual outcome $Y^0$ within this subgroup: \[F(Y^0\mid S=s_1,\:D=d_1),\] where $F$ denotes the cumulative density function. When clear from context, we use the shorthand $F_{Y^0\mid s_1,d_1}$ instead. Throughout, we assume that random variables are defined over a compact domain, and that densities are absolutely continuous with respect to the Lebesgue measure.

One-dimensional Estimator

\import{./figures/}{figure1}

We commence our analysis by considering cases where the outcome of interest, $Y$, is one-dimensional, i.e., a scalar.

In contrast to athey2006identification, we posit that the outcome of an individual can be determined by a combination of the latent variable $U$ (with a common support across groups), and the group to which the individual belongs. This adjustment relaxes the model assumption of the CiC framework, as outlined below.

customass{1}[Model assumption] At each state $S=s$ and treatment group $D=d$, the potential outcomes $Y^0(t)$ can be determined through a production function $h_{s,d}(\cdot\:;t)$ of a latent variable $U$, which models the individual characteristics: \begin{equation*} \forall t\in\{0,1\}:\quad Y^0(t) = h_{s,d}(U; t). \end{equation*}

A key distinction between our setup and the classic CiC framework lies in relaxing (ref) to (ref), allowing the production functions $h_{s,d}(\cdot)$ to be group-specific. In applications involving targeted treatment assignments, (ref) emerges as a more sensible assumption. Assumptions (ref) and (ref) are adapted analogously:

customass{2}[Strict monotonicity] The production functions $h_{s,d}(\,\cdot\,; t)$ are strictly increasing in $U$ for every $t\in\{t_0,t_1\}$, and every $s,d$.
customass{3}[Time invariance] Within every subgroup, the distribution of the latent variable $U$ does not change over time. That is, $\forall u, \forall s,d$, \[F_{U\vert S=s, D=d, T=t_1}(u) = F_{U\vert S=s, D=d, T=t_0}(u).\]

Additionally, we prefer to formulate the overlap assumption in relation to the potential outcomes rather than the latent variable $U$. This preference arises from the broader accessibility and interpretability of the support of the outcome, as opposed to that of latent characteristics.

customass{4}[Outcome support overlap] The potential outcomes $Y^0(t_0)$ and $Y^0(t_1)$ are defined over domains $\mathbb{Y}_0$ and $\mathbb{Y}_1$, respectively, which are common across subgroups $S\in\{s_0,s_1\}, D\in\{d_0,d_1\}$. Moreover, $\mathbb{Y}_0\subseteq\mathbb{Y}_1$.

Assumptions (ref) through (ref) establish the existence of four distinct monotone maps that push forward the density of $Y^0(t_0)$ to that of $Y^0(t_1)$ in each subgroup. In particular, let $T_{s,d}$ denote the monotone map that pushes forward the density of $Y^0(t_0)$ in the group corresponding to $S=s, D=d$ to the density of $Y^0(t_1)$ in the same group (see Figure (ref).) Specifically, $T_{s,d}$ can be expressed in terms of the $h(\cdot)$ functions as\footnote{For the purposes of this section, Eq. (ref) can be considered as the definition of maps $T_{s,d}$.}

equation[equation omitted — 92 chars of source]

or equivalently, in terms of the cumulative density functions,

equation[equation omitted — 106 chars of source]

It is noteworthy that under (ref), $F_{Y^0(t)\mid s,t}=F_{Y(t)\mid s,t}$ for every subgroup except $\{s_1,d_1\}$. Consider the map $T^*_s=T_{s,d_1}\circ T^{-1}_{s,d_0}$. $T^*_s$ can be interpreted as the non-linear drift between the maps from $Y^0(t_0)$ to $Y^0(t_1)$ among control and treatment groups in any state $s$ (see Figure (ref).) Identification in CiC is achieved by assuming that there is no drift, i.e., $T^*_s$ is the identity map (see Eq. (ref)). We shall proceed by relaxing this assumption as follows.

assumption[State-independent drifts] The drift between the mappings of potential outcomes $Y^0(t_0)$ to $Y^0(t_1)$ in the control and treatment groups is independent of the state. Formally, \[T^*\coloneqq T_{s_0,d_1}\circ T^{-1}_{s_0,d_0} \equiv T_{s_1,d_1}\circ T^{-1}_{s_1,d_0}.\]

Specifically, rather than assuming no drift, we have relaxed the assumption to equal drift across the two states. (ref) can also be expressed in terms of the production functions, albeit at the cost of interpretability:

multline[multline omitted — 458 chars of source]

We are now ready to state our identification result. The complete set of proofs for our results can be found in Appendix (ref).

restatable{theorem}{thmid} Under assumptions (ref) - (ref) and (ref) - (ref), the cumulative density function of the missing counterfactual $Y^0(t_1)$ in the group $S=s_1,D=d_1$ is identified as: \begin{equation}\begin{split} F_{Y^0(t_1)}&_{\mid s_1,d_1}(y) =\& F_{Y(t_0)\mid s_1,d_1} \circ T^{-1}_{s_0,d_1} \circ T_{s_0,d_0} \circ T^{-1}_{s_1,d_0} ( y) , \end{split} \end{equation} where $T_{s,d}$ is given by Eq. (ref).
remarkTo avoid unnecessarily heavy notation, we did not discuss the observed covariates. However, an identical analysis can be done after adjusting for the observed covariates, $X$. In particular, production functions may depend on the observed covariates, as long as their monotonicity in $U$ is maintained for every $x$ in the domain of $X$. (ref), (ref) and (ref) need to be valid conditioned on $X$ in this scenario, and Eq. (ref) must hold when every term is conditioned on $X$.
remarkWe articulated our identifiability assumptions in accordance with the original work of athey2006identification. An alternative way of presenting the assumptions would be to do it akin to the quantile-quantile equi-confounding bias assumption proposed by ghassami2022combining. In our context, this would translate to directly assuming the existence of monotone maps $T_{s,d}$ based on Eq. (ref) instead of drawing conclusions from (ref)-(ref) to establish it.
remarkWe formulated production functions to model the potential outcomes under no treatment ($Y^0$), whereas we left the other potential outcome, $Y^1$, unrestricted. Due to symmetry, one could model $Y^1$ using production functions and leave $Y^0$ unrestricted. This scenario might arise for instance, in an study where 3 out of 4 cohorts receive treatment. One should exercise greater caution in such cases however, since in practice, the treatment may have effects that significantly alter the composition of the population under study. Under these circumstances, assuming monotone production functions for $Y^1$ could be a more drastic assumption.

Relaxing monotonicity

As mentioned earlier, (ref) (and its counterpart, (ref) in athey2006identification) is an untestable assumption, and might be drastic to impose in certain applications. The monotonicity of functions $h_{s,d}(\cdot;t)$ in $U$ has two implications: (i) these functions are bijective, establishing a well-defined inverse for them; (ii) the property that $\mathbbm{P}(h_{s,d}(U;t)\leq y)=\mathbbm{P}(U\leq h^{-1}_{s,d}(y;t))$, which is repeatedly utilized in proving the point identification result in Theorem (ref) (see Appendix (ref)). While the bijectivity of $h_{s,d}(\cdot;t)$ appears to be essential for our framework to work, we can relax (ii). Let us first rephrase the monotonicity assumption in the equivalent form:

equation[equation omitted — 122 chars of source]

This assumption can be relaxed as follows.

assumption[$\epsilon$-monotonicity] Functions $h^{-1}_{s,d}(\cdot;t)$ are well-defined. Additionally, for any $u$ in the support of $U$, \begin{equation}\mathbbm{P}\big(\langle U- u, h_{s,d}(U;t)- h_{s,d}(u;t)\rangle<0\mid s,d\big)\leq\frac{\epsilon}{2}.\end{equation}

For example, the monthly income of an individual in terms of her age after adjusting for the other covariates can be a $\epsilon$-monotone function. In general, monthly income increases due to promotions and inflation. However, temporary unemployment and retirement can affect this trend. See Figure (ref) for a visualization.

figure[figure omitted — 199 chars of source]

We expect that under this relaxation, point identification cannot be achieved. However, the following partial identification result holds.

restatable{proposition}{prppartial} Under assumptions (ref), (ref), (ref) and (ref) - (ref), for any $y\in\mathbb{Y}_1$, \begin{equation*} \begin{split} F_{Y(t_0)\mid s_1,d_1}&\circ \mathrm{\Phi}^{-1}_{s_0,d_1} \circ \mathrm{\Phi}_{s_0,d_0} \circ \mathrm{\Phi}^{-1}_{s_1,d_0} ( y) -\epsilon \&\leq F_{Y^0(t_1)\mid s_1,d_1}( y) \leq\\ &F_{Y(t_0)\mid s_1,d_1}\circ \mathrm{\overline{\Phi}}^{-1}_{s_0,d_1} \circ \mathrm{\overline{\Phi}}_{s_0,d_0} \circ \mathrm{\overline{\Phi}}^{-1}_{s_1,d_0} ( y) +\epsilon, \end{split} \end{equation*} where \[ \begin{split} \mathrm{\underline{\Phi}}_{s,d}(y)=F^{-1}_{Y(t_1)\vert s,d}\big(F_{Y(t_0)\vert s,d}(y)-\epsilon\big),\\ \mathrm{\underline{\Phi}}^{-1}_{s,d}(y)=F^{-1}_{Y(t_0)\vert s,d}\big(F_{Y(t_1)\vert s,d}(y)-\epsilon\big),\\ \mathrm{\overline{\Phi}}_{s,d}(y)=F^{-1}_{Y(t_1)\vert s,d}\big(F_{Y(t_0)\vert s,d}(y)+\epsilon\big),\\ \mathrm{\overline{\Phi}}^{-1}_{s,d}(y)=F^{-1}_{Y(t_0)\vert s,d}\big(F_{Y(t_1)\vert s,d}(y)+\epsilon\big). \end{split} \]
remarkNote that strict monotonicity is a special case of (ref), which corresponds to $\epsilon=0$. Accordingly, Proposition (ref) reduces to Theorem (ref) with the choice of $\epsilon=0$.
remarkAn analogous result applies to the original framework of CiC when (ref) is relaxed to $\epsilon$-monotonicity. See Appendix (ref) for details.

Relaxing time invariance

In certain applications, (ref) may also be violated. In particular, in repeated cross-sections, it might be challenging to maintain the same distribution among the cases under study. Even if possible, this may both reduce the number of available samples and result in selection bias, making inference more complicated. As such, we consider a relaxation of (ref) and derive a partial identification result in this setting. We assume that the distribution of the latent variable may change over time, but this change is bounded in Kolmogorov (aka KS) distance kolmogorov1933sulla.

restatable[$\delta$-invariance]{assumption}{asdelta} Within every subgroup, the Kolmogorov distance of the distribution of the latent variable across $T=t_0$ and $T=t_1$ is bounded by $\delta$. That is, \[\sup_u\big\vert F_{U\vert S=s,D=d,T=t_1}(u)-F_{U\vert S=s,D=d,T=t_0}(u)\big\vert\leq\delta.\]

Note again that when $\delta=0$, (ref) reduces to (ref).

restatable{proposition}{prpdelta} Under assumptions (ref), (ref), (ref), (ref), and (ref)-(ref), for any $y\in\mathbb{Y}_1$, \begin{equation*} \begin{split} F_{Y(t_0)\mid s_1,d_1}&\circ \mathrm{\Psi}^{-1}_{s_0,d_1} \circ \mathrm{\Psi}_{s_0,d_0} \circ \mathrm{\Psi}^{-1}_{s_1,d_0} ( y) -\delta \&\leq F_{Y^0(t_1)\mid s_1,d_1}( y) \leq\\ &F_{Y(t_0)\mid s_1,d_1}\circ \mathrm{\overline{\Psi}}^{-1}_{s_0,d_1} \circ \mathrm{\overline{\Psi}}_{s_0,d_0} \circ \mathrm{\overline{\Psi}}^{-1}_{s_1,d_0} ( y) +\delta, \end{split} \end{equation*} where \[ \begin{split} \mathrm{\underline{\Psi}}_{s,d}(y)=F^{-1}_{Y(t_1)\vert s,d}\big(F_{Y(t_0)\vert s,d}(y)-\delta\big),\\ \mathrm{\underline{\Psi}}^{-1}_{s,d}(y)=F^{-1}_{Y(t_0)\vert s,d}\big(F_{Y(t_1)\vert s,d}(y)-\delta\big),\\ \mathrm{\overline{\Psi}}_{s,d}(y)=F^{-1}_{Y(t_1)\vert s,d}\big(F_{Y(t_0)\vert s,d}(y)+\delta\big),\\ \mathrm{\overline{\Psi}}^{-1}_{s,d}(y)=F^{-1}_{Y(t_0)\vert s,d}\big(F_{Y(t_1)\vert s,d}(y)+\delta\big). \end{split} \]
remarkEven more generally, we can simultaneously relax (ref) to (ref) and (ref) to (ref). See Proposition (ref) in Appendix (ref) for the partial identification result for this case.

Marginal contrasts and joint counterfactuals

Theorem (ref) guarantees the identification of the probability density of the missing counterfactual, $Y^0(t_1)$ in the group corresponding to $S=s_1,D=d_1$. Having access to this density, we can compute any marginal contrast estimand franks2019flexible. A marginal contrast estimand is an estimand that can be expressed as a functional of the marginal distribution of the counterfactual outcomes. This includes a vast majority of commonly used estimands, such as average treatment effects, conditional average treatment effects, quantile treatment effects, risk ratios, etc. However, in certain applications, more information is desired. Examples of estimands that are not marginal contrasts include the distribution of the treatment effect, individual-level treatment effects, or the quantiles of the treatment effect. As a concrete example, consider $Z \coloneqq Y^1(t_1) - Y^0(t_1)$. The density of $Z$ is not identifiable from merely the marginal densities of $Y^1(t_1)$ and $Y^0(t_0)$. In particular, the treatment may have non-zero effects on a fraction of the population even if the two potential outcomes have identical distributions. Estimands that are not marginal contrasts require stronger identifiability assumptions in general. In this section, we discuss a stronger version of (ref) that can lead to the identification of such estimands.

customassu{3}[Strong time invariance] The latent variable $U$ does not change over time, and there is no loss to follow-up.

(ref) is stronger than (ref), in the sense that it completely rules out the possibility of any changes in the latent variable itself, or the cohort of the study. In contrast, (ref) would allow for loss to follow-up, as long as similar individuals were recruited for the study. (ref) even accommodates a repeated cross-sections study. (ref) also allows for the evolution of the latent variable, as long as its distribution among the study population remains the same. On the other hand, the stronger assumption (ref) allows for the identification of a wider range of causal estimands. In particular, under (ref), the joint density of the counterfactuals $\big(Y^0(t_1), Y^{1}(t_1)\big)$ is identified.

restatable{proposition}{prpjoint} Under assumptions (ref), (ref), (ref), (ref), and (ref) - (ref), the joint density of $Y^0(t_1)$ and $Y^1(t_1)$ in the group corresponding to $S=s_1, D=d_1$ is identified as \begin{equation*} \begin{split} &F_{Y^0(t_1), Y^1(t_1)\mid s_1,d_1}(y^0,y^1) =\& F_{Y(t_0), Y(t_1)\mid s_1,d_1}\big(F^{-1}_{Y(t_0)\mid s_1,d_1} \circ F_{Y^0(t_1)\mid s_1,d_1}(y^0),y^1\big), \end{split}\end{equation*} where $F_{Y^0(t_1)\mid s_1,d_1}(\cdot)$ is given by Equation (ref).

Inference

In this section, we discuss the estimation aspect of our framework with a focus on the average effect of treatment on the treated. More formally, we consider the estimation of \[\tau\coloneqq\mathbb{E}[Y^1(t_1)-Y^0(t_1)\vert S=s_1,D=d_1].\] Note that under the identifiability assumptions of Section (ref), $\tau$ is identified as

multline[multline omitted — 397 chars of source]

We make the following assumption on the data generating mechanism to render estimation feasible.

assumptionConditioned on $S=s,D=d,T=t$, variables $Y(t)$ are continuous random variables defined on a shared bounded domain $[\underline{y},\overline{y}]$, with continuously differentiable density functions $f_{sdt}$, where $f_{sdt}$ is bounded from above and away from $0$, and $\partial f_{sdt}/\partial y$ is bounded. For all $s,d,t$, $p_{sdt}=\mathbb{P}(S\!=s,D\!=d,T\!=t)>0$, and given $s,d,t$, the samples $Y_i(t)$ are independent draws from $f_{sdt}$.

Let $N$ be the total number of observed outcome samples. Akin to athey2006identification, we build a finite-sample estimator for $\tau$ based on empirical estimators of cumulative density functions. In particular, let $\{Y_{sd,i}(t)\}_{i=1}^{N_{sdt}}$ denote the independent samples of the outcome at time $t$ in group $S=s, D=d$. We define

equation[equation omitted — 153 chars of source]
equation[equation omitted — 139 chars of source]

Finally, the estimator for $\tau$ is built as:

multline[multline omitted — 524 chars of source]

The following theorem establishes the consistency and asymptotic normality of $\hat{\tau}$.

restatable{theorem}{thmconsistency} Under (ref), $\hat{\tau}-\tau = \mathcal{O}_p(N^{-\frac{1}{2}})$, and \begin{multline} \sqrt{N}(\hat{\tau}-\tau)\overset{D}{\rightarrow}\mathcal{N}(0, \dfrac{V_0}{p_{s_1d_1t_1}}+ \dfrac{V_1}{p_{s_0d_1t_1}}+ \dfrac{V_2}{p_{s_0d_1t_0}}+\\ \dfrac{V_3}{p_{s_0d_0t_0}}+ \dfrac{V_4}{p_{s_0d_0t_1}}+ \dfrac{V_5}{p_{s_1d_0t_1}}+ \dfrac{V_6}{p_{s_1d_0t_0}}+ \dfrac{V_7}{p_{s_1d_1t_0}} ), \end{multline} where $\{V_i\}_{i=0}^7$ are given by Eq. (ref).

The expression for the variance and the discussion on its estimation are postponed to Appendix (ref).

Optimal Transport Representation and High-dimensional Extension

So far, we focused on the special case of scalar outcome for ease of presentation. In this section, we discuss the generalization of our results to cover multi-dimensional outcomes. To this end, we first review the one-dimensional case from an optimal transport point of view. Let $\eta_{s,d}$ and $\mu_{s,d}$ denote the probability measures over $Y^0(t_0)$ and $Y^0(t_1)$ in the group corresponding to $S=s,D=d$, respectively. Brenier's theorem brenier1991polar implies that there exists a unique mapping $T$ such that $T_\#\eta_{s,d}=\mu_{s,d}$, i.e., $T$ pushes forward $\eta_{s,d}$ to $\mu_{s,d}$, and $T$ is the gradient of a convex function. Moreover, $T$ is the optimal transport map with quadratic cost. More precisely, let $\Gamma$ denote the space of joint distributions over $\big(Y^0(t_0),Y^0(t_1)\big)$ conditioned on $S=s, D=d$, that agree with the marginal densities $\eta_{s,d}$ and $\mu_{s,d}$. The optimization problem

equation[equation omitted — 92 chars of source]

has a unique solution $\gamma^*$, where $\big(Y^0(t_0),Y^0(t_1)\big)\sim \gamma^*$ if and only if $Y^0(t_0)\sim \eta_{s,d}$ and $Y^0(t_1) = T\big(Y^0(t_0)\big)$, $\eta_{s,d}-a.s$.\footnote{Note that we have omitted the implicit conditioning on $S=s,D=d$ in our notation to improve readability.} It is straightforward to verify that a one-dimensional function is the gradient of a convex function if and only if it is monotone. Therefore, the monotonicity of $h_{s,d}(\cdot)$ and the monotonicity of $T_{s,d}(\cdot)$ (as a consequence of the latter) imply $T_{s,d}\equiv T$. In other words, the mapping $T_{s,d}$ is identified as the optimal transport map with quadratic cost that pushes forward $\eta_{s,d}$ to $\mu_{s,d}$.

Indeed, monotonicity can be slightly relaxed. The identifiability of the map $T_{s,d}$ in one dimension is guaranteed under co-monotonicity of $h_{s,d}(\cdot;t_0)$ and $h_{s,d}(\cdot;t_1)$:

equation[equation omitted — 127 chars of source]

In order to achieve identifiability results in higher dimensions based on Brenier's theorem, we need to make sure that $T_{s,d}$ is the gradient of a convex function. It is known that a function is the gradient of a convex function if and only if it is cyclically monotone rockafellar1970convex. A function $h$ is said to be cyclically monotone if for any sequence $x_0,\dots,x_n$ in its domain,

equation[equation omitted — 84 chars of source]

where $x_{n+1}=x_0$. For $n=1$, Eq. (ref) reduces to Eq. (ref). With this preliminary discussion in place, the identification result in higher dimensions can be stated as follows.

customassu{2}[Co-cyclic monotonocity] Functions $h_{s,d}(\cdot;t_0)$ and $h_{s,d}(\cdot;t_1)$ are co-cyclically monotone. That is, for any sequence $u_0,\dots, u_n$ in their common domain, \[\sum_{i=0}^n\langle h_{s,d}(u_i;t_0), h_{s,d}(u_i;t_1)-h_{s,d}(u_{i+1};t_1)\rangle\geq0,\] where $u_{n+1}=u_0$.
restatable{theorem}{thmhighd} Under assumptions (ref), (ref), (ref), (ref), and (ref) - (ref), the probability measure over the missing counterfactual, i.e., $\mu_{s_1,d_1}$, is identified as \begin{equation} \mu_{s_1,d_1} = (T^*\circ T_{s_1,d_0})_\#\eta_{s_1,d_1}, \end{equation} where $T_{s,d}$ is the Brenier map that pushes forward $\eta_{s,d}$ to $\mu_{s, d}$ for every $s\in\{s_0,s_1\}, d\in\{d_0,d_1\}$, and $T^*$ is the Brenier map that pushes forward $(T_{s_0,d_0\#}\eta_{s_0,d_1})$ to $\mu_{s_0,d_1}$.

Theorem (ref) is in essence a generalization of Theorem (ref). To see this, note that the Brenier map that pushes forward $\eta_{s,d}$ to $\mu_{s,d}$ in one dimension is precisely $F^{-1}_{Y^0(t_1)\vert s,d}\circ F_{Y^0(t_0)\vert s,d}$. However, the co-cyclic monotonicity assumption (ref) is not as easy to interpret as (ref). To address this challenge, we follow the proposition of torous2021optimal based on the following result.

proposition[saks2005weak] Let $K$ and $F$ be a convex and a finite subset of $\mathbb{R}^d$ respectively. A function $T:K\to F$ is cyclically monotone if and only if it is monotone.

Proposition (ref) implies that in our setting, if the densities $\eta_{s,d}$ are supported on convex sets and $Y^0(t_1)$ is finite-valued in every subgroup, then (ref) reduces to co-monotonicity (ref). In other words, restricting $\eta_{s,d}$ and $\mu_{s,d}$ densities to be defined over convex and finite sets, respectively, Assumption (ref) of Theorem (ref) can be replaced by the easily interpretable assumption of Eq. (ref).

figure*[figure* omitted — 628 chars of source]

Simulation Studies

Our empirical analysis is structured into two main sections\footnote{ The code to reproduce the results of this paper are accessible at https://github.com/SinaAkbarii/Triple-Changes.}. In the first part, we compare the estimation error of the triple changes estimator and DiD, triple difference, and CiC estimators using synthetically generated datasets, where we know the ground truth: i.e., the treatment effect on the treated. To evaluate the performance of each estimator, we focus on the relative bias metric, defined as $\varepsilon=\big\vert 1-\frac{\widehat{\tau}}{\tau}\big\vert,$ where $\tau$ is the average treatment effect on the treated, given by (ref), and $\widehat{\tau}$ represents the estimate provided by the respective estimator. In the second part, we apply our estimator to data from National Survey on Children's Health (NSCH)\footnote{https://mchb.hrsa.gov/national-survey-childrens-health-questionnaires-datasets-supporting-documents} to assess the effect of Medicaid expansion under Affordable Care Act (ACA) in the united states on preventive care for children.

Synthetic data

We begin with a simple linear model where the latent variable $U$ given $s,d$ is sampled from a Gaussian distribution with mean $\nu_{s,d}$ and variance $1$ (see Appendix (ref) for a complete table of parameters $\nu_{s,d}$, as well as further details of our experiment setup,) and the production functions are defined as

equation[equation omitted — 102 chars of source]

The actual outcome at time $t_1$ in the treated group $(s_1,d_1)$ is also sampled from a Gaussian distribution with mean $2.75$ and variance $1$.

CiC and triple changes estimators use estimates of the cumulative density functions and their inverses. We implemented two versions of each of these estimators, namely a non-parametric estimator using the empirical estimator of density functions (see Eq. (ref) and (ref),) and a model-based version using the maximum likelihood estimator given a class of distributions. The class of distributions in this section are specified correctly, i.e., the true model lies in the class. See appendix (ref) for a setting where this class is misspecified.

In Figure (ref), the triple difference and triple changes estimators are denoted by DDD and CCC, respectively. As depicted by Figure (ref), DiD and CiC estimators exhibit persistent bias, not converging to zero even with increasing sample size. In contrast, the biases of triple difference and triple changes estimators approach zero with as the sample size grows. Furthermore, since the triple difference estimator only relies on empirical averages, it shows slightly lower bias in small sample sizes.

To add non-linearities to the previous model, we modified the data generating mechanism, including two of the production functions. Specifically, for $(s_0,d_1,t_1)$, and $(s_1,d_1,t_1)$, the production functions were modified to \[h_{s,d}(u;t) = 0.1\exp\Big(2u+\big(\frac{1+s}{4}+\frac{d-0.5}{2}\big) t\Big),\] whereas the other six groups were generated according to Eq. (ref). The outcomes for this model are illustrated in Figure (ref). Notably, under the nonlinear model, the triple-difference estimator yields biased estimates, while the triple-changes estimators remain (asymptotically) unbiased. Furthermore, owing to the complexity of the model, empirical estimators of density functions exhibit more bias compared to their MLE-based counterparts, particularly in low-sample regimes.

Application to NSCH data

The Medicaid expansion under Affordable Care Act was designed to extend Medicaid coverage to more low-income citizens in the US. This expansion raised the income threshold for Medicaid eligibility to 138% of the federal poverty level (FPL). While Louisiana ($s_1$) adopted this expansion in 2016, the neighboring states of Mississippi and Texas ($s_0$) are yet to do so. We utilized publicly available anonymous NSCH data for the years 2016 ($t_0$) and 2017 ($t_1$) to assess the impact of Medicaid expansion on children's access to preventive healthcare, specifically analyzing the change in the frequency of doctor visits. Individuals with FPL $\leq100\%$ were classified as eligible for the expansion ($d_1$), while those with FPL $\geq140\%$ were considered not eligible ($d_0$). Data within the uncertainty margin between these thresholds was excluded.

Applying the triple difference and triple changes estimators on the frequency of doctor visits for preventive purposes with 1000 bootstraps resulted in means of $0.170$ and $0.145$, respectively, with $90\%$ confidence intervals of $[-0.004, 0.344]$ and $[-0.01, 0.331]$, respectively. Our findings suggest that the adoption of Medicaid expansion in Louisiana increased the likelihood of children undergoing an annual preventive care visit, compared to those in non-expansion states. This aligns with the conclusion of roy2020impact.

Concluding Remarks

We formally presented the triple changes estimator for assessing the treatment effect on the treated in observational studies. We derived a set of necessary assumptions for point identification and discussed partial identification results under relaxed versions of these assumptions. Future research avenues could include exploring partial identification under alternative assumptions, investigating extensions to time series settings, and conducting statistical analyses to enhance the robustness and applicability of the estimator.