Applying refutation tests to the Lalonde and IHDP datasets

Import the Dependencies

[3]:

import dowhy
from dowhy import CausalModel
import pandas as pd
import numpy as np

Loading the Dataset

Infant Health and Development Program Dataset (IHDP)

The measurements used are on the child—birth weight, head circumference, weeks bornpreterm, birth order, first born, neonatal health index (see Scott and Bauer 1989), sex, twinstatus—as well as behaviors engaged in during the pregnancy—smoked cigarettes, drankalcohol, took drugs—and measurements on the mother at the time she gave birth—age,marital status, educational attainment (did not graduate from high school, graduated fromhigh school, attended some college but did not graduate, graduated from college), whethershe worked during pregnancy, whether she received prenatal care—and the site (8 total) inwhich the family resided at the start of the intervention. There are 6 continuous covariatesand 19 binary covariates.

Reference

Hill, J. L. (2011). Bayesian nonparametric modeling for causal inference. Journal of Computational and Graphical Statistics, 20(1), 217-240. https://doi.org/10.1198/jcgs.2010.08162

[4]:

data = pd.read_csv("https://raw.githubusercontent.com/AMLab-Amsterdam/CEVAE/master/datasets/IHDP/csv/ihdp_npci_1.csv", header = None)
col =  ["treatment", "y_factual", "y_cfactual", "mu0", "mu1" ,]
for i in range(1,26):
    col.append("x"+str(i))
data.columns = col
data = data.astype({"treatment":'bool'}, copy=False)
data.head()

[4]:
treatment y_factual y_cfactual mu0 mu1 x1 x2 x3 x4 x5 ... x16 x17 x18 x19 x20 x21 x22 x23 x24 x25
0 True 5.599916 4.318780 3.268256 6.854457 -0.528603 -0.343455 1.128554 0.161703 -0.316603 ... 1 1 1 1 0 0 0 0 0 0
1 False 6.875856 7.856495 6.636059 7.562718 -1.736945 -1.802002 0.383828 2.244320 -0.629189 ... 1 1 1 1 0 0 0 0 0 0
2 False 2.996273 6.633952 1.570536 6.121617 -0.807451 -0.202946 -0.360898 -0.879606 0.808706 ... 1 0 1 1 0 0 0 0 0 0
3 False 1.366206 5.697239 1.244738 5.889125 0.390083 0.596582 -1.850350 -0.879606 -0.004017 ... 1 0 1 1 0 0 0 0 0 0
4 False 1.963538 6.202582 1.685048 6.191994 -1.045229 -0.602710 0.011465 0.161703 0.683672 ... 1 1 1 1 0 0 0 0 0 0

5 rows × 30 columns

Lalonde Dataset

A data frame with 445 observations on the following 12 variables.

  • age: age in years.

  • educ: years of schooling.

  • black: indicator variable for blacks.

  • hisp: indicator variable for Hispanics.

  • married: indicator variable for martial status.

  • nodegr: indicator variable for high school diploma.

  • re74: real earnings in 1974.

  • re75: real earnings in 1975.

  • re78: real earnings in 1978.

  • u74: indicator variable for earnings in 1974 being zero.

  • u75: indicator variable for earnings in 1975 being zero.

  • treat: an indicator variable for treatment status.

References

Dehejia, Rajeev and Sadek Wahba. 1999.``Causal Effects in Non-Experimental Studies: Re-Evaluating the Evaluation of Training Programs.’’ Journal of the American Statistical Association 94 (448): 1053-1062.

LaLonde, Robert. 1986. ``Evaluating the Econometric Evaluations of Training Programs.’’ American Economic Review 76:604-620.

[2]:

from rpy2.robjects import r as R

from os.path import expanduser
home = expanduser("~")

%reload_ext rpy2.ipython

# %R install.packages("Matching")
%R library(Matching)

R[write to console]: Loading required package: MASS

R[write to console]: ##
##  Matching (Version 4.9-7, Build Date: 2020-02-05)
##  See http://sekhon.berkeley.edu/matching for additional documentation.
##  Please cite software as:
##   Jasjeet S. Sekhon. 2011. ``Multivariate and Propensity Score Matching
##   Software with Automated Balance Optimization: The Matching package for R.''
##   Journal of Statistical Software, 42(7): 1-52.
##



[2]:

array(['Matching', 'MASS', 'tools', 'stats', 'graphics', 'grDevices',
       'utils', 'datasets', 'methods', 'base'], dtype='<U9')

