Title: | Mechanistic Simulation of Species Range Dynamics |
---|---|
Description: | Species range dynamics simulation toolset. |
Authors: | Katarzyna Markowska [aut, cre], Lechosław Kuczyński [aut], Tad Dallas [rev] (@taddallas), Joanne Potts [rev] (@TheAnalyticalEdge) |
Maintainer: | Katarzyna Markowska <[email protected]> |
License: | MIT + file LICENSE |
Version: | 1.0.5 |
Built: | 2024-11-05 16:17:58 UTC |
Source: | https://github.com/ropensci/rangr |
This function simulates dispersal for each grid cell by calculating the number of individuals dispersing out of the cell and the number of individuals dispersing into the cell.
disp( N_t, id, id_matrix, data_table, kernel, dens_dep, dlist, id_within, within_mask, border, planar, dist_resolution, max_dist, dist_bin, ncells_in_circle, cl = NULL )
disp( N_t, id, id_matrix, data_table, kernel, dens_dep, dlist, id_within, within_mask, border, planar, dist_resolution, max_dist, dist_bin, ncells_in_circle, cl = NULL )
N_t |
integer matrix representing population numbers at a single time step; NA indicates cells outside the study area |
id |
|
id_matrix |
|
data_table |
matrix that contains information about all cells in current time points |
kernel |
function defining dispersal kernel |
dens_dep |
character vector of length 1 specifying if the probability
of settling in a target grid cell is (case-sensitive, default
|
dlist |
list with identifiers of target cells at a specified distance from a focal cell |
id_within |
integer vector with identifiers of cells inside the study area |
within_mask |
logical matrix that specifies boundaries of the study area |
border |
character vector of length 1 defining how to deal
with borders (case-sensitive, default
|
planar |
logical vector of length 1; |
dist_resolution |
integer vector of length 1; dimension of one side of
one cell of |
max_dist |
distance (in the same units as used in the raster |
dist_bin |
numeric vector of length 1 with value |
ncells_in_circle |
numeric vector; number of cells on each distance |
cl |
if simulation is done in parallel, the name of a cluster object
created by |
The function is used by sim
internally and is not intended to be
called by the user. The parameters for this function are passed from
a sim_data
object created by initialise
.
Dispersal distance is expressed in original spatial units of the
SpatRaster
provided to the sim
function
(n1_map
and K_map
). However, it is internally converted to units
of the simulation (i.e. the size of a single cell) by calculating round(distance/resolution)
. If the selected dispersal distance is
smaller than resolution/2
, the individual
does not disperse effectively and remains in the same cell.
The dispersal rate (proportion of dispersing individuals) can be estimated
from the dispersal kernel probability function by calculating the probability
that the dispersal distance is greater than resolution/2
.
The function returns a list that contains two matrices:
em
- emigration matrix with the number of individuals that dispersed
from each cell
im
- immigration matrix with the number of individuals that dispersed
to each cell
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # disp disp_output <- disp( N_t = sim_data$n1_map, id = unwrap(sim_data$id), id_matrix = as.matrix(unwrap(sim_data$id), wide = TRUE), data_table = sim_data$data_table, kernel = sim_data$kernel, dens_dep = sim_data$dens_dep, dlist = sim_data$dlist, id_within = sim_data$id_within, within_mask = sim_data$within_mask, border = sim_data$border, planar = sim_data$planar, dist_resolution = sim_data$dist_resolution, max_dist = sim_data$max_dist, dist_bin = sim_data$dist_bin, ncells_in_circle = sim_data$ncells_in_circle ) # immigration and emigration matrices names(disp_output)
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # disp disp_output <- disp( N_t = sim_data$n1_map, id = unwrap(sim_data$id), id_matrix = as.matrix(unwrap(sim_data$id), wide = TRUE), data_table = sim_data$data_table, kernel = sim_data$kernel, dens_dep = sim_data$dens_dep, dlist = sim_data$dlist, id_within = sim_data$id_within, within_mask = sim_data$within_mask, border = sim_data$border, planar = sim_data$planar, dist_resolution = sim_data$dist_resolution, max_dist = sim_data$max_dist, dist_bin = sim_data$dist_bin, ncells_in_circle = sim_data$ncells_in_circle ) # immigration and emigration matrices names(disp_output)
This function simulates an observation process. It accepts the sim_results
object, which is generated by the sim
function, and applies the virtual
ecologist approach on the N_map
component of the object. The function
returns a data.frame
with the 'observed' abundances.
get_observations( sim_data, sim_results, type = c("random_one_layer", "random_all_layers", "from_data", "monitoring_based"), obs_error = c("rlnorm", "rbinom"), obs_error_param = NULL, ... )
get_observations( sim_data, sim_results, type = c("random_one_layer", "random_all_layers", "from_data", "monitoring_based"), obs_error = c("rlnorm", "rbinom"), obs_error_param = NULL, ... )
sim_data |
|
sim_results |
|
type |
character vector of length 1; describes the sampling type (case-sensitive):
|
obs_error |
character vector of length 1; type of the distribution
that defines the observation process: " |
obs_error_param |
numeric vector of length 1; standard deviation
(on a log scale) of the random noise in observation process generated from
the log-normal distribution ( |
... |
other necessary internal parameters:
|
data.frame
object with geographic coordinates, time steps,
estimated abundance, observation error (if obs_error_param
is
provided), and observer identifiers (if type = "monitoring_based"
). If type = "from_data"
, returned object is sorted in the same order as the input points
.
## Not run: library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) # prepare data sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_1 <- sim(obj = sim_data, time = 110, burn = 10) # 1. random_one_layer sample1 <- get_observations( sim_data, sim_1, type = "random_one_layer", prop = 0.1 ) # 2. random_all_layers sample2 <- get_observations( sim_data, sim_1, type = "random_all_layers", prop = 0.15 ) # 3. from_data sample3 <- get_observations( sim_data, sim_1, type = "from_data", points = observations_points ) # 4. monitoring_based # define observations sites all_points <- xyFromCell(unwrap(sim_data$id), cells(unwrap(sim_data$K_map))) sample_idx <- sample(1:nrow(all_points), size = 20) sample_points <- all_points[sample_idx, ] sample4 <- get_observations( sim_data, sim_1, type = "monitoring_based", cells_coords = sample_points, prob = 0.3, progress_bar = TRUE ) # 5. noise "rlnorm" sample5 <- get_observations(sim_data, sim_1, type = "random_one_layer", obs_error = "rlnorm", obs_error_param = log(1.2) ) # 6. noise "rbinom" sample6 <- get_observations(sim_data, sim_1, type = "random_one_layer", obs_error = "rbinom", obs_error_param = 0.8 ) ## End(Not run)
## Not run: library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) # prepare data sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_1 <- sim(obj = sim_data, time = 110, burn = 10) # 1. random_one_layer sample1 <- get_observations( sim_data, sim_1, type = "random_one_layer", prop = 0.1 ) # 2. random_all_layers sample2 <- get_observations( sim_data, sim_1, type = "random_all_layers", prop = 0.15 ) # 3. from_data sample3 <- get_observations( sim_data, sim_1, type = "from_data", points = observations_points ) # 4. monitoring_based # define observations sites all_points <- xyFromCell(unwrap(sim_data$id), cells(unwrap(sim_data$K_map))) sample_idx <- sample(1:nrow(all_points), size = 20) sample_points <- all_points[sample_idx, ] sample4 <- get_observations( sim_data, sim_1, type = "monitoring_based", cells_coords = sample_points, prob = 0.3, progress_bar = TRUE ) # 5. noise "rlnorm" sample5 <- get_observations(sim_data, sim_1, type = "random_one_layer", obs_error = "rlnorm", obs_error_param = log(1.2) ) # 6. noise "rbinom" sample6 <- get_observations(sim_data, sim_1, type = "random_one_layer", obs_error = "rbinom", obs_error_param = 0.8 ) ## End(Not run)
Population growth functions are used during simulation
conducted by the sim
function.
The user is required to specify the name of a growth function while initialising the
sim_data
object using initialise
.
exponential(x, r, ...) ricker(x, r, K, A = NA) gompertz(x, r, K, A = NA)
exponential(x, r, ...) ricker(x, r, K, A = NA) gompertz(x, r, K, A = NA)
x |
number of individuals |
r |
intrinsic population growth rate |
... |
not used, added for compatibility reasons |
K |
carrying capacity |
A |
coefficient of Allee effect (A <= 0: weak, A > 0: strong, NA: none) |
x
can be a vector, matrix, SpatRaster
or any other R
object for which basic arithmetic operations produce valid results.
These functions are intended to be used in the sim
function, where x
is a matrix of the same dimensions as the SpatRaster
object specified in n1_map
parameter.
Object of the same dimensions as x
that contains expected number
of individuals in the next time step.
Boukal, D. S., & Berec, L. (2002). Single-species models of the Allee effect: extinction boundaries, sex ratios and mate encounters. Journal of Theoretical Biology, 218(3), 375-394. doi:10.1006/jtbi.2002.3084
Gompertz, B. (1825) On the Nature of the Function Expressive of the Law of Human Mortality, and on a New Mode of Determining the Value of Life Contigencies. Philosophical Transactions of the Royal Society of London, 115, 513-583. doi:10.1098/rstl.1825.0026
Ricker, W.E. (1954) Stock and Recruitment. Journal of the Fisheries Research Board of Canada, 11, 559-623. doi:10.1139/f54-039
Hostetler, J.A. and Chandler, R.B. (2015), Improved state-space models for inference about spatial and temporal variation in abundance from count data. Ecology, 96: 1713-1723. doi:10.1890/14-1487.1
Courchamp, F., L. Berec and J. Gascoigne. 2008. Allee Effects in Ecology and Conservation. Oxford University Press, New York. 256 pp. ISBN 978-0-19-857030-1
x <- 1:10 exponential(x, r = 0.4) ricker(x, r = 2, K = 5) ricker(x, r = 2, K = 5, A = -5) gompertz(x, r = 1.2, K = 5) gompertz(x, r = 1.2, K = 5, A = 5)
x <- 1:10 exponential(x, r = 0.4) ricker(x, r = 2, K = 5) ricker(x, r = 2, K = 5, A = -5) gompertz(x, r = 1.2, K = 5) gompertz(x, r = 1.2, K = 5, A = 5)
This function generates a sim_data
object containing all the necessary
information required to run a simulation by the sim
function. The
input maps (n1_map
and K_map
) can be in the Cartesian or longitude/latitude
coordinate system.
initialise( n1_map, K_map, K_sd = 0, r, r_sd = 0, growth = "gompertz", A = NA, dens_dep = c("K2N", "K", "none"), border = c("reprising", "absorbing"), kernel_fun = "rexp", ..., max_dist = NA, calculate_dist = TRUE, dlist = NULL, progress_bar = TRUE, quiet = FALSE )
initialise( n1_map, K_map, K_sd = 0, r, r_sd = 0, growth = "gompertz", A = NA, dens_dep = c("K2N", "K", "none"), border = c("reprising", "absorbing"), kernel_fun = "rexp", ..., max_dist = NA, calculate_dist = TRUE, dlist = NULL, progress_bar = TRUE, quiet = FALSE )
n1_map |
|
K_map |
|
K_sd |
numeric vector of length 1 with value |
r |
numeric vector of length 1; intrinsic population growth rate |
r_sd |
numeric vector of length 1 with value |
growth |
character vector of length 1; the name of a population growth
function, either defined in |
A |
numeric vector of length 1; strength of the Allee effect
(see the |
dens_dep |
character vector of length 1 specifying if the probability
of settling in a target grid cell is (case-sensitive, default
|
border |
character vector of length 1 defining how to deal
with borders (case-sensitive, default
|
kernel_fun |
character vector of length 1; name of a random number
generation function defining a dispersal kernel (case-sensitive,
default |
... |
any parameters required by |
max_dist |
numeric vector of length 1; maximum distance of dispersal to pre-calculate target cells |
calculate_dist |
logical vector of length 1; determines if target cells will be precalculated |
dlist |
list; target cells at a specified distance calculated for every cell within the study area |
progress_bar |
logical vector of length 1; determines if progress bar for calculating distances should be displayed |
quiet |
logical vector of length 1; determines if messages should be displayed |
The most time-consuming part of computations performed by the sim
function is the simulation of dispersal. To speed it up, a list containing
indexes of target cells at a specified distance from a focal cell
is calculated in advance and stored in a dlist
slot.
The max_dist
parameter sets
the maximum distance at which this pre-calculation is performed. If max_dist
is NULL
, it is set to 0.99 quantile from the kernel_fun
.
All distance calculations are always based on metres if the input maps are
latitude/longitude. For planar input maps, distances are calculated in map
units, which are typically metres, but check the crs()
if in doubt.
If the input maps are in the Cartesian coordinate system and the grid cells
are squares, then
the distances between cells are calculated using the distance
function from the terra
package. These distances are later divided by the
resolution of the input maps.
For input maps with grid cells in shapes other than squares (e.g. with
rectangular cells or longitude/latitude coordinate system), the distance
resolution is calculated by finding the shortest distance between each
"queen" type neighbor. All distances calculated by the distance
function are further divided by this distance resolution.
To avoid discontinuities in the distances at which the target cells are located,
an additional parameter dist_bin
is calculated as half of the maximum
distance between each "queen" type neighbour. It is used to expand the
distances at which target cells are located from a single number to a range.
NA in the input maps represents cells outside the study area.
The K_get_interpolation
function can be used to prepare K_map
that changes
over time. This may be useful, when simulating environmental change or
exploring the effects of ecological disturbances.
Object of class sim_data
which inherits from list
. This object
contains all necessary information to perform a simulation using
sim
function.
Hijmans R (2024). terra: Spatial Data Analysis. R package version 1.7-81, https://rspatial.github.io/terra/, https://rspatial.org/
Solymos P, Zawadzki Z (2023). pbapply: Adding Progress Bar to '*apply' Functions. R package version 1.7-2, https://CRAN.R-project.org/package=pbapply.
## Not run: # input maps library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) K_small_changing <- rast(system.file("input_maps/K_small_changing.tif", package = "rangr")) n1_small_lon_lat <- rast(system.file("input_maps/n1_small_lon_lat.tif", package = "rangr")) K_small_lon_lat <- rast(system.file("input_maps/K_small_lon_lat.tif", package = "rangr")) # basic example sim_data_1 <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # example with changing environment K_interpolated <- K_get_interpolation( K_small_changing, K_time_points = c(1, 25, 50) ) sim_data_2 <- initialise( n1_map = n1_small, K_map = K_interpolated, r = log(2), rate = 1 / 1e3 ) # example with lon/lat rasters sim_data_3 <- initialise( n1_map = n1_small_lon_lat, K_map = K_small_lon_lat, r = log(2), rate = 1 / 1e3 ) # example without progress bar and messages sim_data_4 <- initialise( n1_map = n1_small, K_map = K_small, K_sd = 0.1, r = log(5), r_sd = 4, growth = "ricker", rate = 1 / 200, max_dist = 5000, dens_dep = "K2N", progress_bar = FALSE, quiet = TRUE ) ## End(Not run)
## Not run: # input maps library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) K_small_changing <- rast(system.file("input_maps/K_small_changing.tif", package = "rangr")) n1_small_lon_lat <- rast(system.file("input_maps/n1_small_lon_lat.tif", package = "rangr")) K_small_lon_lat <- rast(system.file("input_maps/K_small_lon_lat.tif", package = "rangr")) # basic example sim_data_1 <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # example with changing environment K_interpolated <- K_get_interpolation( K_small_changing, K_time_points = c(1, 25, 50) ) sim_data_2 <- initialise( n1_map = n1_small, K_map = K_interpolated, r = log(2), rate = 1 / 1e3 ) # example with lon/lat rasters sim_data_3 <- initialise( n1_map = n1_small_lon_lat, K_map = K_small_lon_lat, r = log(2), rate = 1 / 1e3 ) # example without progress bar and messages sim_data_4 <- initialise( n1_map = n1_small, K_map = K_small, K_sd = 0.1, r = log(5), r_sd = 4, growth = "ricker", rate = 1 / 200, max_dist = 5000, dens_dep = "K2N", progress_bar = FALSE, quiet = TRUE ) ## End(Not run)
SpatRaster
object representing a carrying
capacity map projected to WGS 84 (CRS84) from the original raster K_big
.
This map can be used as a carrying capacity map to initialise
data necessary
to perform a simulation with the sim
function. It is compatible with the n1_big_lon_lat.tif
raster.
SpatRaster
object with 74 rows
and 125 columns containing integer values 0-25 and NA's
indicating unsuitable areas.
Data generated in-house to serve as an example (using spatial autocorrelation).
system.file("input_maps/K_big_lon_lat.tif", package = "rangr")
system.file("input_maps/K_big_lon_lat.tif", package = "rangr")
SpatRaster
object that can be used as a
carrying capacity map to initialise
data necessary
to perform a simulation with the sim
function. This map is compatible
with n1_big.tif
.
SpatRaster
object with 100 rows
and 100 columns containing integer values 0-25 and NA's indicating
unsuitable areas.
Data generated in-house to serve as an example (using spatial autocorrelation).
system.file("input_maps/K_big.tif", package = "rangr")
system.file("input_maps/K_big.tif", package = "rangr")
This function linearly interpolates values in a series of carrying capacity maps.
K_get_interpolation(K_map, K_time_points = NULL, time = NULL)
K_get_interpolation(K_map, K_time_points = NULL, time = NULL)
K_map |
|
K_time_points |
integer vector; time for each layer in |
time |
integer vector of length 1; number of total time steps required
(this is defined when evoking the function |
To simulate dynamic environmental scenarios (e.g. climate change, land use change, ecological disturbance, etc.) one needs to provide time-varying carrying capacity maps.
Either K_time_points
or the time
parameter is needed to perform
interpolation. If the interpolation should be calculated between two carrying
capacity maps, there is no need to pass the time points, because 1 will
be set as the starting time point and time
will be used as the ending point.
On the other hand, in the absence of the time
argument, the maximum element
of K_time_points
is considered to be the ending point for the interpolation.
SpatRaster
object with number of layers
equal to time
.
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small_changing <- rast(system.file("input_maps/K_small_changing.tif", package = "rangr")) K_interpolated_01 <- K_get_interpolation( K_small_changing, K_time_points = c(1, 10, 15) ) K_two_layers <- subset( K_small_changing, c(1, 2) ) K_interpolated_02 <- K_get_interpolation( K_two_layers, time = 15 ) ## End(Not run)
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small_changing <- rast(system.file("input_maps/K_small_changing.tif", package = "rangr")) K_interpolated_01 <- K_get_interpolation( K_small_changing, K_time_points = c(1, 10, 15) ) K_two_layers <- subset( K_small_changing, c(1, 2) ) K_interpolated_02 <- K_get_interpolation( K_two_layers, time = 15 ) ## End(Not run)
SpatRaster
object representing changing carrying
capacity maps projected to WGS 84 (CRS84) from the original raster
K_small_changing
. These maps can be used as carrying capacity maps to initialise
data necessary to perform a simulation with the sim
function. To utilise
these maps in initialise
the user must first use K_get_interpolation
to generate a map for every time step of the simulation. These maps are
compatible with the n1_small_lon_lat.tif
raster.
SpatRaster
object with 3 layers,
each having 12 rows and 14 columns containing integer values 0-170 and NA's indicating unsuitable areas.
Data generated in-house to serve as an example (using spatial autocorrelation).
system.file("input_maps/K_small_changing_lon_lat.tif", package = "rangr")
system.file("input_maps/K_small_changing_lon_lat.tif", package = "rangr")
SpatRaster
object that can be used as carrying
capacity maps to initialise
data necessary to perform a simulation
with the sim
function.
To utilise these maps in initialise
the user first must
use K_get_interpolation
to generate a map for every time step
of the simulation.
These maps are compatible with n1_small.tif
.
Each subsequent map contains a virtual environment with greater
carrying capacity than the previous one.
SpatRaster
object with 3 layers,
each has 15 rows and 10 columns containing integer values 0-170 and
NA's that indicates unsuitable areas.
Data generated in-house to serve as an example (using spatial autocorrelation).
system.file("input_maps/K_small_changing.tif", package = "rangr")
system.file("input_maps/K_small_changing.tif", package = "rangr")
SpatRaster
object that represents a carrying
capacity map projected to WGS 84 (CRS84) from the original raster K_small
.
This map can be used as a carrying capacity map to initialise
data necessary
to perform a simulation with the sim
function. It is compatible with the n1_small_lon_lat.tif
raster.
SpatRaster
object with 12 rows
and 14 columns containing integer values 0-100 and NA's indicating unsuitable areas.
Data generated in-house to serve as an example (using spatial autocorrelation).
system.file("input_maps/K_small_lon_lat.tif", package = "rangr")
system.file("input_maps/K_small_lon_lat.tif", package = "rangr")
SpatRaster
object that can be used a carrying
capacity map to initialise
data necessary to perform a simulation
with the sim
function.
This map is compatible with n1_small.tif
.
SpatRaster
object with 15 rows
and 10 columns containing integer values 0-100 and NA's indicating
unsuitable areas.
Data generated in-house to serve as an example (using spatial autocorrelation).
system.file("input_maps/K_small.tif", package = "rangr")
system.file("input_maps/K_small.tif", package = "rangr")
SpatRaster
object representing an abundance map
at the first time step of the simulation projected to WGS 84 (CRS84) from the
original raster n1_big
. This map can be used as a simulation starting point
to initialise
data necessary to perform a simulation with the sim
function. It is compatible with the K_big_lon_lat.tif
map.
SpatRaster
object with 74 rows
and 125 columns containing integer values 0-50 and NA's
indicating unsuitable areas.
Data generated in-house to serve as an example.
system.file("input_maps/n1_big_lon_lat.tif", package = "rangr")
system.file("input_maps/n1_big_lon_lat.tif", package = "rangr")
SpatRaster
object that can be used a
as simulation starting point to initialise
data necessary
to perform a simulation with the sim
function. This map is compatible
with K_big.tif
map.
SpatRaster
object with 100 rows
and 100 columns containing integer values 0-50 and NA's
that indicates unsuitable areas.
Data generated in-house to serve as an example.
system.file("input_maps/n1_big.tif", package = "rangr")
system.file("input_maps/n1_big.tif", package = "rangr")
SpatRaster
object representing an abundance map
at the first time step of the simulation projected to WGS 84 (CRS84) from the
original raster n1_small
. This map can be used as a simulation starting point
to initialise
data necessary to perform a simulation with the sim
function. It is compatible with the
K_small_lon_lat.tif
and K_small_changing_lon_lat.tif
maps.
SpatRaster
object with 12 rows
and 14 columns containing integer values 0-10 and NA's
indicating unsuitable areas.
Data generated in-house to serve as an example.
system.file("input_maps/n1_small_lon_lat.tif", package = "rangr")
system.file("input_maps/n1_small_lon_lat.tif", package = "rangr")
SpatRaster
object that can be used a
as simulation starting point to initialise
data necessary
to perform a simulation with the sim
function. This map is compatible
with K_small.tif
and K_small_changing.tif
maps.
SpatRaster
object with 15 rows
and 10 columns containing integer values 0-10 and NA's
indicating unsuitable areas.
Data generated in-house to serve as an example.
system.file("input_maps/n1_small.tif", package = "rangr")
system.file("input_maps/n1_small.tif", package = "rangr")
A data.frame
containing a sample input data to the function
get_observations
when type
argument is set to "from_file".
This data is compatible with n1_small.tif
,
K_small.tif
and K_small_changing.tif
maps.
observations_points
observations_points
A data frame with 1500 rows and 3 variables:
x coordinate
y coordinate
time_step at which the abundances should be observed
Data generated in-house to serve as an example
sim_results
ObjectPlots abundances obtained during simulation.
## S3 method for class 'sim_results' plot(x, template = NULL, time_points = NULL, range, type, ...)
## S3 method for class 'sim_results' plot(x, template = NULL, time_points = NULL, range, type, ...)
x |
|
template |
|
time_points |
numeric vector; specifies points in time from which plots will be generated |
range |
numeric vector of length 2; range of values to be used for the
legend (if |
type |
character vector of length 1; type of map: "continuous" (default), "classes" or "interval" (case-sensitive) |
... |
further arguments passed to |
SpatRaster
object with as many layers
as the length of time_points
parameter
## Not run: library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_res <- sim(sim_data, time = 10) plot(sim_res) plot(sim_res, template = n1_small, time_points = c(1, 10)) # plot specific area plot(sim_res, xlim = c(4, 10), ylim = c(0, 10)) plot(sim_res, ext = c(4, 10, 0, 10)) plot(sim_res, template = n1_small, ext = c(274000, 280000, 610000, 620000)) ## End(Not run)
## Not run: library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_res <- sim(sim_data, time = 10) plot(sim_res) plot(sim_res, template = n1_small, time_points = c(1, 10)) # plot specific area plot(sim_res, xlim = c(4, 10), ylim = c(0, 10)) plot(sim_res, ext = c(4, 10, 0, 10)) plot(sim_res, template = n1_small, ext = c(274000, 280000, 610000, 620000)) ## End(Not run)
sim_data
ObjectPrint sim_data
Object
## S3 method for class 'sim_data' print(x, ...)
## S3 method for class 'sim_data' print(x, ...)
x |
|
... |
further arguments passed to or from other methods; currently none specified |
sim_data
object is invisibly returned (the x
param)
library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) print(sim_data)
library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) print(sim_data)
sim_results
ObjectPrint sim_results
Object
## S3 method for class 'sim_results' print(x, ...)
## S3 method for class 'sim_results' print(x, ...)
x |
|
... |
further arguments passed to or from other methods; none specified |
sim_results
object is invisibly returned (the x
param)
## Not run: library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_res <- sim(obj = sim_data, time = 20, burn = 5) print(sim_res) ## End(Not run)
## Not run: library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_res <- sim(obj = sim_data, time = 20, burn = 5) print(sim_res) ## End(Not run)
summary.sim_data
ObjectPrint summary.sim_data
Object
## S3 method for class 'summary.sim_data' print(x, ...)
## S3 method for class 'summary.sim_data' print(x, ...)
x |
|
... |
further arguments passed to or from other methods; currently none specified |
None
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) summary_sim_data <- summary(sim_data) print(summary_sim_data)
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) summary_sim_data <- summary(sim_data) print(summary_sim_data)
summary.sim_results
ObjectPrint summary.sim_results
Object
## S3 method for class 'summary.sim_results' print(x, ...)
## S3 method for class 'summary.sim_results' print(x, ...)
x |
|
... |
further arguments passed to or from other methods; currently none specified |
None
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_results <- sim(sim_data, time = 10) summary_sim_results <- summary(sim_results) print(summary_sim_results)
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) sim_results <- sim(sim_data, time = 10) summary_sim_results <- summary(sim_results) print(summary_sim_results)
This function simulates population growth and dispersal providing a given
environmental scenario. All parameters for the simulation must be set
in advance using initialise
.
sim( obj, time, burn = 0, return_mu = FALSE, cl = NULL, progress_bar = TRUE, quiet = FALSE )
sim( obj, time, burn = 0, return_mu = FALSE, cl = NULL, progress_bar = TRUE, quiet = FALSE )
obj |
|
time |
positive integer vector of length 1; number of time steps simulated |
burn |
positive integer vector of length 1; the number of burn-in time steps that are discarded from the output |
return_mu |
logical vector of length 1; if |
cl |
an optional cluster object created by
|
progress_bar |
logical vector of length 1 determines if progress bar for simulation should be displayed |
quiet |
logical vector of length 1; determines if warnings should be displayed |
This is the main simulation module. It takes the sim_data
object prepared
by initialise
and runs simulation for a given number of time steps.
The initial (specified by the burn
parameter) time steps are skipped
and discarded from the output. Computations can be done in parallel if
the name of a cluster created by makeCluster
is provided.
Generally, at each time step, simulation consists of two phases: local dynamics and dispersal.
Local dynamics (which connects habitat patches in time) is defined by
the function growth
. This parameter is specified while creating
the obj
using initialise
, but can be later modified by using
the update
function. Population growth can be either exponential
or density-dependent, and the regulation is implemented by the use of
Gompertz or Ricker models (with a possibility of providing any other,
user defined function). For every cell, the expected population density
during the next time step is calculated from the corresponding number
of individuals currently present in this cell, and the actual number
of individuals is set by drawing a random number from the Poisson
distribution using this expected value. This procedure introduces a realistic
randomness, however additional levels of random variability can be
incorporated by providing parameters of both demographic and environmental
stochasticity while specifying the sim_data
object using the initialise
function (parameters r_sd
and K_sd
, respectively).
To simulate dispersal (which connects habitat patches in space), for each
individual in a given cell, a dispersal distance is randomly drawn from
the dispersal kernel density function. Then, each individual is translocated
to a randomly chosen cell at this distance apart from the current location.
For more details, see the disp
function.
This function returns an object of class sim_results
which is
a list containing the following components:
extinction
- TRUE
if population is extinct or FALSE
otherwise
simulated_time
- number of simulated time steps without
the burn-in ones
N_map
- 3-dimensional array representing spatio-temporal
variability in population numbers. The first two dimensions correspond to
the spatial aspect of the output and the third dimension represents time.
In case of a total extinction, a simulation is stopped before reaching
the specified number of time steps. If the population died out before reaching
the burn
threshold, then nothing can be returned and an error occurs.
Solymos P, Zawadzki Z (2023). pbapply: Adding Progress Bar to '*apply' Functions. R package version 1.7-2, https://CRAN.R-project.org/package=pbapply.
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # simulation sim_1 <- sim(obj = sim_data, time = 20) # simulation with burned time steps sim_2 <- sim(obj = sim_data, time = 20, burn = 10) # example with parallelization library(parallel) cl <- makeCluster(detectCores()) # parallelized simulation sim_3 <- sim(obj = sim_data, time = 20, cl = cl) stopCluster(cl) # visualisation plot( sim_1, time_points = 20, template = sim_data$K_map ) plot( sim_1, time_points = c(1, 5, 10, 20), template = sim_data$K_map ) plot( sim_1, template = sim_data$K_map ) ## End(Not run)
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # simulation sim_1 <- sim(obj = sim_data, time = 20) # simulation with burned time steps sim_2 <- sim(obj = sim_data, time = 20, burn = 10) # example with parallelization library(parallel) cl <- makeCluster(detectCores()) # parallelized simulation sim_3 <- sim(obj = sim_data, time = 20, cl = cl) stopCluster(cl) # visualisation plot( sim_1, time_points = 20, template = sim_data$K_map ) plot( sim_1, time_points = c(1, 5, 10, 20), template = sim_data$K_map ) plot( sim_1, template = sim_data$K_map ) ## End(Not run)
sim_results
ObjectThis function creates a subset of given time points from the sim_results
object.
## S3 method for class 'sim_results' subset(x, from = NULL, time_points = NULL, ...)
## S3 method for class 'sim_results' subset(x, from = NULL, time_points = NULL, ...)
x |
|
from |
numeric vector of length 1; indicates the starting time point from which all time point should be kept |
time_points |
numeric vector; indicates all time points to keep |
... |
further arguments to be passed to or from other methods |
Either from
or time_points
argument has to be specified.
Time point passed by the from
argument will be set as a cutoff point
and all abundances for previous time points will be discarded.
sim_results
object with only selected time_points
present
in the N_map
slot
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n = n1_small, r = log(2), K_map = K_small, max_dist = 1000, rate = 1 / 1e3 ) sim_results <- sim(sim_data, time = 10) summary(sim_results) sim_results_cropped <- subset(sim_results, time_points = c(1:2)) summary(sim_results_cropped) ## End(Not run)
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n = n1_small, r = log(2), K_map = K_small, max_dist = 1000, rate = 1 / 1e3 ) sim_results <- sim(sim_data, time = 10) summary(sim_results) sim_results_cropped <- subset(sim_results, time_points = c(1:2)) summary(sim_results_cropped) ## End(Not run)
sim_data
ObjectSummary Of sim_data
Object
## S3 method for class 'sim_data' summary(object, ...)
## S3 method for class 'sim_data' summary(object, ...)
object |
|
... |
further arguments passed to or from other methods; currently none specified |
summary.sim_data
object
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) summary(sim_data)
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) summary(sim_data)
sim_results
ObjectSummary Of sim_results
Object
## S3 method for class 'sim_results' summary(object, ...)
## S3 method for class 'sim_results' summary(object, ...)
object |
|
... |
further arguments passed to or from other methods; none specified |
summary.sim_results
object
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # simulation sim_results <- sim(sim_data, time = 10) summary(sim_results)
# data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # simulation sim_results <- sim(sim_data, time = 10) summary(sim_results)
sim_results
To RasterThis function transforms selected subset of abundance matrices from
sim_results
into SpatRaster
object. Layers are
specified by time_points
, which can be one or multiple points in time.
to_rast(sim_results, time_points = sim_results$simulated_time, template = NULL)
to_rast(sim_results, time_points = sim_results$simulated_time, template = NULL)
sim_results |
|
time_points |
numeric vector of length 1 or more; specifies points in
time from which |
template |
|
SpatRaster
based on sim_results
object
with layers corresponding to time_points
.
Hijmans R (2024). terra: Spatial Data Analysis. R package version 1.7-81, https://rspatial.github.io/terra/, https://rspatial.org/
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # simulation sim_1 <- sim(obj = sim_data, time = 100) # raster construction my_rast <- to_rast( sim_1, time_points = c(1, 10, 20, 100), template = sim_data$K_map ) # visualization plot(my_rast, range = range(sim_1$N_map, na.rm = TRUE)) ## End(Not run)
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) # simulation sim_1 <- sim(obj = sim_data, time = 100) # raster construction my_rast <- to_rast( sim_1, time_points = c(1, 10, 20, 100), template = sim_data$K_map ) # visualization plot(my_rast, range = range(sim_1$N_map, na.rm = TRUE)) ## End(Not run)
sim_data
ObjectThis function updates a sim_data
object.
## S3 method for class 'sim_data' update(object, ..., evaluate = TRUE)
## S3 method for class 'sim_data' update(object, ..., evaluate = TRUE)
object |
|
... |
further arguments passed to or from other methods; currently none specified |
evaluate |
logical vector of length 1; if |
If evaluate = TRUE
then the updated sim_data
object,
otherwise the updated call.
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data_1 <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) summary(sim_data_1) sim_data_2 <- update(sim_data_1, max_dist = 3000) summary(sim_data_2) ## End(Not run)
## Not run: # data preparation library(terra) n1_small <- rast(system.file("input_maps/n1_small.tif", package = "rangr")) K_small <- rast(system.file("input_maps/K_small.tif", package = "rangr")) sim_data_1 <- initialise( n1_map = n1_small, K_map = K_small, r = log(2), rate = 1 / 1e3 ) summary(sim_data_1) sim_data_2 <- update(sim_data_1, max_dist = 3000) summary(sim_data_2) ## End(Not run)