2021年3月19日星期五

R - Vectorize nested for loops to assign new values to a matrix

I'm currently trying to vectorize this nested for loop to save time during execution but it doesnt seem to work. What I want is to go through every cell of my matrix and check if the value is 0 or 1 then change the value based on a condition. This is the algorithm for the Forest-fire model

for (i in 1:nrow(X)) {    for (j in 1:ncol(X)) {            if (X[i, j] == 2) {        if (runif(1) > (1 - a)^neighbours(X, i, j)) {          B[i, j] <- 1        }      }       else if (X[i, j] == 1) {        burning <- TRUE        if (runif(1) < b) {          B[i, j] <- 0        }      }          }  }  

Here is the neighbours function :

neighbours <- function(A, i, j) {    # calculate number of neighbours of A[i,j] that are infected    # we have to check for the edge of the grid    nbrs <- 0    # sum across row i - 1    if (i > 1) {      if (j > 1) nbrs <- nbrs + (A[i-1, j-1] == 1)      nbrs <- nbrs + (A[i-1, j] == 1)      if (j < ncol(A)) nbrs <- nbrs + (A[i-1, j+1] == 1)    }    # sum across row i    if (j > 1) nbrs <- nbrs + (A[i, j-1] == 1)    nbrs <- nbrs + (A[i, j] == 1)    if (j < ncol(A)) nbrs <- nbrs + (A[i, j+1] == 1)    # sum across row i + 1    if (i < nrow(A)) {      if (j > 1) nbrs <- nbrs + (A[i+1, j-1] == 1)      nbrs <- nbrs + (A[i+1, j] == 1)      if (j < ncol(A)) nbrs <- nbrs + (A[i+1, j+1] == 1)    }    return(nbrs)  }  

And some code to make it work :

set.seed(3)  X <- matrix(2, 21, 21)  X[11, 11:13] <- 1  burning <- FALSE  a= 0.2  b = 0.4  B <- X  

I've started by trying with sapply but couldn't get results back into the matrix and for the past hour I've been trying to use nested foreach loops

library(foreach)  B <-  foreach(i=1:nrow(X), .combine='cbind') %:%    foreach(j=1:ncol(X), .combine='c') %do% {      if (X[i, j] == 2) {        if (runif(1) > (1 - a)^neighbours(X, i, j)) {          1        }      }       else if (X[i, j] == 1) {        burning <- TRUE        if (runif(1) < b) {          0          print(i)          print(j)        }      }    }  

But I'm only getting back the lines I need to change I'm not familiar with vectorization so perhaps I'm missing some basic steps !

https://stackoverflow.com/questions/66714556/r-vectorize-nested-for-loops-to-assign-new-values-to-a-matrix March 20, 2021 at 03:33AM

没有评论:

发表评论