[7]:

%R data(lalonde)
%R -o lalonde
lalonde = lalonde.astype({'treat':'bool'}, copy=False)
lalonde.head()

[7]:
age educ black hisp married nodegr re74 re75 re78 u74 u75 treat
1 37 11 1 0 1 1 0.0 0.0 9930.05 1 1 True
2 22 9 0 1 0 1 0.0 0.0 3595.89 1 1 True
3 30 12 1 0 0 0 0.0 0.0 24909.50 1 1 True
4 27 11 1 0 0 1 0.0 0.0 7506.15 1 1 True
5 33 8 1 0 0 1 0.0 0.0 289.79 1 1 True

Step 1: Building the model

IHDP

[12]:

# Create a causal model from the data and given common causes

common_causes = []

for i in range(1, 26):
    common_causes += ["x"+str(i)]

ihdp_model = CausalModel(
                data=data,
                treatment='treatment',
                outcome='y_factual',
                common_causes=common_causes
            )
ihdp_model

WARNING:dowhy.causal_model:Causal Graph not provided. DoWhy will construct a graph based on data inputs.
INFO:dowhy.causal_graph:If this is observed data (not from a randomized experiment), there might always be missing confounders. Adding a node named "Unobserved Confounders" to reflect this.
INFO:dowhy.causal_model:Model to find the causal effect of treatment ['treatment'] on outcome ['y_factual']

[12]:

<dowhy.causal_model.CausalModel at 0x7fafe92321d0>

Lalonde

[15]:

lalonde_model = CausalModel(
                data=lalonde,
                treatment='treat',
                outcome='re78',
                common_causes='nodegr+black+hisp+age+educ+married'.split('+')
            )
lalonde_model

WARNING:dowhy.causal_model:Causal Graph not provided. DoWhy will construct a graph based on data inputs.
INFO:dowhy.causal_graph:If this is observed data (not from a randomized experiment), there might always be missing confounders. Adding a node named "Unobserved Confounders" to reflect this.
INFO:dowhy.causal_model:Model to find the causal effect of treatment ['treat'] on outcome ['re78']

[15]:

<dowhy.causal_model.CausalModel at 0x7fafe9170610>

Step 2: Identification

IHDP

[25]:

#Identify the causal effect for the ihdp dataset
ihdp_identified_estimand = ihdp_model.identify_effect()

INFO:dowhy.causal_identifier:Common causes of treatment and outcome:['x25', 'x12', 'U', 'x4', 'x14', 'x5', 'x17', 'x18', 'x16', 'x9', 'x6', 'x1', 'x15', 'x19', 'x20', 'x10', 'x23', 'x21', 'x2', 'x22', 'x13', 'x11', 'x24', 'x8', 'x3', 'x7']
WARNING:dowhy.causal_identifier:If this is observed data (not from a randomized experiment), there might always be missing confounders. Causal effect cannot be identified perfectly.

WARN: Do you want to continue by ignoring any unobserved confounders? (use proceed_when_unidentifiable=True to disable this prompt) [y/n] y

INFO:dowhy.causal_identifier:Instrumental variables for treatment and outcome:[]

Lalonde

[26]:

#Identify the causal effect for the lalonde dataset
lalonde_identified_estimand = lalonde_model.identify_effect()

INFO:dowhy.causal_identifier:Common causes of treatment and outcome:['hisp', 'black', 'U', 'age', 'married', 'educ', 'nodegr']
WARNING:dowhy.causal_identifier:If this is observed data (not from a randomized experiment), there might always be missing confounders. Causal effect cannot be identified perfectly.

WARN: Do you want to continue by ignoring any unobserved confounders? (use proceed_when_unidentifiable=True to disable this prompt) [y/n] y

INFO:dowhy.causal_identifier:Instrumental variables for treatment and outcome:[]

Step 3: Estimation (using propensity score weighting)

IHDP

[32]:

ihdp_estimate = ihdp_model.estimate_effect(
                    ihdp_identified_estimand,
                    method_name="backdoor.propensity_score_weighting"
                )

print("The Causal Estimate is " + str(ihdp_estimate.value))

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: y_factual~treatment+x25+x12+x4+x14+x5+x17+x18+x16+x9+x6+x1+x15+x19+x20+x10+x23+x21+x2+x22+x13+x11+x24+x8+x3+x7

