Mediation analysis

[kshedden]

This is an implementation of mediation analysis following the approach of Imai and collaborators.

[coveralls]

Coverage increased (+0.05%) to 83.95% when pulling 00a8a43 on kshedden:mediation into a1f53b0 on statsmodels:master.

[coveralls]

Coverage increased (+0.05%) to 83.95% when pulling 00a8a43 on kshedden:mediation into a1f53b0 on statsmodels:master.

[coveralls]

Coverage increased (+0.04%) to 83.95% when pulling 836396e on kshedden:mediation into a1f53b0 on statsmodels:master.

[kshedden]

A notebook illustrating the mediation procedure:

http://nbviewer.ipython.org/urls/umich.box.com/shared/static/jpmd9y99259u6dv0rj6p46993981m7zm.ipynb

[josef-pkt]

just some general comments:

I looked at the reference briefly after you submitted the PR. I understand the basic idea but haven't figured out the assumptions yet. Overall I'm far behind in my readings about the differences of "causality" between statistics and econometrics, mainly Imbens and Rubin articles (and a new book).

Essentially, the direction that we have been and are going is that additional functionality in statistics are your area in statsmodels and I try to keep up enough to help with the integration into statsmodels.
(It took me more than half a year to understand the basic idea of GEE and it's relationship to GMM, and several years to understand most of GLM/LEF/QMLE)

One question here is where to put it. Should we create a new subdirectory causality or treatment or with a similar category? Or can we add it to stats or another existing subpackage.

Similar question for propensity score matching, regression discontinuity design, and whatever else will come into statsmodels in future.
I have no idea whether we will get SEM and graphical models anytime soon, but I doubt it.

[kshedden]

Thanks for the comments. I don't have a strong opinion about where to put it, your suggestions are all fine with me.

I've been using this in a project and it holds up pretty well. I used to just use the simple method of multiplying coefficients but was never sure if that was meaningful outside of linear models. This approach seems more correct for generic use.

[coveralls]

Coverage increased (+0.04%) to 83.95% when pulling 511f602 on kshedden:mediation into a1f53b0 on statsmodels:master.

[sjgiorgi]

Line 351 of mediation.py: index = ["ACME (control)", "ACME (treated)", "ADE (control)", "ADE (treated)",
Line 356 of mediation.py: for i, vec in enumerate([self.ACME_ctrl, self.ACME_tx, self.ADE_tx, self.ADE_ctrl,

Seems like self.ADE_tx and self.ADE_ctrl should be switched on line 356?

[kshedden]

@sjgiorgi yes, you appear to be correct. I have fixed it. Thanks. Any other comments are very welcome.

[josef-pkt]

since these are user facing attributes, longer more descriptive names would be usefull

[josef-pkt]

Looks good based on skimming the code

incomplete docstrings, needs eventually some standardized naming with treatment effects, and can go into new folder treatment (if we stick with that name)
related PR #2455 which doesn't allow different types of outcome or treatment models yet, and might be able to follow the same or similar pattern as here.

I did a bit of background reading and opened #2627 for a possible second approach.

[josef-pkt]

I'm merging this, naming convention might and location will still change before 0.8 release

[josef-pkt]

Thanks @kshedden

As I mentioned, I started to read to catch up with the topic but I'm still far behind.