Nothing
#' @title Simulate new points on a specific sf object
#'
#' @description Simulate new points on a specific \link[sf]{sf} object. See \link[sf]{st_sample} for details.
#'
#'
#' @param sf object of class sf or sfc
#' @param size sample size(s) requested; either total size, or a numeric vector with sample sizes for each feature geometry.
#' When sampling polygons, the returned sampling size may differ from the requested size, as the bounding box is
#' sampled, and sampled points intersecting the polygon are returned.
#'
#' @export
#'
brk_newPoints <- function(sf, size = 100){
ls_points <- st_sf( geometry = st_sample(sf, size = size, type = "random", exact = TRUE))
return(ls_points)
}
#' @title Add sequence to element of a \link[sf]{sf} object
#'
#' @description Add column with sequence for each row to a \link[sf]{sf} object
#'
#' @param sf object of class sf or sfc
#' @param key a character string used as name of the new column
#' @param from the starting value of the sequence. Of length 1 unless just from is supplied as an unnamed argument.
#' @param to the end (maximal) value of the sequence. Of length 1 unless just from is supplied as an unnamed argument.
#' @param by number: increment of the sequence.
#'
#'
#' @export
#'
brk_timeline <- function(sf, key, from, to, by){
if(length(from)==1) from = rep(from, nrow(sf))
if(length(to)==1) to = rep(to, nrow(sf))
sf[[key]] <- lapply(1:nrow(sf), function(i) return(seq(from[i], to[i], by)))
return(sf)
}
#' @title Add raster value to element of sf object
#'
#' @description Add the raster value(s) of a Raster* object to element of sf object.
#'
#' @param stackRaster_exposure The Raster* object
#' @param sf the sf object
#' @param stackTimeline sequence with matching to elements of RasterStack (length of stackTimeline
#' must be the same as length of a list of stackRaster_exposure)
#' @param keyTime name of the column to match exposure timeline from stackRaster_exposure
#' @param key name of the new column
#'
#' @export
#'
brk_exposureMatch <- function(stackRaster_exposure, sf, stackTimeline, keyTime = "TIMELINE", key = "EXPOSURE"){
# 1. Get the list of exposure of each point over the timeline of stackexposure
list_raster_exposure <- raster::as.list(stackRaster_exposure)
if(length(list_raster_exposure) != length(stackTimeline)){
stop("number of Raster in stackRaster_exposure and length of stackTimeline mismatch")
}
ls_EXPOSURE <-
lapply(1:length(list_raster_exposure),
function(i){
raster::extract(x = list_raster_exposure[[i]], y = sf)
})
# 2. get index matching between stackTimeline each individuals timeline
index_matching <- lapply(1:nrow(sf),
function(i){
match(sf[[keyTime]][[i]], stackTimeline) # ATTENTION, ORDER MATTER !!!
})
# Note: Since `NA` may appear in index_matching when timeline of individuals is longer than stackTimeline
# we have to clean timeline element where exposure is modified !
if(any(lapply(index_matching, function(x) any(is.na(x))) == TRUE)){
warning("Note that some elements in column ", keyTime, " does not match with stackTimeline. So the column ", keyTime, " has changed.")
sf[[keyTime]] = lapply(1:nrow(sf), function(i) sf[[keyTime]][[i]][which(!is.na(index_matching[[i]]))] )
}
# 3. add exposure to individual timeline
sf_EXPOSURE <-
lapply(1:nrow(sf),
function(i){
unlist(lapply(ls_EXPOSURE, `[[`, i)[index_matching[[i]]])
})
sf[[key]] <- sf_EXPOSURE
return(sf)
}
#' @title Compute exposure for a `RasterStack` class object from package `raster`
#'
#'
#' @details Local intensity depends of \code{beta} and \code{alpha} parameters.
#' Beta represents the toxic adherence between [0,1].
#' Alpha represents a list of parameters of the lost of toxic particules
#' due to covariates (precipitation).
#' There are two configurations to integrate the loss in the function :
#' (i) simulating covariate (simulate=TRUE) or (ii) uploading covariate (simulate=FALSE).
#' The covariate is linked to the loss by a linear regression with paramaters minalpha,
#' maxalpha, covariate_threshold.
#'
#' @importFrom utils flush.console
#'
#' @param RasterStack_dispersal RasterStack. An object of class\code{RasterStack} (see package \code{raster} for details).
#' @param sf sf. And object of class `sf` on which exposure is computed from the previous
#' list of raster by patch `RasterStack_dispersal`. See \link[sf]{sf} for details.
#' @param key name of the column in `sf` object providing emission amount
#' @param keyTime name o the column of sf for time
#' @param loss numeric. scalar or vector (of the same length as the number the timeline include is argument sf) to apply a loss on exposure cells.
#' @param beta numeric. toxic adherence parameter between 0 and 1 (default = 0.4).
#' @param nbr_cores integer. Set the number of cores to used for parallel computing.
#' @param quiet boolean. Set `TRUE` to remove progress bar.
#' @param unit default is meter \code{"m"}. but should be more generic: "any".
#'
#' @export
#'
brk_exposure <- function(RasterStack_dispersal,
sf,
key,
keyTime,
loss, ### to be removed once well integrated !!!
beta,
nbr_cores,
quiet,
unit){
UseMethod("brk_exposure")
}
#' @export
#'
brk_exposure.RasterStack <- function(RasterStack_dispersal,
sf,
key,
keyTime = NULL,
loss = NULL,
beta = 0.4,
nbr_cores = 1,
quiet=FALSE,
unit = "m"){
# should also check ID?
nbr_RasterStack <- length(RasterStack_dispersal@layers)
nbr_patch <- nrow(sf)
if(nbr_RasterStack != nbr_patch){
stop("Number of stack layers and sf objects differ.")
}
if(is.null(keyTime)){stop("please provide a column name of 'sf' refering to the timeline.")}
if(colListEqual(sf[[keyTime]])){
timeline <- sf[[keyTime]][[1]] # all patch have the same timeline!
} else{
stop(paste("element of column", keyTime, " in object 'sf' are not all equal."))
}
#-------------------------- BEGIN FUNCTION
# perte sur les feuilles de plante hôte : dépend du temps et des conditions météo
# fonction simulateur stochastique de climat :(attention chemin vers data de climat dans la fonction)
if(is.null(loss)){
loss = rep(0, length(timeline)) # all patch have the same timeline!
}
# ------------------------- END LOSS FUNCTION
# --- cell surface
bbox <- RasterStack_dispersal@extent
size_raster <- RasterStack_dispersal@ncols # should be equal to RasterStack_dispersal@ncols
if(unit == "m"){
cell_area <- ((bbox@xmax - bbox@xmin) * (bbox@ymax - bbox@ymin)) / ( RasterStack_dispersal@ncols * RasterStack_dispersal@nrows)
} else{
cell_area <- raster::area(RasterStack_dispersal)[1] * 10^6 # convert km^2 to m^2
}
message("Step 1/2: Compute global spatio-temporal exposure profile... ", appendLF = T)
#sf <- base::as.list(dplyr::select(sf, -patch_id))
list_matrix_convol <- lapply(1:nbr_RasterStack, function(i){ raster::as.matrix(RasterStack_dispersal[[i]]) })
out_matrix_exposure <- sptp_exposure(sf = sf,
list_matrix_convol = list_matrix_convol,
key = key,
timeline = timeline,
size_raster = size_raster,
loss = loss,
cell_area = cell_area,
beta = beta,
quiet = quiet)
message("done ", appendLF = T)
message("Step 2/2: Convert all matrices to rasters... ", appendLF = F)
RasterStack_exposure <- lapply(1:length(out_matrix_exposure), function(i){
raster::raster(x = out_matrix_exposure[[i]],
crs = RasterStack_dispersal@crs,
xmn = bbox@xmin,
xmx = bbox@xmax,
ymn = bbox@ymin,
ymx = bbox@ymax)
})
RasterStack_exposure <- raster::stack(RasterStack_exposure)
#names(RasterStack_exposure) <- paste0(timeline)
names(RasterStack_exposure) <- as.factor(timeline)
message("done ", appendLF = TRUE)
return(RasterStack_exposure)
}
###############################################################################
###############################################################################
###############################################################################
###############################################################################
############## PRIVATE ########################################################
###############################################################################
# local exposure
sptp_exposure <- function(sf,
list_matrix_convol,
key,
timeline,
size_raster,
loss,
cell_area,
beta,
quiet=FALSE){
ls_matrix = list() # list is faster, and it allows to name each matrix
size_raster <- ncol(list_matrix_convol[[1]]) # all matrix have the same size!
if(length(loss) == 1){
loss = rep(loss, length(timeline))
}
if(length(loss) != length(timeline)){
stop("length of 'loss' vector is different of the length of timeline")
}
for (time_id in 1:length(timeline)) {
matrix_toxic <- matrix(data = 0, nrow = size_raster, ncol = size_raster)
for (patch_id in 1:length(list_matrix_convol)) {
if (sf[[key]][[patch_id]][time_id] != 0) {
matrix_toxic <- matrix_toxic + (cell_area * sf[[key]][[patch_id]][time_id] * beta * list_matrix_convol[[patch_id]] )
}
}
if (time_id == 1) {
ls_matrix[[time_id]] <- matrix_toxic
} else {
ls_matrix[[time_id]] <- matrix_toxic + ls_matrix[[time_id - 1]] * (1 - loss[time_id])
}
if(quiet != TRUE){
seq_print <- 1:length(timeline)/length(timeline) * 10
# cat(paste0(time_id, "/", length(timeline), "..."))
cat('\r', "|", rep("*", round(seq_print[time_id])), rep("_", 10-round(seq_print[time_id])))
flush.console()
}
}
names(ls_matrix) <- timeline
return(ls_matrix)
}
# -------- internal
colListEqual = function(colList){
all(sapply(1:length(colList), function(i) all.equal(colList[[1]], colList[[i]])))
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.