I'm relatively new to programming in R. I want to simulate the movement of one individual across a 5x5 grid given that the grid varies in its environmental conditions and that movement from one cell to another is based on the environmental conditions of their immediate neighbours. The final result of this simulation I want is the location of the individual after x number of time steps.
First, I made a data frame that contained the x,y coordinates of the grid and their environmental conditions. I then calculated the resistance to movement and its inverse based off my random environmental conditions (v1,v2).
env_cond<-data.frame(x=rep(1:5,5),y=rep(1:5,each=5),v1=rnorm(25),v2=rnorm(25))
env_cond$resistance<- res_surf<- (env_cond [1,3] - env_cond [,3])^2 + (env_cond [1,4]- env_cond [,4])^2
env_cond$inv_res <- 1/env_cond$resistance #where movement is based on inverse resistance
env_cond$cell_num <- 1:25
head (env_cond)
x y v1 v2 resistance inv_res cell_num
1 1 1 1.233266019 0.3554372 0.0000000 Inf 1
2 2 1 0.499331993 0.3780565 0.5391708 1.8546999 2
3 3 1 1.633103368 0.7464020 0.3127234 3.1977142 3
4 4 1 -0.583125893 0.6591043 3.3914933 0.2948554 4
5 5 1 0.929743728 -0.7338991 1.2787793 0.7819958 5
6 1 2 0.009317203 0.2060074 1.5203800 0.6577303 6
>
Next, I created a neighbours matrix. I'm assuming that an individual can move to only its 4 immediate neighbours and nowhere else on the grid. It shows the cell numbers of the grid which correspond to the immediate 4 neighbours of a cell. For example, cell 1 (which corresponds to x=1,y=1) gives NA for North because it cannot move above the grid space.
north <- ifelse (env_cond$y==1, NA, env_cond$cell_num-5) #y+1
south <- ifelse (env_cond$y==5, NA, env_cond$cell_num+5) #y-1
west <- ifelse (env_cond$x==1, NA, env_cond$cell_num-1) #x-1
east <- ifelse (env_cond$x==5, NA, env_cond$cell_num+1) #x+1
neighbours <- data.frame(north, south, west, east)
head (neighbours)
north south west east
1 NA 6 NA 2
2 NA 7 1 3
3 NA 8 2 4
4 NA 9 3 5
5 NA 10 4 NA
6 1 11 NA 7
>
I created a probability matrix by first assigning the inverse resistance values of the neighbours to the cell numbers. I replaced the NA's by 0 to illustrate the impossibility of movement and any values for infinity arbitrarily by 10. I then converted the values to probability:
prob_mat <- cbind (env_cond$inv_res [neighbours$north], env_cond$inv_res [neighbours$south],env_cond$inv_res [neighbours$west], env_cond$inv_res [neighbours$east])
colnames(prob_mat) <- c("y+1", "y-1", "x-1", "x+1") #renamed the columns of prob matrix
#changing NA to O
prob_mat[is.na(prob_mat)]<-0
#changing inf to 10
prob_mat [6, 1] <- 10
prob_mat [2, 3] <- 10
prob_mat1 <- matrix (nrow = nrow(prob_mat), ncol=4)
for (i in 1:nrow (prob_mat)) {
prob_mat1 [i,]<- prob_mat[i,]/sum(prob_mat[i,])
head (prob_mat1)
[,1] [,2] [,3] [,4]
[1,] 0.0000000 0.26179048 0.0000000 0.73820952
[2,] 0.0000000 0.01556767 0.7459112 0.23852109
[3,] 0.0000000 0.06208574 0.8092602 0.12865408
[4,] 0.0000000 0.10119069 0.7221972 0.17661214
[5,] 0.0000000 0.39156264 0.6084374 0.00000000
[6,] 0.9246218 0.05608074 0.0000000 0.01929748
This probability matrix shows for each cell number the probability of moving to its neighbouring cell (without showing the cell number of those individual neighbours). From here, I'm kinda stuck. I don't know how to actually simulate the movement of an individual from cell 1 (given that each choice is made independently of the previous step, sort of like a Markov chain & where there are different probabilities of moving based on your current step). I suspect it has something to do with indexing but I haven't figured out how to manage the different probabilities for each cell yet. This is my first time posting here, so hopefully this makes sense/is reproducible. Any help is much appreciated!
The best way would probably be to write code to convert your matrix into a 25x25 transition matrix and the use a Markov chain library, but it is reasonably straightforward to use your set up as is:
rand_walk <- function(start,steps){
walk = numeric(steps)
walker = start
for(i in 1:steps){
walk[i] <- walker
walker <- walker + sample(c(-5,5,-1,1),1,prob = prob_mat1[walker,])
}
walk
}
The basic idea is that moving up or down is adding or subtracting 5 from the current cell number and moving right or left is adding or subtracting 1, so it is enough to sample from the vector c(-5,5,-1,1) with the probabilities of those 4 steps given by the appropriate row of the probability matrix.
Typical output:
> rand_walk(1,100)
[1] 1 2 1 6 1 2 1 2 1 2 1 2 1 6 1 2 1 2 1 6 1 6
[23] 1 6 1 2 1 2 3 8 9 8 13 12 7 8 7 8 3 8 7 8 7 8
[45] 7 8 7 8 7 8 7 8 7 8 7 12 17 22 21 22 17 12 7 12 7 8
[67] 3 8 13 8 7 12 7 8 9 8 9 8 7 6 7 8 7 2 1 6 1 2
[89] 1 6 1 2 1 2 1 2 1 2 1 6
In this code I gave the complete walk (which is useful for debugging purposes) but you could of course drop the accumlating matrix walk completely and just return the final walker. Also, note that in this code, I recorded steps positions so only steps - 1 transitions.