inst/snippet/lm-simulate2.R

In fastR: Data sets and utilities for Foundations and Applications of Statistics by R Pruim

sim <- function(b0=3,b1=5,sigma=2,
                x=rep(1:5,each=4)    # 4 observations at each of 5 values
    ){
    e <-rnorm(length(x),sd=sigma)
    y <- b0 + b1*x + e
    model <- lm(y~x)  
    ci <- confint(model,2)
    return(b1 > ci[1] && b1 < ci[2])
}
t <- table(replicate(10000,sim())); t / 10000
prop.test(t[2],sum(t),p=0.95)