-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathplants.R
169 lines (159 loc) · 7.94 KB
/
plants.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
#' setup.plants is a function that creates a list of plant info
#' @param repro reproduction probability, must be a repro prob for every plant species, a vector with each value between 0 and 1
#' @param survival survival probability, must be a survival prob for every plant species, a vector with each value between 0 and 1
#' @param comp.matrix competition matrix of survival probabilities between species, rows and columns must be equal to the number of plant species, each cell is populated with a value between 0 and 1
#' @param names names of the plants being simulated, if null plants are "a", "b", "c", etc.
#' @return a list with the repro probability, survival probability, comp matrix probabilities, and names
setup.plants <- function(repro, survival, comp.matrix, names=NULL){
if(is.null(names)){
names <- letters[seq_along(repro)]
}
if(length(repro) != length(survival)){
stop("Reproduction and Survival parameters must be the same length")
}
if(class(comp.matrix) != "matrix"){
stop("Comp.matrix parameter must be a matrix class")
}
if(nrow(comp.matrix) != length(repro) & ncol(comp.matrix) != length(repro)){
stop("Comp.matrix dimensions must be same length as reproduction and survival")
}
repro <- setNames(repro, names)
survival <- setNames(survival, names)
rownames(comp.matrix) <- names
colnames(comp.matrix) <- names
return(list(repro=repro, survival=survival, comp.matrix=comp.matrix, names=names))
}
#creating some initial parameters to test
#r <- c(0.5,0.5)
#s <- c(0.7, 0.6)
#cm <- matrix(data=c(0.9, 0.7, 0.3, 0.5), 2, 2)
#names <- c("Bill", "Ted")
#info <- setup.plants(repro = r, survival = s, comp.matrix = cm, names = names)
#t <- terrain.fun(65)
#' init.plants is a function that creates an array and initiates plants across a terrain
#' @param terrain a matrix and the output of the terrain.fun function
#' @param timesteps the number of times the simulation will run
#' @param names names of the plants being simulated, if null plants are "a", "b", "c", etc.
#' @return an array where columns and rows equal columns and rows of terrain and the third dimension equals with timestep + 1 (an initial matrix, plus one matrix for each timestep)
init.plants <- function(terrain, timesteps, names){
if(is.null(names)){
names <- letters[seq_along(repro)]
}
plants <- array("", dim=c(dim(terrain), timesteps+1))
for(i in 1:nrow(terrain)){
for(j in 1:ncol(terrain)){
plants[i,j,1] <- sample(names,1)
}
}
#puts NAs into plants array where there is water in terrain
for(i in seq_len(dim(plants)[timesteps+1])){
plants[,,i][is.na(terrain)] <- NA
}
return(plants)
}
#creating a plant array to work with
#plants <- init.plants(terrain = t, timesteps = 2, names=names)
#' survival function determins the if a plant in a given cell survives each timestep
#' @param cell the cell (row, col, k) in the array (the output from init.plants) that is being tested
#' @param info a list and the output from setup.plants
#' @return returns "" if the plant died and the name of the plant if it survived
survival <- function(plants, row, column, time, info){
#checks to see if the cell is underwater
if(is.na(plants[row, column, time])){
plants[row,column, time] <- NA
}
for(n in 1:length(info$names)){
#if the cell has a plant in it
if(isTRUE(plants[row, column, time] == info$names[n])){
#checks to see if the plant in the cell survived
if(runif(1) > info$survival[plants[row, column, time]]){
#plant died, cell is empty
plants[row, column, time] <- ""
}
}
}
return(plants[row, column, time])
}
#' reproduce determines if and where a plant will reproduce. Limited to adjacent cells. Also includes a competition function if a cell is already occupied
#' @param timesteps the number of times the simulation will run
#' @param row the number of rows in the matrix. Determined from the number of rows in the Plants array
#' @param column the number of columns in the matrix. Determined from the number of columns in the Plants array
#' @param plants an array and output of init.plants function
#' @param info a list and output of the setup.plants function
#' @return returns updated array (timestep) with plant information depending on outcome of simulated repro and competition
reproduce <- function(time, row, column, plants, info){
#loops through all rows in plants array
for(p in 1:nrow(plants)){
#creates a matrix of possible repro locations
pos.loc <- as.matrix(expand.grid(row+c(-1,0,1), column+c(-1,0,1)))
#If pos repro loc is out of bounds, populates pos.loc with NA
pos.loc[pos.loc < 1] <- NA
pos.loc[pos.loc > dim(plants)[1]] <- NA
#Deletes all rows containing NAs from pos.loc
for(i in 9:1){
if(is.na(pos.loc[i,1]) | is.na(pos.loc[i,2])){
pos.loc <- pos.loc[-i,]
}
}
#loops through all the plant species
for(n in 1:length(info$names)){
#chooses a cell to reproduce into
pr <- sample(pos.loc[,1],1)
pc <- sample(pos.loc[,2],1)
#if a cell has a plant in it
if(isTRUE(plants[row,column,time] == info$names[n])){
#if the plant reproduced
if(runif(1) <= info$repro[plants[row, column, time]]){
#if the repro location does not contain NA (i.e., isn't under water)
if(!is.na(plants[pr, pc, time])){
#if an adjacent cell is a possible repro location and already has a plant in it
if(isTRUE(plants[pr, pc, time] == info$names[n])){
#if the propagule outcompetes the original plant
if(runif(1) <= info$comp.matrix[plants[row,column,time], plants[pr, pc, time]]){
#then the cell gets populated with the new plant
plants[pr, pc, time] <- plants[row,column,time]
}
}
}
}
}
}
}
return(plants[row, column, time])
}
#' plant.timestep is a wrapper around the survival and reproduce functions. Loops through every row, column, and timestep in the plants array
#' @param plants an array and output of init.plants function
#' @param terrain a matrix and the output of the terrain.fun function
#' @param info a list and output of the setup.plants function
#' @param timesteps the number of times the simulation will run
#' @return a new matrix with the updated plant location info based on the outcome of reproduction, survival, and competition
plant.timestep <- function(plants, terrain, info, timesteps){
#loops through the timesteps
for(t in 1:(timesteps+1)){
#loops through the rows
for(i in 1:nrow(plants)){
#loops through the columns
for(j in 1:ncol(plants)){
#applies the survival function to each cell in the matrix
survival(plants, i, j, t, info)
reproduce(t, i, j, plants, info)
}
}
}
return(plants)
}
#' run.plant.eco is a function that wraps around the setup.plants, init.plants, and plant.timestep functions
#' @param timesteps the number of times the simulation will run
#' @param terrain a matrix and the output of the terrain.fun function
#' @param repro reproduction probability, must be a repro prob for every plant species, a vector with each value between 0 and 1
#' @param survival survival probability, must be a survival prob for every plant species, a vector with each value between 0 and 1
#' @param comp.matrix competition matrix of survival probabilities between species, rows and columns must be equal to the number of plant species, each cell is populated with a value between 0 and 1
#' @param names names of the plants being simulated, if null plants are "a", "b", "c", etc.
#' @return returns an array with plant locations overtime
#' @export
run.plant.eco <- function(timesteps, terrain, repro, survival, comp.matrix, names=NULL){
info <- setup.plants(repro=repro, survival=survival, comp.matrix=comp.matrix, names = names)
plants <- init.plants(terrain = terrain, timesteps = timesteps, names = names)
output <- plant.timestep(plants=plants, terrain=terrain, info=info, timesteps=timesteps)
return(output)
}