AZCausal Integration
Amazon's AZCausal library provides the
functionality to fit synthetic control and difference-in-difference models to your
data. Integrating the synthetic data generating process of causal_validation with
AZCausal is trivial, as we show in this notebook. To start, we'll simulate a toy
dataset.
from azcausal.estimators.panel.sdid import SDID
import scipy.stats as st
@@ -469,6 +515,14 @@ AZCausal Integration
)
from causal_validation.transforms.parameter import UnitVaryingParameter
cfg = Config(
n_control_units=10,
n_pre_intervention_timepoints=60,
@@ -480,18 +534,65 @@ AZCausal Integration
data = linear_trend(simulate(cfg))
plot(data)
<Axes: xlabel='Time', ylabel='Observed'>
-
+will inflate the treated group's observations in the post-intervention window.
+TRUE_EFFECT = 0.05
effect = StaticEffect(effect=TRUE_EFFECT)
inflated_data = effect(data)
plot(inflated_data)
<Axes: xlabel='Time', ylabel='Observed'>
-
+Fitting a model
We now have some very toy data on which we may apply a model. For this demonstration we shall use the Synthetic Difference-in-Differences model implemented in AZCausal; @@ -500,13 +601,28 @@
Fitting a model
AZCausal. Through the.to_azcausal() method implemented here, this is
straightforward to achieve. Once we have a AZCausal compatible dataset, the modelling
is very simple by virtue of the clean design of AZCausal.
+panel = inflated_data.to_azcausal()
model = SDID()
result = model.fit(panel)
print(f"Delta: {TRUE_EFFECT - result.effect.percentage().value / 100}")
print(result.summary(title="Synthetic Data Experiment"))
Delta: -2.3592239273284576e-16
+
+
+effect. However, given the simplicity of the data, this is not surprising. With the functionality within this package though we can easily construct more complex datasets in effort to fully stress-test any new model and identify its limitations.
@@ -545,6 +670,13 @@Fitting a model
$\delta_{t, n}$ is 5% when $n=1$ and $t\geq 60$ and 0 otherwise. Meanwhile, $\mathbf{Y}$ is the matrix of observations, long in the number of time points and wide in the number of units. +cfg = Config(
n_control_units=10,
n_pre_intervention_timepoints=60,
@@ -565,18 +697,53 @@ Fitting a model
data = effect(periodic(linear_trend(simulate(cfg))))
plot(data)
<Axes: xlabel='Time', ylabel='Observed'>
-
+time we see that the delta between the estaimted and true effect is much larger than before.
+panel = data.to_azcausal()
model = SDID()
result = model.fit(panel)
print(f"Delta: {100*(TRUE_EFFECT - result.effect.percentage().value / 100): .2f}%")
print(result.summary(title="Synthetic Data Experiment"))
Delta: 1.71%
+
+
+Data Synthesis
In this notebook we'll demonstrate how causal-validation can be used to simulate
synthetic datasets. We'll start with very simple data to which a static treatment
effect may be applied. From there, we'll build up to complex datasets. Along the way,
we'll show how reproducibility can be ensured, plots can be generated, and unit-level
parameters may be specified.
from itertools import product
import matplotlib.pyplot as plt
@@ -597,9 +643,29 @@ Data Synthesis
)
from causal_validation.transforms.parameter import UnitVaryingParameter
Simulating a Dataset
+then invoking the simulate function. Once simulated, we may visualise the data
through the plot function.
cfg = Config(
n_control_units=10,
n_pre_intervention_timepoints=60,
@@ -610,14 +676,41 @@ Simulating a Dataset
data = simulate(cfg)
plot(data)
<Axes: xlabel='Time', ylabel='Observed'>
-
+Controlling baseline behaviour
We observe that we have 10 control units, each of which were sampled from a Gaussian distribution with mean 20 and scale 0.2. Had we wished for our underlying observations to have more or less noise, or to have a different global mean, then we can simply specify that through the config file.
+means = [10, 50]
scales = [0.1, 0.5]
@@ -633,12 +726,34 @@ Controlling baseline behaviour
data = simulate(cfg)
plot(data, ax=ax, title=f"Mean: {m}, Scale: {s}")
+Reproducibility
In the above four panels, we can see that whilst the mean and scale of the underlying data generating process is varying, the functional form of the data is the same. This is by design to ensure that data sampling is reproducible. To sample a new dataset, you may either change the underlying seed in the config file.
+cfg = Config(
n_control_units=10,
n_pre_intervention_timepoints=60,
@@ -646,14 +761,50 @@ Reproducibility
seed=42,
)
Reusing the same config file across simulations
+fig, axes = plt.subplots(ncols=2, figsize=(10, 3))
for ax in axes:
data = simulate(cfg)
plot(data, ax=ax)
+Or manually specifying and passing your own pseudorandom number generator key
+
rng = np.random.RandomState(42)
@@ -662,67 +813,232 @@ Reproducibility
data = simulate(cfg, key=rng)
plot(data, ax=ax)
+Simulating an effect
In the data we have seen up until now, the treated unit has been drawn from the same
data generating process as the control units. However, it can be helpful to also
inflate the treated unit to observe how well our model can recover the the true
treatment effect. To do this, we simply compose our dataset with an Effect object.
In the below, we shall inflate our data by 2%.
effect = StaticEffect(effect=0.02)
inflated_data = effect(data)
fig, (ax0, ax1) = plt.subplots(ncols=2, figsize=(10, 3))
plot(data, ax=ax0, title="Original data")
plot(inflated_data, ax=ax1, title="Inflated data")
<Axes: title={'center': 'Inflated data'}, xlabel='Time', ylabel='Observed'>
-
+More complex generation processes
The example presented above shows a very simple stationary data generation process. However, we may make our example more complex by including a non-stationary trend to the data.
+trend_term = Trend(degree=1, coefficient=0.1)
data_with_trend = effect(trend_term(data))
plot(data_with_trend)
<Axes: xlabel='Time', ylabel='Observed'>
-
+trend_term = Trend(degree=2, coefficient=0.0025)
data_with_trend = effect(trend_term(data))
plot(data_with_trend)
<Axes: xlabel='Time', ylabel='Observed'>
-
+We may also include periodic components in our data
+periodicity = Periodic(amplitude=2, frequency=6)
perioidic_data = effect(periodicity(trend_term(data)))
plot(perioidic_data)
<Axes: xlabel='Time', ylabel='Observed'>
-
+Unit-level parameterisation
+sampling_dist = norm(0.0, 1.0)
intercept = UnitVaryingParameter(sampling_dist=sampling_dist)
trend_term = Trend(degree=1, intercept=intercept, coefficient=0.1)
data_with_trend = effect(trend_term(data))
plot(data_with_trend)
<Axes: xlabel='Time', ylabel='Observed'>
-
+sampling_dist = poisson(2)
frequency = UnitVaryingParameter(sampling_dist=sampling_dist)
p = Periodic(frequency=frequency)
plot(p(data))
<Axes: xlabel='Time', ylabel='Observed'>
-
+Conclusions
In this notebook we have shown how one can define their model's true underlying data generating process, starting from simple white-noise samples through to more complex @@ -731,6 +1047,9 @@
Conclusions
Amazon's own AZCausal library to compare the effect estimated by a model with the true effect of the underlying data generating process. A link to this notebook is here. +Placebo Testing
A placebo test is an approach to assess the validity of a causal model by checking if the effect can truly be attributed to the treatment, or to other spurious factors. A @@ -518,6 +557,13 @@
Placebo Testing
reliable. Placebo testing is thus a critical step to ensure the robustness of findings in RCTs. In this notebook, we demonstrate how a placebo test can be conducted incausal-validation.
+from azcausal.core.error import JackKnife
from azcausal.estimators.panel.did import DID
from azcausal.estimators.panel.sdid import SDID
@@ -531,13 +577,37 @@ Placebo Testing
from causal_validation.plotters import plot
from causal_validation.validation.placebo import PlaceboTest
/home/runner/.local/share/hatch/env/virtual/causal-validation/CYBYs5D-/docs/lib/python3.10/site-packages/pandera/engines/pandas_engine.py:67: UserWarning: Using typeguard < 3. Generic types like List[TYPE], Dict[TYPE, TYPE] will only validate the first element in the collection.
+
+
+Data simulation
To demonstrate a placebo test, we must first simulate some data. For the purposes of illustration, we'll simulate a very simple dataset containing 10 control units where each unit has 60 pre-intervention observations, and 30 post-intervention observations.
+cfg = Config(
n_control_units=10,
n_pre_intervention_timepoints=60,
@@ -550,16 +620,50 @@ Data simulation
data = effect(simulate(cfg))
plot(data)
<Axes: xlabel='Time', ylabel='Observed'>
-
+Model
We'll now define our model. To do this, we'll use the synthetic
difference-in-differences implementation of AZCausal. This implementation, along with
any other model from AZCausal, can be neatly wrapped up in our AZCausalWrapper to
make fitting and effect estimation simpler.
model = AZCausalWrapper(model=SDID(), error_estimator=JackKnife())
Placebo Test Results
Now that we have a dataset and model defined, we may conduct our placebo test. With 10 control units, the test will estimate 10 individual effects; 1 per control unit when @@ -571,32 +675,91 @@
Placebo Test Results
Accordingly, the p-value attains a value of 0.5, indicating that we have insufficient evidence to reject the null hypothesis and we, therefore, have no evidence to suggest that there is bias within this particular setup. +result = PlaceboTest(model, data).execute()
result.summary()
Output()
-Model Comparison
We can also use the results of a placebo test to compare two or more models. Using
causal-validation, this is as simple as supplying a series of models to the placebo
test and comparing their outputs. To demonstrate this, we will compare the previously
used synthetic difference-in-differences model with regular difference-in-differences.
did_model = AZCausalWrapper(model=DID())
PlaceboTest([model, did_model], data).execute().summary()
Output()
-Causal Validation is a library designed to validate and test your causal models. To achieve this, we provide functionality to simulate causal data, and vaildate your model through a placebo test.
"},{"location":"#data-synthesis","title":"Data Synthesis","text":"Data Synthesis in Causal Validation is a fully composable process whereby a set of functions are sequentially applied to a dataset. At some point in this process we also induce a treatment effect. Any of these functions can be parameterised to either have constant parameter values across all control units, or a value that varies across parameters. To see this, consider the below example where we simulate a dataset whose trend varies across each of the 10 control units.
from causal_validation import Config, simulate\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.plotters import plot\nfrom causal_validation.transforms import Trend, Periodic\nfrom causal_validation.transforms.parameter import UnitVaryingParameter\nfrom scipy.stats import norm\n\ncfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n)\n\n# Simulate the base observation\nbase_data = simulate(cfg)\n\n# Apply a linear trend with unit-varying intercept\nintercept = UnitVaryingParameter(sampling_dist = norm(0, 1))\ntrend_component = Trend(degree=1, coefficient=0.1, intercept=intercept)\ntrended_data = trend_component(base_data)\n\n# Simulate a 5% lift in the treated unit's post-intervention data\neffect = StaticEffect(0.05)\ninflated_data = effect(trended_data)\nOnce a dataset has been synthesised, we may wish to validate our model using a placebo test. In Causal Validation this is straightforward and can be accomplished in combination with AZCausal by the following.
from azcausal.estimators.panel.sdid import SDID\nfrom causal_validation.validation.placebo import PlaceboTest\n\nmodel = AZCausalWrapper(model=SDID())\nresult = PlaceboTest(model, inflated_data).execute()\nresult.summary()\nAmazon's AZCausal library provides the functionality to fit synthetic control and difference-in-difference models to your data. Integrating the synthetic data generating process of causal_validation with AZCausal is trivial, as we show in this notebook. To start, we'll simulate a toy dataset.
from azcausal.estimators.panel.sdid import SDID\nimport scipy.stats as st\n\nfrom causal_validation import (\n Config,\n simulate,\n)\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.plotters import plot\nfrom causal_validation.transforms import (\n Periodic,\n Trend,\n)\nfrom causal_validation.transforms.parameter import UnitVaryingParameter\ncfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=123,\n)\n\nlinear_trend = Trend(degree=1, coefficient=0.05)\ndata = linear_trend(simulate(cfg))\nplot(data)\n<Axes: xlabel='Time', ylabel='Observed'>\nwill inflate the treated group's observations in the post-intervention window.
TRUE_EFFECT = 0.05\neffect = StaticEffect(effect=TRUE_EFFECT)\ninflated_data = effect(data)\nplot(inflated_data)\n<Axes: xlabel='Time', ylabel='Observed'>\nWe now have some very toy data on which we may apply a model. For this demonstration we shall use the Synthetic Difference-in-Differences model implemented in AZCausal; however, the approach shown here will work for any model implemented in AZCausal. To achieve this, we must first coerce the data into a format that is digestible for AZCausal. Through the .to_azcausal() method implemented here, this is straightforward to achieve. Once we have a AZCausal compatible dataset, the modelling is very simple by virtue of the clean design of AZCausal.
panel = inflated_data.to_azcausal()\nmodel = SDID()\nresult = model.fit(panel)\nprint(f\"Delta: {TRUE_EFFECT - result.effect.percentage().value / 100}\")\nprint(result.summary(title=\"Synthetic Data Experiment\"))\nDelta: -2.3592239273284576e-16\n\u256d\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256e\n| Synthetic Data Experiment |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Panel |\n| Time Periods: 90 (60/30) total (pre/post) |\n| Units: 11 (10/1) total (contr/treat) |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| ATT |\n| Effect: 1.1858 |\n| Observed: 24.90 |\n| Counter Factual: 23.72 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Percentage |\n| Effect: 5.0000 |\n| Observed: 105.00 |\n| Counter Factual: 100.00 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Cumulative |\n| Effect: 35.57 |\n| Observed: 747.03 |\n| Counter Factual: 711.46 |\n\u2570\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256f\neffect. However, given the simplicity of the data, this is not surprising. With the functionality within this package though we can easily construct more complex datasets in effort to fully stress-test any new model and identify its limitations.
To achieve this, we'll simulate 10 control units, 60 pre-intervention time points, and 30 post-intervention time points according to the following process: $$ \\begin{align} \\mu_{n, t} & \\sim\\mathcal{N}(20, 0.5^2)\\ \\alpha_{n} & \\sim \\mathcal{N}(0, 1^2)\\ \\beta_{n} & \\sim \\mathcal{N}(0.05, 0.01^2)\\ \\nu_n & \\sim \\mathcal{N}(1, 1^2)\\ \\gamma_n & \\sim \\operatorname{Student-t}{10}(1, 1^2)\\ \\mathbf{Y}{n, t} & = \\mu_{n, t} + \\alpha_{n} + \\beta_{n}t + \\nu_n\\sin\\left(3\\times 2\\pi t + \\gamma\\right) + \\delta_{t, n} \\end{align} $$ where the true treatment effect $\\delta_{t, n}$ is 5% when $n=1$ and $t\\geq 60$ and 0 otherwise. Meanwhile, $\\mathbf{Y}$ is the matrix of observations, long in the number of time points and wide in the number of units.
cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n global_mean=20,\n global_scale=1,\n seed=123,\n)\n\nintercept = UnitVaryingParameter(sampling_dist=st.norm(loc=0.0, scale=1))\ncoefficient = UnitVaryingParameter(sampling_dist=st.norm(loc=0.05, scale=0.01))\nlinear_trend = Trend(degree=1, coefficient=coefficient, intercept=intercept)\n\namplitude = UnitVaryingParameter(sampling_dist=st.norm(loc=1.0, scale=2))\nshift = UnitVaryingParameter(sampling_dist=st.t(df=10))\nperiodic = Periodic(amplitude=amplitude, shift=shift, frequency=3)\n\ndata = effect(periodic(linear_trend(simulate(cfg))))\nplot(data)\n<Axes: xlabel='Time', ylabel='Observed'>\ntime we see that the delta between the estaimted and true effect is much larger than before.
panel = data.to_azcausal()\nmodel = SDID()\nresult = model.fit(panel)\nprint(f\"Delta: {100*(TRUE_EFFECT - result.effect.percentage().value / 100): .2f}%\")\nprint(result.summary(title=\"Synthetic Data Experiment\"))\nDelta: 1.71%\n\u256d\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256e\n| Synthetic Data Experiment |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Panel |\n| Time Periods: 90 (60/30) total (pre/post) |\n| Units: 11 (10/1) total (contr/treat) |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| ATT |\n| Effect: 0.728265 |\n| Observed: 22.88 |\n| Counter Factual: 22.15 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Percentage |\n| Effect: 3.2874 |\n| Observed: 103.29 |\n| Counter Factual: 100.00 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Cumulative |\n| Effect: 21.85 |\n| Observed: 686.44 |\n| Counter Factual: 664.59 |\n\u2570\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256f\nIn this notebook we'll demonstrate how causal-validation can be used to simulate synthetic datasets. We'll start with very simple data to which a static treatment effect may be applied. From there, we'll build up to complex datasets. Along the way, we'll show how reproducibility can be ensured, plots can be generated, and unit-level parameters may be specified.
from itertools import product\n\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom scipy.stats import (\n norm,\n poisson,\n)\n\nfrom causal_validation import (\n Config,\n simulate,\n)\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.plotters import plot\nfrom causal_validation.transforms import (\n Periodic,\n Trend,\n)\nfrom causal_validation.transforms.parameter import UnitVaryingParameter\nthen invoking the simulate function. Once simulated, we may visualise the data through the plot function.
cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=123,\n)\n\ndata = simulate(cfg)\nplot(data)\n<Axes: xlabel='Time', ylabel='Observed'>\nWe observe that we have 10 control units, each of which were sampled from a Gaussian distribution with mean 20 and scale 0.2. Had we wished for our underlying observations to have more or less noise, or to have a different global mean, then we can simply specify that through the config file.
means = [10, 50]\nscales = [0.1, 0.5]\n\nfig, axes = plt.subplots(ncols=2, nrows=2, figsize=(10, 6), tight_layout=True)\nfor (m, s), ax in zip(product(means, scales), axes.ravel(), strict=False):\n cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n global_mean=m,\n global_scale=s,\n )\n data = simulate(cfg)\n plot(data, ax=ax, title=f\"Mean: {m}, Scale: {s}\")\nIn the above four panels, we can see that whilst the mean and scale of the underlying data generating process is varying, the functional form of the data is the same. This is by design to ensure that data sampling is reproducible. To sample a new dataset, you may either change the underlying seed in the config file.
cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=42,\n)\nReusing the same config file across simulations
fig, axes = plt.subplots(ncols=2, figsize=(10, 3))\nfor ax in axes:\n data = simulate(cfg)\n plot(data, ax=ax)\nOr manually specifying and passing your own pseudorandom number generator key
\nrng = np.random.RandomState(42)\n\nfig, axes = plt.subplots(ncols=2, figsize=(10, 3))\nfor ax in axes:\n data = simulate(cfg, key=rng)\n plot(data, ax=ax)\nIn the data we have seen up until now, the treated unit has been drawn from the same data generating process as the control units. However, it can be helpful to also inflate the treated unit to observe how well our model can recover the the true treatment effect. To do this, we simply compose our dataset with an Effect object. In the below, we shall inflate our data by 2%.
effect = StaticEffect(effect=0.02)\ninflated_data = effect(data)\nfig, (ax0, ax1) = plt.subplots(ncols=2, figsize=(10, 3))\nplot(data, ax=ax0, title=\"Original data\")\nplot(inflated_data, ax=ax1, title=\"Inflated data\")\n<Axes: title={'center': 'Inflated data'}, xlabel='Time', ylabel='Observed'>\nThe example presented above shows a very simple stationary data generation process. However, we may make our example more complex by including a non-stationary trend to the data.
trend_term = Trend(degree=1, coefficient=0.1)\ndata_with_trend = effect(trend_term(data))\nplot(data_with_trend)\n<Axes: xlabel='Time', ylabel='Observed'>\ntrend_term = Trend(degree=2, coefficient=0.0025)\ndata_with_trend = effect(trend_term(data))\nplot(data_with_trend)\n<Axes: xlabel='Time', ylabel='Observed'>\nWe may also include periodic components in our data
periodicity = Periodic(amplitude=2, frequency=6)\nperioidic_data = effect(periodicity(trend_term(data)))\nplot(perioidic_data)\n<Axes: xlabel='Time', ylabel='Observed'>\nsampling_dist = norm(0.0, 1.0)\nintercept = UnitVaryingParameter(sampling_dist=sampling_dist)\ntrend_term = Trend(degree=1, intercept=intercept, coefficient=0.1)\ndata_with_trend = effect(trend_term(data))\nplot(data_with_trend)\n<Axes: xlabel='Time', ylabel='Observed'>\nsampling_dist = poisson(2)\nfrequency = UnitVaryingParameter(sampling_dist=sampling_dist)\n\np = Periodic(frequency=frequency)\nplot(p(data))\n<Axes: xlabel='Time', ylabel='Observed'>\nIn this notebook we have shown how one can define their model's true underlying data generating process, starting from simple white-noise samples through to more complex example with periodic and temporal components, perhaps containing unit-level variation. In a follow-up notebook, we show how these datasets may be integrated with Amazon's own AZCausal library to compare the effect estimated by a model with the true effect of the underlying data generating process. A link to this notebook is here.
"},{"location":"_examples/placebo_test/","title":"Placebo Testing","text":"A placebo test is an approach to assess the validity of a causal model by checking if the effect can truly be attributed to the treatment, or to other spurious factors. A placebo test is conducted by iterating through the set of control units and at each iteration, replacing the treated unit by one of the control units and measuring the effect. If the model detects a significant effect, then it suggests potential bias or omitted variables in the analysis, indicating that the causal inference is flawed.
A successful placebo test will show no statistically significant results and we may then conclude that the estimated effect can be attributed to the treatment and not driven by confounding factors. Conversely, a failed placebo test, which shows significant results, suggests that the identified treatment effect may not be reliable. Placebo testing is thus a critical step to ensure the robustness of findings in RCTs. In this notebook, we demonstrate how a placebo test can be conducted in causal-validation.
from azcausal.core.error import JackKnife\nfrom azcausal.estimators.panel.did import DID\nfrom azcausal.estimators.panel.sdid import SDID\n\nfrom causal_validation import (\n Config,\n simulate,\n)\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.models import AZCausalWrapper\nfrom causal_validation.plotters import plot\nfrom causal_validation.validation.placebo import PlaceboTest\n/home/runner/.local/share/hatch/env/virtual/causal-validation/CYBYs5D-/docs/lib/python3.10/site-packages/pandera/engines/pandas_engine.py:67: UserWarning: Using typeguard < 3. Generic types like List[TYPE], Dict[TYPE, TYPE] will only validate the first element in the collection.\n warnings.warn(\nTo demonstrate a placebo test, we must first simulate some data. For the purposes of illustration, we'll simulate a very simple dataset containing 10 control units where each unit has 60 pre-intervention observations, and 30 post-intervention observations.
cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=123,\n)\n\nTRUE_EFFECT = 0.05\neffect = StaticEffect(effect=TRUE_EFFECT)\ndata = effect(simulate(cfg))\nplot(data)\n<Axes: xlabel='Time', ylabel='Observed'>\nWe'll now define our model. To do this, we'll use the synthetic difference-in-differences implementation of AZCausal. This implementation, along with any other model from AZCausal, can be neatly wrapped up in our AZCausalWrapper to make fitting and effect estimation simpler.
model = AZCausalWrapper(model=SDID(), error_estimator=JackKnife())\nNow that we have a dataset and model defined, we may conduct our placebo test. With 10 control units, the test will estimate 10 individual effects; 1 per control unit when it is mocked as the treated group. With those 10 effects, the routine will then produce the mean estimated effect, along with the standard deviation across the estimated effect, the effect's standard error, and the p-value that corresponds to the null-hypothesis test that the effect is 0.
In the below, we see that expected estimated effect is small at just 0.08. Accordingly, the p-value attains a value of 0.5, indicating that we have insufficient evidence to reject the null hypothesis and we, therefore, have no evidence to suggest that there is bias within this particular setup.
result = PlaceboTest(model, data).execute()\nresult.summary()\nOutput()\n\n| Model | Effect | Standard Deviation | Standard Error | p-value |\n|-------|--------|--------------------|----------------|---------|\n| SDID | 0.0851 | 0.4079 | 0.129 | 0.5472 |\n\n"},{"location":"_examples/placebo_test/#model-comparison","title":"Model Comparison","text":"
We can also use the results of a placebo test to compare two or more models. Using causal-validation, this is as simple as supplying a series of models to the placebo test and comparing their outputs. To demonstrate this, we will compare the previously used synthetic difference-in-differences model with regular difference-in-differences.
did_model = AZCausalWrapper(model=DID())\nPlaceboTest([model, did_model], data).execute().summary()\nOutput()\n\n| Model | Effect | Standard Deviation | Standard Error | p-value |\n|-------|--------|--------------------|----------------|---------|\n| SDID | 0.0851 | 0.4079 | 0.129 | 0.5472 |\n| DID | 0.0002 | 0.2818 | 0.0891 | 0.9982 |\n\n"}]} \ No newline at end of file +{"config":{"lang":["en"],"separator":"[\\s\\-]+","pipeline":["stopWordFilter"]},"docs":[{"location":"","title":"Welcome to Causal Validation","text":"
Causal Validation is a library designed to validate and test your causal models. To achieve this, we provide functionality to simulate causal data, and vaildate your model through a placebo test.
"},{"location":"#data-synthesis","title":"Data Synthesis","text":"Data Synthesis in Causal Validation is a fully composable process whereby a set of functions are sequentially applied to a dataset. At some point in this process we also induce a treatment effect. Any of these functions can be parameterised to either have constant parameter values across all control units, or a value that varies across parameters. To see this, consider the below example where we simulate a dataset whose trend varies across each of the 10 control units.
from causal_validation import Config, simulate\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.plotters import plot\nfrom causal_validation.transforms import Trend, Periodic\nfrom causal_validation.transforms.parameter import UnitVaryingParameter\nfrom scipy.stats import norm\n\ncfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n)\n\n# Simulate the base observation\nbase_data = simulate(cfg)\n\n# Apply a linear trend with unit-varying intercept\nintercept = UnitVaryingParameter(sampling_dist = norm(0, 1))\ntrend_component = Trend(degree=1, coefficient=0.1, intercept=intercept)\ntrended_data = trend_component(base_data)\n\n# Simulate a 5% lift in the treated unit's post-intervention data\neffect = StaticEffect(0.05)\ninflated_data = effect(trended_data)\nOnce a dataset has been synthesised, we may wish to validate our model using a placebo test. In Causal Validation this is straightforward and can be accomplished in combination with AZCausal by the following.
from azcausal.estimators.panel.sdid import SDID\nfrom causal_validation.validation.placebo import PlaceboTest\n\nmodel = AZCausalWrapper(model=SDID())\nresult = PlaceboTest(model, inflated_data).execute()\nresult.summary()\nfrom azcausal.estimators.panel.sdid import SDID\nimport scipy.stats as st\n\nfrom causal_validation import (\n Config,\n simulate,\n)\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.plotters import plot\nfrom causal_validation.transforms import (\n Periodic,\n Trend,\n)\nfrom causal_validation.transforms.parameter import UnitVaryingParameter\ncfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=123,\n)\n\nlinear_trend = Trend(degree=1, coefficient=0.05)\ndata = linear_trend(simulate(cfg))\nplot(data)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n will inflate the treated group's observations in the post-intervention window.
TRUE_EFFECT = 0.05\neffect = StaticEffect(effect=TRUE_EFFECT)\ninflated_data = effect(data)\nplot(inflated_data)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n panel = inflated_data.to_azcausal()\nmodel = SDID()\nresult = model.fit(panel)\nprint(f\"Delta: {TRUE_EFFECT - result.effect.percentage().value / 100}\")\nprint(result.summary(title=\"Synthetic Data Experiment\"))\n\nDelta: -2.3592239273284576e-16\n\u256d\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256e\n| Synthetic Data Experiment |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Panel |\n| Time Periods: 90 (60/30) total (pre/post) |\n| Units: 11 (10/1) total (contr/treat) |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| ATT |\n| Effect: 1.1858 |\n| Observed: 24.90 |\n| Counter Factual: 23.72 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Percentage |\n| Effect: 5.0000 |\n| Observed: 105.00 |\n| Counter Factual: 100.00 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Cumulative |\n| Effect: 35.57 |\n| Observed: 747.03 |\n| Counter Factual: 711.46 |\n\u2570\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256f\n\n effect. However, given the simplicity of the data, this is not surprising. With the functionality within this package though we can easily construct more complex datasets in effort to fully stress-test any new model and identify its limitations.
To achieve this, we'll simulate 10 control units, 60 pre-intervention time points, and 30 post-intervention time points according to the following process: $$ \\begin{align} \\mu_{n, t} & \\sim\\mathcal{N}(20, 0.5^2)\\ \\alpha_{n} & \\sim \\mathcal{N}(0, 1^2)\\ \\beta_{n} & \\sim \\mathcal{N}(0.05, 0.01^2)\\ \\nu_n & \\sim \\mathcal{N}(1, 1^2)\\ \\gamma_n & \\sim \\operatorname{Student-t}{10}(1, 1^2)\\ \\mathbf{Y}{n, t} & = \\mu_{n, t} + \\alpha_{n} + \\beta_{n}t + \\nu_n\\sin\\left(3\\times 2\\pi t + \\gamma\\right) + \\delta_{t, n} \\end{align} $$ where the true treatment effect $\\delta_{t, n}$ is 5% when $n=1$ and $t\\geq 60$ and 0 otherwise. Meanwhile, $\\mathbf{Y}$ is the matrix of observations, long in the number of time points and wide in the number of units.
cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n global_mean=20,\n global_scale=1,\n seed=123,\n)\n\nintercept = UnitVaryingParameter(sampling_dist=st.norm(loc=0.0, scale=1))\ncoefficient = UnitVaryingParameter(sampling_dist=st.norm(loc=0.05, scale=0.01))\nlinear_trend = Trend(degree=1, coefficient=coefficient, intercept=intercept)\n\namplitude = UnitVaryingParameter(sampling_dist=st.norm(loc=1.0, scale=2))\nshift = UnitVaryingParameter(sampling_dist=st.t(df=10))\nperiodic = Periodic(amplitude=amplitude, shift=shift, frequency=3)\n\ndata = effect(periodic(linear_trend(simulate(cfg))))\nplot(data)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n time we see that the delta between the estaimted and true effect is much larger than before.
panel = data.to_azcausal()\nmodel = SDID()\nresult = model.fit(panel)\nprint(f\"Delta: {100*(TRUE_EFFECT - result.effect.percentage().value / 100): .2f}%\")\nprint(result.summary(title=\"Synthetic Data Experiment\"))\n\nDelta: 1.71%\n\u256d\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256e\n| Synthetic Data Experiment |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Panel |\n| Time Periods: 90 (60/30) total (pre/post) |\n| Units: 11 (10/1) total (contr/treat) |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| ATT |\n| Effect: 0.728265 |\n| Observed: 22.88 |\n| Counter Factual: 22.15 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Percentage |\n| Effect: 3.2874 |\n| Observed: 103.29 |\n| Counter Factual: 100.00 |\n\u251c\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2524\n| Cumulative |\n| Effect: 21.85 |\n| Observed: 686.44 |\n| Counter Factual: 664.59 |\n\u2570\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u256f\n\n"},{"location":"examples/azcausal/#azcausal-integration","title":"AZCausal Integration","text":"Amazon's AZCausal library provides the functionality to fit synthetic control and difference-in-difference models to your data. Integrating the synthetic data generating process of causal_validation with AZCausal is trivial, as we show in this notebook. To start, we'll simulate a toy dataset.
We now have some very toy data on which we may apply a model. For this demonstration we shall use the Synthetic Difference-in-Differences model implemented in AZCausal; however, the approach shown here will work for any model implemented in AZCausal. To achieve this, we must first coerce the data into a format that is digestible for AZCausal. Through the .to_azcausal() method implemented here, this is straightforward to achieve. Once we have a AZCausal compatible dataset, the modelling is very simple by virtue of the clean design of AZCausal.
from itertools import product\n\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom scipy.stats import (\n norm,\n poisson,\n)\n\nfrom causal_validation import (\n Config,\n simulate,\n)\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.plotters import plot\nfrom causal_validation.transforms import (\n Periodic,\n Trend,\n)\nfrom causal_validation.transforms.parameter import UnitVaryingParameter\nthen invoking the simulate function. Once simulated, we may visualise the data through the plot function.
cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=123,\n)\n\ndata = simulate(cfg)\nplot(data)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n means = [10, 50]\nscales = [0.1, 0.5]\n\nfig, axes = plt.subplots(ncols=2, nrows=2, figsize=(10, 6), tight_layout=True)\nfor (m, s), ax in zip(product(means, scales), axes.ravel(), strict=False):\n cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n global_mean=m,\n global_scale=s,\n )\n data = simulate(cfg)\n plot(data, ax=ax, title=f\"Mean: {m}, Scale: {s}\")\ncfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=42,\n)\nReusing the same config file across simulations
fig, axes = plt.subplots(ncols=2, figsize=(10, 3))\nfor ax in axes:\n data = simulate(cfg)\n plot(data, ax=ax)\nOr manually specifying and passing your own pseudorandom number generator key
\nrng = np.random.RandomState(42)\n\nfig, axes = plt.subplots(ncols=2, figsize=(10, 3))\nfor ax in axes:\n data = simulate(cfg, key=rng)\n plot(data, ax=ax)\neffect = StaticEffect(effect=0.02)\ninflated_data = effect(data)\nfig, (ax0, ax1) = plt.subplots(ncols=2, figsize=(10, 3))\nplot(data, ax=ax0, title=\"Original data\")\nplot(inflated_data, ax=ax1, title=\"Inflated data\")\n\n<Axes: title={'center': 'Inflated data'}, xlabel='Time', ylabel='Observed'>\n trend_term = Trend(degree=1, coefficient=0.1)\ndata_with_trend = effect(trend_term(data))\nplot(data_with_trend)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n trend_term = Trend(degree=2, coefficient=0.0025)\ndata_with_trend = effect(trend_term(data))\nplot(data_with_trend)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n We may also include periodic components in our data
periodicity = Periodic(amplitude=2, frequency=6)\nperioidic_data = effect(periodicity(trend_term(data)))\nplot(perioidic_data)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n sampling_dist = norm(0.0, 1.0)\nintercept = UnitVaryingParameter(sampling_dist=sampling_dist)\ntrend_term = Trend(degree=1, intercept=intercept, coefficient=0.1)\ndata_with_trend = effect(trend_term(data))\nplot(data_with_trend)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n sampling_dist = poisson(2)\nfrequency = UnitVaryingParameter(sampling_dist=sampling_dist)\n\np = Periodic(frequency=frequency)\nplot(p(data))\n\n<Axes: xlabel='Time', ylabel='Observed'>\n"},{"location":"examples/basic/#data-synthesis","title":"Data Synthesis","text":"In this notebook we'll demonstrate how causal-validation can be used to simulate synthetic datasets. We'll start with very simple data to which a static treatment effect may be applied. From there, we'll build up to complex datasets. Along the way, we'll show how reproducibility can be ensured, plots can be generated, and unit-level parameters may be specified.
We observe that we have 10 control units, each of which were sampled from a Gaussian distribution with mean 20 and scale 0.2. Had we wished for our underlying observations to have more or less noise, or to have a different global mean, then we can simply specify that through the config file.
"},{"location":"examples/basic/#reproducibility","title":"Reproducibility","text":"In the above four panels, we can see that whilst the mean and scale of the underlying data generating process is varying, the functional form of the data is the same. This is by design to ensure that data sampling is reproducible. To sample a new dataset, you may either change the underlying seed in the config file.
"},{"location":"examples/basic/#simulating-an-effect","title":"Simulating an effect","text":"In the data we have seen up until now, the treated unit has been drawn from the same data generating process as the control units. However, it can be helpful to also inflate the treated unit to observe how well our model can recover the the true treatment effect. To do this, we simply compose our dataset with an Effect object. In the below, we shall inflate our data by 2%.
The example presented above shows a very simple stationary data generation process. However, we may make our example more complex by including a non-stationary trend to the data.
"},{"location":"examples/basic/#unit-level-parameterisation","title":"Unit-level parameterisation","text":""},{"location":"examples/basic/#conclusions","title":"Conclusions","text":"In this notebook we have shown how one can define their model's true underlying data generating process, starting from simple white-noise samples through to more complex example with periodic and temporal components, perhaps containing unit-level variation. In a follow-up notebook, we show how these datasets may be integrated with Amazon's own AZCausal library to compare the effect estimated by a model with the true effect of the underlying data generating process. A link to this notebook is here.
"},{"location":"examples/placebo_test/","title":"Placebo Testing","text":"from azcausal.core.error import JackKnife\nfrom azcausal.estimators.panel.did import DID\nfrom azcausal.estimators.panel.sdid import SDID\n\nfrom causal_validation import (\n Config,\n simulate,\n)\nfrom causal_validation.effects import StaticEffect\nfrom causal_validation.models import AZCausalWrapper\nfrom causal_validation.plotters import plot\nfrom causal_validation.validation.placebo import PlaceboTest\n\n/home/runner/.local/share/hatch/env/virtual/causal-validation/CYBYs5D-/docs/lib/python3.10/site-packages/pandera/engines/pandas_engine.py:67: UserWarning: Using typeguard < 3. Generic types like List[TYPE], Dict[TYPE, TYPE] will only validate the first element in the collection.\n warnings.warn(\n\n cfg = Config(\n n_control_units=10,\n n_pre_intervention_timepoints=60,\n n_post_intervention_timepoints=30,\n seed=123,\n)\n\nTRUE_EFFECT = 0.05\neffect = StaticEffect(effect=TRUE_EFFECT)\ndata = effect(simulate(cfg))\nplot(data)\n\n<Axes: xlabel='Time', ylabel='Observed'>\n model = AZCausalWrapper(model=SDID(), error_estimator=JackKnife())\nresult = PlaceboTest(model, data).execute()\nresult.summary()\n\n| Model | Effect | Standard Deviation | Standard Error | p-value |\n|-------|--------|--------------------|----------------|---------|\n| SDID | 0.0851 | 0.4079 | 0.129 | 0.5472 |\n\n
did_model = AZCausalWrapper(model=DID())\nPlaceboTest([model, did_model], data).execute().summary()\n\n| Model | Effect | Standard Deviation | Standard Error | p-value |\n|-------|--------|--------------------|----------------|---------|\n| SDID | 0.0851 | 0.4079 | 0.129 | 0.5472 |\n| DID | 0.0002 | 0.2818 | 0.0891 | 0.9982 |\n\n"},{"location":"examples/placebo_test/#placebo-testing","title":"Placebo Testing","text":"
A placebo test is an approach to assess the validity of a causal model by checking if the effect can truly be attributed to the treatment, or to other spurious factors. A placebo test is conducted by iterating through the set of control units and at each iteration, replacing the treated unit by one of the control units and measuring the effect. If the model detects a significant effect, then it suggests potential bias or omitted variables in the analysis, indicating that the causal inference is flawed.
A successful placebo test will show no statistically significant results and we may then conclude that the estimated effect can be attributed to the treatment and not driven by confounding factors. Conversely, a failed placebo test, which shows significant results, suggests that the identified treatment effect may not be reliable. Placebo testing is thus a critical step to ensure the robustness of findings in RCTs. In this notebook, we demonstrate how a placebo test can be conducted in causal-validation.
To demonstrate a placebo test, we must first simulate some data. For the purposes of illustration, we'll simulate a very simple dataset containing 10 control units where each unit has 60 pre-intervention observations, and 30 post-intervention observations.
"},{"location":"examples/placebo_test/#model","title":"Model","text":"We'll now define our model. To do this, we'll use the synthetic difference-in-differences implementation of AZCausal. This implementation, along with any other model from AZCausal, can be neatly wrapped up in our AZCausalWrapper to make fitting and effect estimation simpler.
Now that we have a dataset and model defined, we may conduct our placebo test. With 10 control units, the test will estimate 10 individual effects; 1 per control unit when it is mocked as the treated group. With those 10 effects, the routine will then produce the mean estimated effect, along with the standard deviation across the estimated effect, the effect's standard error, and the p-value that corresponds to the null-hypothesis test that the effect is 0.
In the below, we see that expected estimated effect is small at just 0.08. Accordingly, the p-value attains a value of 0.5, indicating that we have insufficient evidence to reject the null hypothesis and we, therefore, have no evidence to suggest that there is bias within this particular setup.
"},{"location":"examples/placebo_test/#model-comparison","title":"Model Comparison","text":"We can also use the results of a placebo test to compare two or more models. Using causal-validation, this is as simple as supplying a series of models to the placebo test and comparing their outputs. To demonstrate this, we will compare the previously used synthetic difference-in-differences model with regular difference-in-differences.