Note
-
Download Jupyter notebook:
https://docs.doubleml.org/stable/examples/did/py_panel_simple.ipynb.
Python: Panel Data Introduction#
In this example, we replicate the results from the guide Getting Started with the did Package of the did-R-package.
As the did-R-package the implementation of DoubleML is based on Callaway and Sant’Anna(2021).
The notebook requires the following packages:
[1]:
import pandas as pd
import numpy as np
from sklearn.linear_model import LinearRegression, LogisticRegression
from doubleml.data import DoubleMLPanelData
from doubleml.did import DoubleMLDIDMulti
Data#
The data we will use is simulated and part of the CSDID-Python-Package.
A description of the data generating process can be found at the CSDID-documentation.
[2]:
dta = pd.read_csv("https://raw.githubusercontent.com/d2cml-ai/csdid/main/data/sim_data.csv")
dta.head()
[2]:
| G | X | id | cluster | period | Y | treat | |
|---|---|---|---|---|---|---|---|
| 0 | 3 | -0.876233 | 1 | 5 | 1 | 5.562556 | 1 |
| 1 | 3 | -0.876233 | 1 | 5 | 2 | 4.349213 | 1 |
| 2 | 3 | -0.876233 | 1 | 5 | 3 | 7.134037 | 1 |
| 3 | 3 | -0.876233 | 1 | 5 | 4 | 6.243056 | 1 |
| 4 | 2 | -0.873848 | 2 | 36 | 1 | -3.659387 | 1 |
To work with the DoubleML-package, we initialize a DoubleMLPanelData object.
Therefore, we set the never-treated units in group column G to np.inf (we have to change the datatype to float).
[3]:
# set dtype for G to float
dta["G"] = dta["G"].astype(float)
dta.loc[dta["G"] == 0, "G"] = np.inf
dta.head()
[3]:
| G | X | id | cluster | period | Y | treat | |
|---|---|---|---|---|---|---|---|
| 0 | 3.0 | -0.876233 | 1 | 5 | 1 | 5.562556 | 1 |
| 1 | 3.0 | -0.876233 | 1 | 5 | 2 | 4.349213 | 1 |
| 2 | 3.0 | -0.876233 | 1 | 5 | 3 | 7.134037 | 1 |
| 3 | 3.0 | -0.876233 | 1 | 5 | 4 | 6.243056 | 1 |
| 4 | 2.0 | -0.873848 | 2 | 36 | 1 | -3.659387 | 1 |
Now, we can initialize the DoubleMLPanelData object, specifying
y_col: the outcomed_cols: the group variable indicating the first treated period for each unitid_col: the unique identification column for each unitt_col: the time columnx_cols: the additional pre-treatment controls
[4]:
dml_data = DoubleMLPanelData(
data=dta,
y_col="Y",
d_cols="G",
id_col="id",
t_col="period",
x_cols=["X"]
)
print(dml_data)
================== DoubleMLPanelData Object ==================
------------------ Data summary ------------------
Outcome variable: Y
Treatment variable(s): ['G']
Covariates: ['X']
Instrument variable(s): None
Time variable: period
Id variable: id
No. Observations: 3979
------------------ DataFrame info ------------------
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 15916 entries, 0 to 15915
Columns: 7 entries, G to treat
dtypes: float64(3), int64(4)
memory usage: 870.5 KB
ATT Estimation#
The DoubleML-package implements estimation of group-time average treatment effect via the DoubleMLDIDMulti class (see model documentation).
The class basically behaves like other DoubleML classes and requires the specification of two learners (for more details on the regression elements, see score documentation). The model will be estimated using the fit() method.
[5]:
dml_obj = DoubleMLDIDMulti(
obj_dml_data=dml_data,
ml_g=LinearRegression(),
ml_m=LogisticRegression(),
control_group="never_treated",
)
dml_obj.fit()
print(dml_obj)
================== DoubleMLDIDMulti Object ==================
------------------ Data summary ------------------
Outcome variable: Y
Treatment variable(s): ['G']
Covariates: ['X']
Instrument variable(s): None
Time variable: period
Id variable: id
No. Observations: 3979
------------------ Score & algorithm ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
------------------ Machine learner ------------------
Learner ml_g: LinearRegression()
Learner ml_m: LogisticRegression()
Out-of-sample Performance:
Regression:
Learner ml_g0 RMSE: [[1.42394271 1.40890149 1.3969726 1.42604003 1.40773671 1.4208955
1.42472087 1.40404173 1.42393828]]
Learner ml_g1 RMSE: [[1.40221074 1.43404007 1.40090569 1.4127083 1.42644506 1.38411456
1.45765855 1.41725837 1.40807185]]
Classification:
Learner ml_m Log Loss: [[0.69072347 0.69065051 0.6907234 0.67936955 0.67947252 0.67945108
0.66203388 0.66239911 0.66228564]]
------------------ Resampling ------------------
No. folds: 5
No. repeated sample splits: 1
------------------ Fit summary ------------------
coef std err t P>|t| 2.5 % 97.5 %
ATT(2.0,1,2) 0.919500 0.064018 14.363176 0.000000 0.794027 1.044972
ATT(2.0,1,3) 1.986508 0.064573 30.763943 0.000000 1.859948 2.113068
ATT(2.0,1,4) 2.956141 0.063195 46.778017 0.000000 2.832281 3.080001
ATT(3.0,1,2) -0.040777 0.065960 -0.618205 0.536440 -0.170056 0.088503
ATT(3.0,2,3) 1.103168 0.065626 16.810028 0.000000 0.974544 1.231792
ATT(3.0,2,4) 2.058363 0.065425 31.461331 0.000000 1.930132 2.186594
ATT(4.0,1,2) 0.002906 0.068334 0.042531 0.966075 -0.131025 0.136838
ATT(4.0,2,3) 0.062507 0.066358 0.941974 0.346206 -0.067552 0.192567
ATT(4.0,3,4) 0.951614 0.067411 14.116704 0.000000 0.819492 1.083736
The summary displays estimates of the \(ATT(g,t_\text{eval})\) effects for different combinations of \((g,t_\text{eval})\) via \(\widehat{ATT}(\mathrm{g},t_\text{pre},t_\text{eval})\), where
\(\mathrm{g}\) specifies the group
\(t_\text{pre}\) specifies the corresponding pre-treatment period
\(t_\text{eval}\) specifies the evaluation period
This corresponds to the estimates given in att_gt function in the did-R-package, where the standard choice is \(t_\text{pre} = \min(\mathrm{g}, t_\text{eval}) - 1\) (without anticipation).
Remark that this includes pre-tests effects if \(\mathrm{g} > t_{eval}\), e.g. \(ATT(4,2)\).
As usual for the DoubleML-package, you can obtain joint confidence intervals via bootstrap.
[6]:
level = 0.95
ci = dml_obj.confint(level=level)
dml_obj.bootstrap(n_rep_boot=5000)
ci_joint = dml_obj.confint(level=level, joint=True)
ci_joint
[6]:
| 2.5 % | 97.5 % | |
|---|---|---|
| ATT(2.0,1,2) | 0.745522 | 1.093477 |
| ATT(2.0,1,3) | 1.811023 | 2.161994 |
| ATT(2.0,1,4) | 2.784399 | 3.127883 |
| ATT(3.0,1,2) | -0.220033 | 0.138479 |
| ATT(3.0,2,3) | 0.924821 | 1.281515 |
| ATT(3.0,2,4) | 1.880560 | 2.236165 |
| ATT(4.0,1,2) | -0.182801 | 0.188613 |
| ATT(4.0,2,3) | -0.117830 | 0.242845 |
| ATT(4.0,3,4) | 0.768416 | 1.134812 |
A visualization of the effects can be obtained via the plot_effects() method.
Remark that the plot used joint confidence intervals per default.
[7]:
fig, ax = dml_obj.plot_effects()
/opt/hostedtoolcache/Python/3.12.10/x64/lib/python3.12/site-packages/matplotlib/cbook.py:1719: FutureWarning: Calling float on a single element Series is deprecated and will raise a TypeError in the future. Use float(ser.iloc[0]) instead
return math.isfinite(val)
/opt/hostedtoolcache/Python/3.12.10/x64/lib/python3.12/site-packages/matplotlib/cbook.py:1719: FutureWarning: Calling float on a single element Series is deprecated and will raise a TypeError in the future. Use float(ser.iloc[0]) instead
return math.isfinite(val)
Effect Aggregation#
As the did-R-package, the \(ATT\)’s can be aggregated to summarize multiple effects. For details on different aggregations and details on their interpretations see Callaway and Sant’Anna(2021).
The aggregations are implemented via the aggregate() method.
Group Aggregation#
To obtain group-specific effects it is possible to aggregate several \(\widehat{ATT}(\mathrm{g},t_\text{pre},t_\text{eval})\) values based on the group \(\mathrm{g}\) by setting the aggregation="group" argument.
[8]:
aggregated = dml_obj.aggregate(aggregation="group")
print(aggregated)
_ = aggregated.plot_effects()
================== DoubleMLDIDAggregation Object ==================
Group Aggregation
------------------ Overall Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
1.486739 0.03422 43.446375 0.0 1.419669 1.553809
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
2.0 1.954050 0.052197 37.436204 0.0 1.851746 2.056353
3.0 1.580765 0.056337 28.059076 0.0 1.470347 1.691184
4.0 0.951614 0.067411 14.116704 0.0 0.819492 1.083736
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
/home/runner/work/doubleml-docs/doubleml-docs/doubleml-for-py/doubleml/did/did_aggregation.py:368: UserWarning: Joint confidence intervals require bootstrapping which hasn't been performed yet. Automatically applying '.aggregated_frameworks.bootstrap(method="normal", n_rep_boot=500)' with default values. For different bootstrap settings, call bootstrap() explicitly before plotting.
warnings.warn(
The output is a DoubleMLDIDAggregation object which includes an overall aggregation summary based on group size.
Time Aggregation#
This aggregates \(\widehat{ATT}(\mathrm{g},t_\text{pre},t_\text{eval})\), based on \(t_\text{eval}\), but weighted with respect to group size. Corresponds to Calendar Time Effects from the did-R-package.
For calendar time effects set aggregation="time".
[9]:
aggregated_time = dml_obj.aggregate("time")
print(aggregated_time)
fig, ax = aggregated_time.plot_effects()
================== DoubleMLDIDAggregation Object ==================
Time Aggregation
------------------ Overall Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
1.479308 0.03508 42.169142 0.0 1.410552 1.548064
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
2 0.919500 0.064018 14.363176 0.0 0.794027 1.044972
3 1.546193 0.051421 30.069143 0.0 1.445409 1.646977
4 1.972231 0.046577 42.343103 0.0 1.880941 2.063521
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
/home/runner/work/doubleml-docs/doubleml-docs/doubleml-for-py/doubleml/did/did_aggregation.py:368: UserWarning: Joint confidence intervals require bootstrapping which hasn't been performed yet. Automatically applying '.aggregated_frameworks.bootstrap(method="normal", n_rep_boot=500)' with default values. For different bootstrap settings, call bootstrap() explicitly before plotting.
warnings.warn(
Event Study Aggregation#
Finally, aggregation="eventstudy" aggregates \(\widehat{ATT}(\mathrm{g},t_\text{pre},t_\text{eval})\) based on exposure time \(e = t_\text{eval} - \mathrm{g}\) (respecting group size).
[10]:
aggregated_eventstudy = dml_obj.aggregate("eventstudy")
print(aggregated_eventstudy)
fig, ax = aggregated_eventstudy.plot_effects()
================== DoubleMLDIDAggregation Object ==================
Event Study Aggregation
------------------ Overall Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
1.989681 0.038712 51.397039 0.0 1.913806 2.065555
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
-2.0 0.002906 0.068334 0.042531 0.966075 -0.131025 0.136838
-1.0 0.012232 0.040396 0.302799 0.762043 -0.066942 0.091406
0.0 0.990576 0.030739 32.225762 0.000000 0.930329 1.050822
1.0 2.022325 0.045645 44.305123 0.000000 1.932862 2.111789
2.0 2.956141 0.063195 46.778017 0.000000 2.832281 3.080001
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
/home/runner/work/doubleml-docs/doubleml-docs/doubleml-for-py/doubleml/did/did_aggregation.py:368: UserWarning: Joint confidence intervals require bootstrapping which hasn't been performed yet. Automatically applying '.aggregated_frameworks.bootstrap(method="normal", n_rep_boot=500)' with default values. For different bootstrap settings, call bootstrap() explicitly before plotting.
warnings.warn(
Aggregation Details#
The DoubleMLDIDAggregation objects include several DoubleMLFrameworks which support methods like bootstrap() or confint(). Further, the weights can be accessed via the properties
overall_aggregation_weights: weights for the overall aggregationaggregation_weights: weights for the aggregation
To clarify, e.g. for the eventstudy aggregation
[11]:
print(aggregated_eventstudy)
================== DoubleMLDIDAggregation Object ==================
Event Study Aggregation
------------------ Overall Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
1.989681 0.038712 51.397039 0.0 1.913806 2.065555
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
-2.0 0.002906 0.068334 0.042531 0.966075 -0.131025 0.136838
-1.0 0.012232 0.040396 0.302799 0.762043 -0.066942 0.091406
0.0 0.990576 0.030739 32.225762 0.000000 0.930329 1.050822
1.0 2.022325 0.045645 44.305123 0.000000 1.932862 2.111789
2.0 2.956141 0.063195 46.778017 0.000000 2.832281 3.080001
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
Here, the overall effect aggregation aggregates each effect with positive exposure
[12]:
print(aggregated_eventstudy.overall_aggregation_weights)
[0. 0. 0.33333333 0.33333333 0.33333333]
If one would like to consider how the aggregated effect with \(e=0\) is computed, one would have to look at the third set of weights within the aggregation_weights property
[13]:
aggregated_eventstudy.aggregation_weights[2]
[13]:
array([0.32875335, 0. , 0. , 0. , 0.32674263,
0. , 0. , 0. , 0.34450402])
Taking a look at the original dml_obj, one can see that this combines the following estimates:
\(\widehat{ATT}(2,1,2)\)
\(\widehat{ATT}(3,2,3)\)
\(\widehat{ATT}(4,3,4)\)
[14]:
print(dml_obj.summary)
coef std err t P>|t| 2.5 % 97.5 %
ATT(2.0,1,2) 0.919500 0.064018 14.363176 0.000000 0.794027 1.044972
ATT(2.0,1,3) 1.986508 0.064573 30.763943 0.000000 1.859948 2.113068
ATT(2.0,1,4) 2.956141 0.063195 46.778017 0.000000 2.832281 3.080001
ATT(3.0,1,2) -0.040777 0.065960 -0.618205 0.536440 -0.170056 0.088503
ATT(3.0,2,3) 1.103168 0.065626 16.810028 0.000000 0.974544 1.231792
ATT(3.0,2,4) 2.058363 0.065425 31.461331 0.000000 1.930132 2.186594
ATT(4.0,1,2) 0.002906 0.068334 0.042531 0.966075 -0.131025 0.136838
ATT(4.0,2,3) 0.062507 0.066358 0.941974 0.346206 -0.067552 0.192567
ATT(4.0,3,4) 0.951614 0.067411 14.116704 0.000000 0.819492 1.083736