# R  CMD BATCH mc_asydist.R mc_asydist.out

postscript(file="xyplot_asy2.eps",append=F,print.it=F)

make.data <- function( b, N, sseed=NULL) {
  if (!is.null(sseed)) set.seed(sseed)
  x  <- rnorm(N)

# assume errors are uniformly distributed from -4 to 4

  e  <- runif(N,-4,4) 
  xb <- b[1] + b[2]* x
  y <- xb + e
  dta <- as.data.frame(cbind( y=y, x=x, e=e))
  return(dta)
}

## set some parameters

b <- c(1,3)

sim <- 1000

## number of cases in simulated data; choose values to put in here

Nvec <- c(n1, n2, n3, n4, n5, n6)

## set mfrow paramater to produce multiple plots on same page for 
## purposes of comparison

par(mfrow=c(3,2))

for(i in 1:nrow(as.matrix(Nvec))){

# define vector to hold OLS estimates

 bsim <- rep(0,sim)

 for(j in 1:sim){

  dta <- make.data( b, Nvec[i])
    Y <- dta$y
    X <- dta$x

  xyreg<-lm(Y~X)

# save the coefficient values to bhat

  bhat <- coef(xyreg)

  bsim[j] <- bhat[2]

 }

# plot a histogram

  hist(bsim,xlim=c(0,5))

}