siftah

set.seed(1000)
r.steps <- seq(0,1,length.out = 7)
r.sq.change <- matrix(NA, nrow = 500, ncol=7)
empirical = F

for (j in c(1:length(r.steps))) {
  for (i in c(1:500)) {
    Sigma <- matrix(r.steps[j],ncol = 3,nrow = 3)
    Sigma[1,] <- 0.25
    Sigma[,1] <- 0.25
    diag(Sigma) <- 1
    samp <- mvrnorm(n = 100,rep(0,3),Sigma,empirical)

    model <- lm(samp[,1] ~ samp[,2] + samp[,3])
    r <- summary(model)$r.squared

    r.sq.change[i,j] <- r
  }
}

colMeans(r.sq.change)
# [1] 0.14161429 0.12446803 0.10742202 0.10532141 0.09171425 0.08773507 0.07180824