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. Unique Ids: 3979
No. Observations: 15916
------------------ 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. Unique Ids: 3979
No. Observations: 15916
------------------ 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.42399836 1.41039895 1.39637598 1.42444937 1.40335355 1.42067316
1.42312791 1.40455858 1.42352007]]
Learner ml_g1 RMSE: [[1.40251639 1.43584983 1.39490406 1.41473112 1.42591492 1.38378077
1.45549857 1.41437203 1.4083904 ]]
Classification:
Learner ml_m Log Loss: [[0.69084345 0.69096141 0.6906615 0.67925681 0.67972123 0.67971272
0.66256474 0.66215298 0.66240443]]
------------------ 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.918098 0.064019 14.340978 0.000000 0.792623 1.043574
ATT(2.0,1,3) 1.992117 0.064648 30.814992 0.000000 1.865410 2.118824
ATT(2.0,1,4) 2.954547 0.063175 46.767521 0.000000 2.830726 3.078368
ATT(3.0,1,2) -0.042733 0.065989 -0.647582 0.517255 -0.172070 0.086603
ATT(3.0,2,3) 1.107601 0.065447 16.923713 0.000000 0.979328 1.235874
ATT(3.0,2,4) 2.056854 0.065678 31.317373 0.000000 1.928128 2.185580
ATT(4.0,1,2) 0.003125 0.068599 0.045558 0.963663 -0.131327 0.137577
ATT(4.0,2,3) 0.063969 0.066386 0.963581 0.335256 -0.066146 0.194084
ATT(4.0,3,4) 0.949730 0.067450 14.080530 0.000000 0.817531 1.081929
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.746717 | 1.089479 |
| ATT(2.0,1,3) | 1.819054 | 2.165180 |
| ATT(2.0,1,4) | 2.785426 | 3.123669 |
| ATT(3.0,1,2) | -0.219388 | 0.133921 |
| ATT(3.0,2,3) | 0.932399 | 1.282803 |
| ATT(3.0,2,4) | 1.881033 | 2.232675 |
| ATT(4.0,1,2) | -0.180517 | 0.186767 |
| ATT(4.0,2,3) | -0.113749 | 0.241686 |
| ATT(4.0,3,4) | 0.769165 | 1.130295 |
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.11/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.11/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.486854 0.034271 43.385353 0.0 1.419684 1.554023
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
2.0 1.954921 0.052232 37.427685 0.0 1.852548 2.057294
3.0 1.582227 0.056300 28.103346 0.0 1.471881 1.692574
4.0 0.949730 0.067450 14.080530 0.0 0.817531 1.081929
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
/opt/hostedtoolcache/Python/3.12.11/x64/lib/python3.12/site-packages/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.47996 0.035096 42.168944 0.0 1.411173 1.548746
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
2 0.918098 0.064019 14.340978 0.0 0.792623 1.043574
3 1.551215 0.051339 30.215353 0.0 1.450594 1.651837
4 1.970565 0.046712 42.185370 0.0 1.879011 2.062119
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
/opt/hostedtoolcache/Python/3.12.11/x64/lib/python3.12/site-packages/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.989949 0.038747 51.357273 0.0 1.914006 2.065892
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
-2.0 0.003125 0.068599 0.045558 0.963663 -0.131327 0.137577
-1.0 0.012029 0.040454 0.297359 0.766193 -0.067259 0.091317
0.0 0.990914 0.030721 32.254904 0.000000 0.930701 1.051127
1.0 2.024386 0.045737 44.261574 0.000000 1.934744 2.114029
2.0 2.954547 0.063175 46.767521 0.000000 2.830726 3.078368
------------------ Additional Information ------------------
Score function: observational
Control group: never_treated
Anticipation periods: 0
/opt/hostedtoolcache/Python/3.12.11/x64/lib/python3.12/site-packages/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.989949 0.038747 51.357273 0.0 1.914006 2.065892
------------------ Aggregated Effects ------------------
coef std err t P>|t| 2.5 % 97.5 %
-2.0 0.003125 0.068599 0.045558 0.963663 -0.131327 0.137577
-1.0 0.012029 0.040454 0.297359 0.766193 -0.067259 0.091317
0.0 0.990914 0.030721 32.254904 0.000000 0.930701 1.051127
1.0 2.024386 0.045737 44.261574 0.000000 1.934744 2.114029
2.0 2.954547 0.063175 46.767521 0.000000 2.830726 3.078368
------------------ 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.918098 0.064019 14.340978 0.000000 0.792623 1.043574
ATT(2.0,1,3) 1.992117 0.064648 30.814992 0.000000 1.865410 2.118824
ATT(2.0,1,4) 2.954547 0.063175 46.767521 0.000000 2.830726 3.078368
ATT(3.0,1,2) -0.042733 0.065989 -0.647582 0.517255 -0.172070 0.086603
ATT(3.0,2,3) 1.107601 0.065447 16.923713 0.000000 0.979328 1.235874
ATT(3.0,2,4) 2.056854 0.065678 31.317373 0.000000 1.928128 2.185580
ATT(4.0,1,2) 0.003125 0.068599 0.045558 0.963663 -0.131327 0.137577
ATT(4.0,2,3) 0.063969 0.066386 0.963581 0.335256 -0.066146 0.194084
ATT(4.0,3,4) 0.949730 0.067450 14.080530 0.000000 0.817531 1.081929