*
*  PoorT.PRG
*
*

seed 878765
comp ndraws=10000
comp useobs=300
comp endobs=1+useobs
all endobs

comp sigma = 3.0
comp alpha = 0.0
comp beta = 1.0
comp gamma = 0.5


set x1 = %ran(1.0)
set x2 = %ran(1.0)


frml(variance=sigma**2) vdef v = 0.0
frml(identity) ydef y = alpha + beta*x1 +gamma*x2 + v
group emodel vdef>>v ydef>>y

set v 1 1 = 0.0
set y 1 1 = 0.0
set tsimu 1 ndraws = 0.0

decl vect h(3)

do draws=1,ndraws
simulate(model=emodel) * endobs-1 2


linreg(noprint,robusterrors) y 2 endobs resids
# constant x1 x2  

comp h1 = 0.0
comp h2 = 1.0/%beta(3)
comp h3 = -%beta(2)/(%beta(3)**2)
comp h(1) = h1
comp h(2) = h2
comp h(3) = h3
comp cvm = %xx
comp Sn = tr(h)*cvm*h
comp Sn1 = Sn(1,1)

comp theta = (%beta(2)/%beta(3))-(beta/gamma)
comp Sn2 = sqrt(Sn1)
comp tsimu(draws) = theta/Sn2

end do draws


stat tsimu
density(type=gaussian) tsimu 1 ndraws sgrid sdensity

set zv 1 10000 = %ran(1.0)
density(type=gaussian) zv 1 10000 zgrid zdensity

scatter(style=line,overlay=line,ovsame,patterns) 2
# sgrid sdensity
# zgrid zdensity







set tsize 1 ndraws = 0.0
do jj=1,ndraws
comp tsize(jj) = %if(abs(tsimu(jj))>1.96,1.0,0.0) 
end do draws

stat tsize

/*

comp c01=fix(0.01*ndraws)
comp c05=fix(0.05*ndraws)
comp c10=fix(0.10*ndraws)
comp c50=fix(0.50*ndraws)
comp c90=fix(0.90*ndraws)
comp c95=fix(0.95*ndraws)
comp c975=fix(0.975*ndraws)
comp c99=fix(0.99*ndraws)

comp sorttsimu=%sort(tsimu)

comp t01=sorttsimu(c01,1) 
comp t05=sorttsimu(c05,1) 
comp t10=sorttsimu(c10,1) 
comp t50=sorttsimu(c50,1) 
comp t90=sorttsimu(c90,1) 
comp t95=sorttsimu(c95,1) 
comp t99=sorttsimu(c99,1) 


disp "t-stat"
disp $
"      1%            5%            10%           50%           90%          95%           99%"
disp t01 t05 t10 t50 t90 t95 t99


*/