par(mfrow=c(1,2))
sigma <- matrix(c(1,.99,.99,1), nrow=2, ncol=2)
out <- MASS::mvrnorm(n=500, mu=c(0,0), Sigma = sigma)

plot(out[,1], out[,2], ylim=c(-3,3), xlim=c(-3,3),
     ylab="Y_i", xlab="mu_i",
     main="High Validity")
abline(h=0, lty=2)
abline(v=0, lty=2)

sigma <- matrix(c(.6,.3,.3,.6), nrow=2, ncol=2)
out <- MASS::mvrnorm(n=500, mu=c(0,0), Sigma = sigma)
plot(out[,1], out[,2], ylim=c(-3,3), xlim=c(-3,3),
     ylab="Y_i", xlab="mu_i",
     main="Low Validity")
abline(h=0, lty=2)
abline(v=0, lty=2)