## Friday, May 25, 2012

### Estimating Random Coefficients in a panel data

* Imagine that you are concerned that each of your individual panel blocks has a unique random coefficient.

* This might be reasonable if say you are wondering what the effect on employment will be if you increase the federal minimum wage.

* States that already have a minimum wage is higher than the expected increase will probably not find much change while states that are much below the intended increase might find a large change.

* Note that if you define your x variable reasonably well in this case you might not even worry about random coefficients.

* Let us start by setting up the simulation.

clear
set seed 11
set obs 50
* Generate 50 states

gen state=_n

* Now let us speficy a random number of observations per state

gen nobs = rpoisson(15)*4

* Now lets specify a random coefficient per state for the policy change.

gen b=-1+rnormal()

* Now let's expand our sample
expand nobs

gen x=rnormal()*3 + 10

gen y= b*x + rnormal()*3

* Some states have 1 observation while others have 12

* In Econometrics terms this is what we are worrying about:

* Y=Xb + e    e is a random variable
* b = B + u   u is a random variable

* Therefore:
* Y=X(B+u) + e
* Y=XB + Xu + e
* Y=XB + v ;  v=(Xu + e)

xtmixed y x || state: x, nocon

reg y x

* We can see that both xtmixed and normal OLS work pretty well.

* This is primarily because v is uncorrelated with x

* Image now that b is correlated.

gen b2=-1+rnormal()*2+x-10

sum b2
* b2 still has an average effect of -1

gen y2= b2*x + rnormal()*3

reg y2 x
* We can see that OLS is biased positively because:

corr b2 x

xtmixed y2 x || state: x, nocon

* xtmixed does seems to improve the outcome because xtmixed state: x  allows v to be correlation with x.

* However let's test this theory with a Monte Carlo estiamtion:

* First we need to define a program for stata to repeat
cap program drop xtmixed_monte_carlo
program define xtmixed_monte_carlo, rclass

clear
* remove the set see option
set obs 50
* Generate 50 states

gen state=_n

* Now let us speficy a random number of observations per state

gen nobs = rpoisson(15)*4

* Now lets specify a random coefficient per state for the policy change.

gen b=-1+rnormal()

* Now let's expand our sample
expand nobs

gen x=rnormal()*3 + 10

gen y= b*x + rnormal()*3

* Some states have 1 observation while others have 12

* In Econometrics terms this is what we are worrying about:

* Y=Xb + e    e is a random variable
* b = B + u   u is a random variable

* Therefore:
* Y=X(B+u) + e
* Y=XB + Xu + e
* Y=XB + v ;  v=(Xu + e)

xtmixed y x || state: x, nocon
* Now we define a return command that will return the cofficients of interest.
return scalar XTM1 = _b[x]

reg y x
* Now we define a return command that will return the cofficients of interest.
return scalar OLS1 = _b[x]

* We can see that both xtmixed and normal OLS work pretty well.

* This is primarily because v is uncorrelated with x

* Imagine now that b is correlated.

gen b2=-1+rnormal()*2+(x-10)

sum b2
* b2 still has an average effect of -1

gen y2= b2*x + rnormal()*3

reg y2 x
* We can see that OLS is biased positively because:
* Now we define a return command that will return the cofficients of interest.
return scalar OLS2 = _b[x]

corr b2 x

xtmixed y x || x: R.state, nocon
* Now we define a return command that will return the cofficients of interest.
return scalar XTM2 = _b[x]

end

xtmixed_monte_carlo
* We can see the coefficient estimates of our first simulation run.

* This might take a little while.  This command will replate the entire simulation 50 times with the estimates as well.
simulate OLS1=r(OLS1) OLS2=r(OLS2) XTM1=r(XTM1) XTM2=r(XTM2), reps(50): xtmixed_monte_carlo

sum
* We can see that clearly OLS and Xtmixed are both biased because of the correlation between u and x.

* However, they appear to work quite well at identifying the average coefficient when the correlation does not work.

* What xtmixed helps do is to identify the variance of u as well.

simulate OLS1=r(OLS1) OLS2=r(OLS2) XTM1=r(XTM1) XTM2=r(XTM2), reps(500): xtmixed_monte_carlo