/*
**	Monte Carlo Program
*/

@   Load Data @

load data[100,5] = exer.txt ;

@ Data generation under Strong Ideal Conditions @

/*
** The regressors are common to individual data sets.
** The errors are different across different data sets
** That is, regressors are nonstochastic but the errors are
*/

seed   = 1;
tt     = 100; @ # of observations @
kk     = 5;  @ # of betas @
iter   = 5000; @ # of sets of different data @
xx  = ones(tt,1)~data[.,2:5]; @ Regressors are fixed @
tb = {0,1,2,3,4} ; @ y = x(2)*1 + x(3)*2 + x(4)*3 + x(5)*4 + e @

storb  = zeros(iter,1);
storse = zeros(iter,1);

i = 1; do while i <= iter;

@ Generating y @
yy  = xx*tb + 2*rndns(tt,1,seed);

@	OLS using yy and xx  @

b  = invpd(xx'xx)*(xx'yy);
e  = yy - xx*b ;

s2 = (e'e)/(tt-kk);
v  = s2*invpd(xx'xx);

se = sqrt(diag(v));

storb[i,1]  = b[2,1];
storse[i,1] = se[2,1];

i = i + 1; endo;

@ Reporting Monte Carlo results @

output file = olsmonte.out reset;

format /rd 12,3;

"Monte Carlo results";
"-----------";
"Mean of OLS b(2)                   ="  meanc(storb);
"s.e. of OLS b(2)                   ="  stdc(storb);
"mean of estimated s.e. of OLS b(2) =" meanc(storse) ;

library pgraph;
graphset;
{a1,a2,a3}=hist(storb,50);
output off ;