The Causal Estimate is 3.409737824404964

/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Lalonde

[31]:

lalonde_estimate = lalonde_model.estimate_effect(
                        lalonde_identified_estimand,
                        method_name="backdoor.propensity_score_weighting"
                    )

print("The Causal Estimate is " + str(lalonde_estimate.value))

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: re78~treat+hisp+black+age+married+educ+nodegr

The Causal Estimate is 1614.1203743033784

/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Step 4: Refutation

IHDP

Add Random Common Cause

[33]:

ihdp_refute_random_common_cause = ihdp_model.refute_estimate(
                                        ihdp_identified_estimand,
                                        ihdp_estimate,
                                        method_name="random_common_cause"
                                    )

print(ihdp_refute_random_common_cause)

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: y_factual~treatment+x25+x12+x4+x14+x5+x17+x18+x16+x9+x6+x1+x15+x19+x20+x10+x23+x21+x2+x22+x13+x11+x24+x8+x3+x7+w_random

Refute: Add a Random Common Cause
Estimated effect:(3.409737824404964,)
New effect:(3.4082516978396793,)


/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Replace Treatment with Placebo

[35]:

ihdp_refute_placebo_treatment = ihdp_model.refute_estimate(
                                    ihdp_identified_estimand,
                                    ihdp_estimate,
                                    method_name="placebo_treatment_refuter",
                                    placebo_type="permute"
                                )

print(ihdp_refute_placebo_treatment)

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: y_factual~placebo+x25+x12+x4+x14+x5+x17+x18+x16+x9+x6+x1+x15+x19+x20+x10+x23+x21+x2+x22+x13+x11+x24+x8+x3+x7

Refute: Use a Placebo Treatment
Estimated effect:(3.409737824404964,)
New effect:(-0.030843695817419192,)


/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Remove Random Subset of Data

[36]:

ihdp_refute_random_subset = ihdp_model.refute_estimate(
                                    ihdp_identified_estimand,
                                    ihdp_estimate,
                                    method_name="data_subset_refuter",
                                    subset_fraction=0.8
                            )
print(ihdp_refute_random_subset)

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: y_factual~treatment+x25+x12+x4+x14+x5+x17+x18+x16+x9+x6+x1+x15+x19+x20+x10+x23+x21+x2+x22+x13+x11+x24+x8+x3+x7

Refute: Use a subset of data
Estimated effect:(3.409737824404964,)
New effect:(3.490166414737983,)


/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Lalonde

Add Random Common Cause

[38]:

lalonde_refute_random_common_cause = lalonde_model.refute_estimate(
                                            lalonde_identified_estimand,
                                            lalonde_estimate,
                                            method_name="random_common_cause"
                                        )

print(lalonde_refute_random_common_cause)

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: re78~treat+hisp+black+age+married+educ+nodegr+w_random

Refute: Add a Random Common Cause
Estimated effect:(1614.1203743033784,)
New effect:(1623.995451276467,)


/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Replace Treatment with Placebo

[39]:

lalonde_refute_placebo_treatment = lalonde_model.refute_estimate(
                                        lalonde_identified_estimand,
                                        lalonde_estimate,
                                        method_name="placebo_treatment_refuter",
                                        placebo_type="permute"
                                    )

print(lalonde_refute_placebo_treatment)

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: re78~placebo+hisp+black+age+married+educ+nodegr

Refute: Use a Placebo Treatment
Estimated effect:(1614.1203743033784,)
New effect:(913.2564092536604,)


/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)

Remove Random Subset of Data

[40]:

lalonde_refute_random_subset = lalonde_model.refute_estimate(
                                    lalonde_identified_estimand,
                                    lalonde_estimate,
                                    method_name="data_subset_refuter",
                                    subset_fraction=0.9
                                )

print(lalonde_refute_random_subset)

INFO:dowhy.causal_estimator:INFO: Using Propensity Score Weighting Estimator
INFO:dowhy.causal_estimator:b: re78~treat+hisp+black+age+married+educ+nodegr

Refute: Use a subset of data
Estimated effect:(1614.1203743033784,)
New effect:(1425.6467153243575,)


/home/tanmay/model_env/lib/python3.7/site-packages/sklearn/utils/validation.py:760: DataConversionWarning: A column-vector y was passed when a 1d array was expected. Please change the shape of y to (n_samples, ), for example using ravel().
  y = column_or_1d(y, warn=True)