## Thursday, January 17, 2013

### model selection tests

do file

* It is typically simple to test for superiority of one model over another model when using nested models.

* This can be accomplished in linear regression through use of a Wald test.

* For example:

* Let's say you have one hypothesis that H0: y = x1*B1 and another hypothesis that H1: y = x1*B1 + x2*B2 + x3*B3.

* We can see that H0 is nested within H1.  That is H0 = H1 when B2=B3=0

* Let's generate some sample data:

clear
set obs 10000

gen x1 = rnormal()
gen x2 = rnormal()
gen x3 = rnormal()

gen u = rnormal()*2

gen y1 = x1*2 + u

* Under this case the null is true.

reg y1 x?

test x2=x3=0
* This test should reject the null about %5 of the time when the null is at an alpha = .05

gen y2 = x1*2 - x2*1 + x3*1 + u

reg y2 x?
* In Stata this test is very easy to perform.

test x2=x3=0
* This test should almost always reject the null since the effects of x2 and x3 on y is large.

* For example.

gen y_p1 = normal(x1*.5)
gen y_p2 = normal(x1*.5+x2*1-x3*.2)

gen Y1 = rbinomial(1,y_p1)

probit Y1 x1 x2 x3

test x2=x3=0

* When looking at nested maximum likelihood estimators the likelihood ratio test can be used to test if the difference in explanatory powers of the models is statistically significant.

* The estimator is 2(lk0-lkA) distributed chi-squared with degrees of freedom equal to the difference in degrees of freedom of the two models.

* Let's see this in action:

probit Y1 x1 x2 x3
* Now let's save the degrees of freedom and log likelihood.

global df0 = e(df_m)
global ll0 = e(ll)

* Now for the other model

probit Y1 x1
* Now let's save the degrees of freedom and log likelihood.

global df1 = e(df_m)
global ll1 = e(ll)

* Now let's get the chi-squared statistic.

di "Chi-squared =" 2*(${ll0}-${ll1}) " with DF = " ${df0}-${df1}

di "P value = " 1-chi2(${df0}-${df1},2*(${ll0}-${ll1})) " probability"

* Notice how very similar this LR test is to the previous Wald test in this case.

* Alternatively when the models are actually different:
gen Y2 = rbinomial(1,y_p2)

probit Y2 x1 x2 x3

test x2=x3=0

probit Y2 x1 x2 x3

global df0 = e(df_m)
global ll0 = e(ll)

probit Y2 x1

global df1 = e(df_m)
global ll1 = e(ll)

di "Chi-squared =" 2*(${ll0}-${ll1}) " with DF = " ${df0}-${df1}
di "P value = " 1-chi2(${df0}-${df1},2*(${ll0}-${ll1})) " probability"

* Notice that the Chi-squared value for the Wald test and the LR test is quite different yet the conclusion is identical when it comes to rejecting the null.

* When models are not nested the problem becomes a bit more challenging.

* The Vuong test is a useful test of the goodness of fit of non-nested models.

* Imagine the above data Y2 which is generated based on the assumptions for the probit model.

* We are not sure that the probit is the better model or the logit.

* With the Vuong (1989) test of non-nested models we can construct a test statistic based on the log likelihoods of each individual observation.

* This is done by treating the difference in log likelihoods as a random variable and constructing a chi-squared estimator based on that difference (once again referencing Wooldridge class notes from EC 821B).

* First we will estimate the probit.

probit Y2 x1 x2 x3

* We will predict the fitted values which are the predicted probabilities of a draw of 1.

predict Y2_probit

* The log likelihood can be found from this by taking the log of those probabilities when the draw is 1 and the log of 1 less those probabilities when the draw is 0.

gen ll_probit = Y2*log(Y2_probit) + (1-Y2)*log(1-Y2_probit)

* In order to test that this is actually the case:
sum ll_probit

di "The log likelihood from the MLE should be equal to " r(N)*r(mean)

* Bingo!

* Now we will do the same with the logit.

logit Y2 x1 x2 x3

* We will predict the fitted values which are the predicted probabilities of a draw of 1.

predict Y2_logit

* The log likelihood can be found from this by taking the log of those probabilities when the draw is 1 and the log of 1 less those probabilities when the draw is 0.

gen ll_logit = Y2*log(Y2_logit) + (1-Y2)*log(1-Y2_logit)

* In order to test that this is actually the case:
sum ll_logit

di "The log likelihood from the MLE should be equal to " r(N)*r(mean)

* We next construct a variable called dif, which is a measure of the individual differences of the likelihoods.

gen dif = ll_probit-ll_logit

* Now this variable should only be statistically different from 0 about 5% of the time (at the 5% level) if the two models have equal explanatory power.

reg dif

* The constant dif only seems to be statistically significant a small proportion of the times indicating that as discussed previous the difference between a probit and a logit model is extremely small.

* It is possible to compare many other non-nested models in this way.

* Binary response models are particularly convenient examples because the log-likelihood statistic is so easy to construct.

* However, log likelihoods are easy to recover and or necessary to construct for many maximum likelihood procedures.