size <- 6 mat <- matrix(seq_len(size ^ 2), ncol = size) low <- 0 high <- 3 delta <- rep(seq_len(ncol(mat)), nrow(mat)) - rep(seq_len(nrow(mat)), each = ncol(mat)) mat[delta < low | delta > high] <- NA mat 

this works with 5000 x 5000 matrices on my machine

size <- 6 mat <- matrix(seq_len(size ^ 2), ncol = size) low <- 0 high <- 3 delta <- rep(seq_len(ncol(mat)), nrow(mat)) - rep(seq_len(nrow(mat)), each = ncol(mat)) #or Ben Bolker's better alternative delta <- row(mat) - col(mat) mat[delta < low | delta > high] <- NA mat 

this works with 5000 x 5000 matrices on my machine