proc iml; /* Set up the simulation */ x = j(10000,1,.); call randgen(x,"Normal",2,2); e = j(10000,1,.); call randgen(e,"Normal",0,1); y = 2 + 0.5*x + e; /* Determine the population regression line */ xx = j(nrow(x),1,1)||x; b = inv(xx`*xx)*xx`*y; print b; yxx = y||xx; /* Look at the bias properties */ /* Specify module that will sample with replacement */ start samplereplace(A, nSamples); results = j(nSamples,1); /** allocate result matrix **/ call randgen(results, "Uniform"); /** fill with random U(0,1) **/ if ncol(A)>nrow(A) then n=ncol(A); else n=nrow(A); results = ceil(n*results); /** convert to integers 1,2,...n **/ ss = A[results[1,],]; do i=2 to nrow(results); ss = ss//A[results[i,],]; end; return(ss); /** reshape and return from A **/ finish; iter = 1000; b_ = j(iter,2,.); /* Perform simulation; use many samples of 200 observations and keep the parameter values */ do i=1 to iter; s = samplereplace(yxx,100); y_ = s[,1]; xx_ = s[,2:3]; b_[i,] = t(inv(xx_`*xx_)*xx_`*y_); end; /* Check whether you have bias */ b_mean = b_[:,]; print b_mean;