*
* Dougherty, Introduction to Econometrics, 4th ed
* Example from Section 2.4
* Monte Carlo experiment
* (It's probably not expected that you would be able to do this at this point).
*
set x 1 20 = t
compute b1=2,b2=0.5
set y 1 20 = b1+b2*x+%ran(1.0)
*
* Because we are generating our own set of random numbers, we won't get the same
* result as the sample in the book.
*
linreg y
# constant x
*
* 10,000,000 draws takes quite a while. This is set to do 100,000 instead. You can
* experiment with a higher number of you have patience.
*
compute ndraws=100000
*
* This saves both the estimated coefficients and standard errors from the
* simulated regressions.
*
set b2draws 1 ndraws = 0.0
set b2stderrs 1 ndraws = 0.0
do i=1,ndraws
set y 1 20 = b1+b2*x+%ran(1.0)
linreg(noprint) y
# constant x
compute b2draws(i)=%beta(2)
compute b2stderrs(i)=%stderrs(2)
end do i
*
density(grid=automatic,maxgrid=100,smoothing=1.5) b2draws / bx fxn
scatter(style=line,vmin=0.0,$
footer="Figure 2.2 Distribution of b2 in the monte carlo experiment")
# bx fxn
*
* With log normal draws. The following compute the standard deviation of the
* N(0,1) necessary to give a variance 1 log normal draw.
*
compute omega=(1+sqrt(5))/2
compute sigma=sqrt(log(omega))
*
set b2drawsln 1 ndraws = 0.0
do i=1,ndraws
set y 1 20 = b1+b2*x+(exp(%ran(sigma))-sqrt(omega))
linreg(noprint) y
# constant x
compute b2drawsln(i)=%beta(2)
end do i
density(grid=input,maxgrid=100,smoothing=1.5) b2drawsln / bx fxln
*
* With uniform draws. The following computes the upper bound for drawing a uniform
* with unit variance. (The lower bound is - the upper bound).
*
compute ubounds=sqrt(3)
set b2drawsun 1 ndraws = 0.0
do i=1,ndraws
set y 1 20 = b1+b2*x+%uniform(-ubounds,ubounds)
linreg(noprint) y
# constant x
compute b2drawsun(i)=%beta(2)
end do i
*
density(grid=input,maxgrid=100,smoothing=1.5) b2drawsun / bx fxun
scatter(style=line,vmin=0.0,$
footer="Figure 2.3 Distributions of b2 with uniform and lognormal distributions for the disturbance term") 2
# bx fxun
# bx fxln
*
* Examination of standard errors (from page 131)
*
compute theoretical=1.0*sqrt(%xx(2,2))
*
density(grid=automatic,maxgrid=100,smoothing=1.5) b2stderrs / bse fse
scatter(style=line,vmin=0.0,hgrid=theoretical,$
footer="Figure 2.7 Distribution of the standard error of b2")
# bse fse