Simulating fake data:
> x <- 10 * runif(10)
> y <- rnorm(10, mean = 2 + 3 * x, 3)Visualizing:
> plot(y ~ x, col = "blue", pch = 19)
Estimating a linear model:
> model <- lm(y ~ x)Calculating and plotting 95 % confidence interval based on simulations:
> xr <- 100
> nb <- 1000
> ci <- .95
> eps <- (1 - ci) / 2
> xs <- seq(min(x), max(x), length = xr)
> coef_val <- MASS::mvrnorm(nb, coef(model), vcov(model))
> ys <- t(coef_val %*% rbind(1, xs))
> predconf <- t(apply(ys, 1, quantile, c(eps, 1 - eps)))Let’s plot all this:
> plot(y ~ x, type = "n")
> matlines(xs, ys, col = adjustcolor("black", .05), lty = 1)
> matlines(xs, predconf, col = "red", lty = 1, lwd = 2)
> points(y ~ x, col = "blue", pch = 19)