I’ve been inspired by two recent talks describing the challenges of using instrumental variable (IV) methods. IV methods are used to estimate the causal effects of an exposure or intervention when there is unmeasured confounding. This estimated causal effect is very specific: the complier average causal effect (CACE). But, the CACE is just one of several possible causal estimands that we might be interested in. For example, there’s the average causal effect (ACE) that represents a population average (not just based the subset of compliers). Or there’s the average causal effect for the exposed or treated (ACT) that allows for the fact that the exposed could be different from the unexposed.
I thought it would be illuminating to analyze a single data set using different causal inference methods, including IV as well as propensity score matching and inverse probability weighting. Each of these methods targets different causal estimands, which may or may not be equivalent depending on the subgroup-level causal effects and underlying population distribution of those subgroups.
This is the first of a two-part post. In this first part, I am focusing entirely on the data generation process (DGP). In the follow-up, I will get to the model estimation.
Underlying assumptions of the DGP
Since the motivation here is instrumental variable analysis, it seems natural that the data generation process include a possible instrument. (Once again, I am going to refer to elsewhere in case you want more details on the theory and estimation of IV models. Here is an excellent in-depth tutorial by Baiocchi et al that provides great background. I’ve even touched on the topic of CACE in an earlier series of posts. Certainly, there is no lack of discussion on this topic, as a quick search around the internet will make readily obvious.)
The figure below is a variation on the directed acyclic graph (DAG) that is often very useful in laying out causal assumptions of a DGP. This particular figure is a type of SWIG: single world intervention graph. SWIGs, developed by Robins and Richardson, fuse the worlds of potential outcomes and DAGs.
Important things to note here:
There is an instrumental variable \(A\) that has a direct causal relationship only to the exposure of interest, \(T\). If the exposure is a particular medical intervention, think of the instrument as some kind of encouragement to get that treatment. Some people get the encouragement, others don’t - though on average folks who are encouraged are no different from folks who are not (at least not in ways that relate to the outcome.)
There is a confounder \(U\), possibly unmeasured, that is related both to potential outcomes and the exposure, but not to the encouragement (the instrument)! In the example below, we conceive of \(U\) as an underlying health status.
Exposure variable \(T\) (that, in this case, is binary, just to keep things simpler) indicates whether a person gets the treatment or not.
Each individual will have two potential treatments \(T^0\) and \(T^1\), where \(T^0\) is the treatment when there is no encouragement (i.e. A = 0), and \(T^1\) is the treatment when \(A = 1\). For any individual, we actually only observe one of these treatments (depending on the actual value of \(A\). The population of interest consists of always-takers, compliers, and never-takers. Never-takers always reject the treatment regardless of whether or not they get encouragement - that is, \(T^0 = T^1 = 0\). Compliers only seek out the treatment when they are encouraged, otherwise they don’t: \(T^0 = 0\) and \(T^1 = 1\). And always-takers always (of course) seek out the treatment: \(T^0 = T^1 = 1\). (In order for the model to be identifiable, we need to make a not-so-crazy assumption that there are no so-called deniers, where \(T^0 = 1\) and \(T^1 = 0\).) An individual may have a different complier status depending on the instrument and exposure (i.e., one person might be a never-taker in one scenario but a complier in another). In this simulation, larger values of the confounder \(U\) will increase \(P(T^a = 1)\) for both \(a \in (0,1)\).
Each individual will have two potential outcomes, only one of which is observed. \(Y_i^0\) is the outcome for person \(i\) when they are unexposed or do not receive the treatment. \(Y_i^1\) is the outcome for that same person when they are exposed or do receive the treatment. In this case, the confounder \(U\) can affect the potential outcomes. (This diagram is technically a SWIT, which is template, since I have generically referred to the potential treatment \(T^a\) and potential outcome \(Y^t\).)
Not shown in this diagram are the observed \(T_i\) and \(Y_i\); we assume that \(T_i = (T_i^a | A = a)\) and \(Y_i = (Y_i^t | T = t)\)
Also not shown on the graph is the causal estimand of an exposure for individual \(i\), which can be defined as \(CE_i \equiv Y^1_i - Y^0_i\). We can calculate the average causal effect, \(E[CE]\), for the sample as a whole as well as for subgroups.
DGP for potential outcomes
The workhorse of this data generating process is a logistic sigmoid function that represents the mean potential outcome \(Y^t\) at each value of \(u\). This allows us to easily generate homogeneous or heterogeneous causal effects. The function has four parameters, \(M\), \(\gamma\), \(\delta\), and \(\alpha\):
\[
Y^t = f(u) = M/[1 + exp(-\gamma(u - \delta))] + \alpha,
\]
where \(M\) is the maximum of the function (assuming the minimum is \(0\)), \(\gamma\) is the steepness of the curve, \(\delta\) is the inflection point of the curve, and \(\alpha\) is a vertical shift of the entire curve. This function is easily implemented in R:
fYt