*! version 1.5 Nov 2016
help psestimate

Title

psestimate -- Estimate the propensity score proposed by Imbens and Rubin (2015)

Syntax

psestimate depvar [indepvars] [if] [in]
[, options]

totry(indepvars) specify list of covariates to try; default is all

notry(varlist) specify list of covariates to exclude; default is none

nolin prevent algorithm of testing linear terms

noquad prevent algorithm of testing quadratic terms

clinear(real) threshold value for likelihood ratio test of first order covariates; default is 1

cquadratic(real) threshold value for likelihood ratio test of second order covariates; default is 2.71

iterate(#) perform maximum of # iterations in each logit; default is 16000

genpscore(newvar) generate new variable with propensity score estimation

genlor(newvar) generate new variable with log odds ratio

indepvars may contain factor variables; see fvvarlist.

Description

The psestimate command estimates the propensity score proposed by Imbens and Rubin (2015).
In particular, it implements the algorithm outlined by Imbens (2015), which estimates the propensity score for a binary dependent variable indicating treatment status.
The main purpose of the program is to select a linear or quadratic function of covariates to include in the estimation function of the propensity score.

A binary treatment variable must be specified as dependent variable.
A subset of covariates may be explicitly included in the linear part of the specification as independent variables.
This initial configuration corresponds to the base model.
All specifications are fitted with a [R] logit model by maximum likelihood.

The algorithm selects first order terms from all remaining variables of the dataset (i.e. excluding variables of the base model),
unless a subset of variables is specified with the totry(indepvars) option.
Specific variables may be excluded using the notry(varlist) option.

The selection of first order terms is performed in a stepwise fashion, comparing the base (nested) model to a model with one single additional covariate.
A likelihood ratio test (LRT) on the null hypothesis of non-significance of the additional coefficient is performed (see [R] lrtest).
All covariates that have not been included are tested and the algorithm selects the one associated with the highest LRT statistic,
unless no covariate meets the LRT statistic threshold specified in clinear(real).
This covariate is then included in the model.

The process of selecting and including additional linear terms is carried out until the highest LRT statistic is less than clinear(real) or there are no remaining covariates to add.

The algorithm chooses second order terms only from covariates selected for the linear specification, performing analogue tests for selection among all interactions and quadratic terms.
This second process of selecting and including additional quadratic terms is carried out until the highest LRT statistic is less than cquadratic(real) or there are no remaining quadratic terms to add.

Options

totry(varlist) specifies the vector of covariates from which the first (and potentially second) order terms are going to be selected.
The default is to include all variables in the dataset, exluding the depvar and other base model covariates indicated in indepvars (if any).

notry(varlist) specifies a vector of covariates to be excluded from the selection of terms.

nolin prevents the program from testing linear terms, choosing quadratic terms from covariates specified as independent variables.
It can be useful to speed up the algorithm if the linear part is already chosen.
If specified, option clinear is irrelevant.
This option may not be combined with noquad.

noquad prevents the program from testing quadratic terms, ending when all linear terms have been added.
It can be useful to speed up the algorithm if the quadratic part is not desired.
If specified, option cquadratic is irrelevant.
This option may not be combined with nolin.

clinear(real) specifies the threshold value used for the addition of first order (linear) terms.
The decision is based on the likelihood ratio test statistic for the null hypothesis that the coefficient of the additional first order term is equal to zero.
See [R] lrtest for additional information.
If the nolin option is specified, then this option is irrelevant.
Default value is 1.

cquadratic(real) specifies the threshold value used for the addition of second order (quadratic) terms.
The decision is based on the likelihood ratio test statistic for the null hypothesis that the coefficient of the additional second order term is equal to zero.
See [R] lrtest for additional information.
If the noquad option is specified, then this option is irrelevant.
Default value is 2.71.

iterate(#) specifies the maximum number of iterations in each logit estimation.
Stata's default value is 1600.
See [R] logit and [R] maximize for additional information.

genpscore(newvar) specifies that a new variable with the estimated propensity scores is generated, named newvar.

genlor(newvar) specifies that a new variable with the log odds ratio of the estimated propensity score is generated, named newvar.

Remarks on executing time

The algorithm implemented by psestimate may take a (very) long time executing.
A progress indicator is displayed while the program selects first and second order terms, to monitor progress.
The number in parenthesis corresponds to the upper bound of iterations the algorithm could perform before running out of covariates (or its combinations, if applicable) to try.

The selection of linear terms is usually faster than that of quadratic ones.
It is a good idea to start using the command with the noquad option and then, when linear terms are chosen, include them explicitely as indepvars and use nolin to skip the first stage.

Examples

For these examples I use the "Lalonde Experimental Data (Dehejia-Wahba Sample)", corresponding to the data analyzed by Dehejia and Wahba (1999) and available on Dehejia's website.
The dataset contains 445 observations with information on treatment status and various other characteristics.

To install the ancillary files (nswre74.dta and replicate_lalonde.do), remember to use net gate after ssc install:

. ssc install psestimate

. net get psestimate

Setup

. use nswre74

Select PS model for treatment variable

. psestimate treat

Select PS model from restricted list of covariates and lowered quadratic threshold

. psestimate treat, totry(age-nodeg re*) cquad(.8)

Select PS model with income and unemployment dummies as basic covariates

. foreach y in 74 75 78 {

.	gen u`y' = (re`y'==0)

. }

. psestimate treat re* u*

Estimate propensity score with no quadratic terms

. psestimate treat, genpscore(ps) noquad

Estimate log odds ratio with explicit selection of linear terms

. psestimate treat age-nodeg, nolin genlor(logodds)

Stored results

psestimate stores the following in r():

Scalars

r(C_l) threshold value of linear terms

r(C_q) threshold value of quadratic terms

Macros

r(tvar) treatment variable

r(K_b) basic terms (explicitly included)

r(K_l) linear terms

r(K_q) quadratic terms

r(h) full model

Additionally, psestimate stores all results stored in e() by [R] logit after fitting the final model with all selected terms.

Author

Alvaro Carril

Research Analyst at J-PAL LAC

acarril@fen.uchile.cl

Acknowledgements

This program started as a refinement of Juan Ignacio Elorrieta's work, whose code provided a significant headstart.
I'm also indebted to Diego Escobar, who helped me to fine tune several features and to squash some bugs.
All remaining errors are my own.

References

Dehejia, Rajeev H. and Sadek Wahba. 1999.
"Causal Effects in Nonexperimental Studies".
Journal of the American Statistical Association 94(448): 1053-1062.

Imbens, Guido W. and Donald B. Rubin. 2015.
Causal Inference in Statistics, Social, and Biomedical Sciences.
New York: Cambridge University Press.

Imbens, Guido W. 2015.
"Matching Methods in Practice: Three Examples".
Journal of Human Resources 50(2): 373-419.