Commit cedf3ede authored by loup.rimbaud@csiro.au's avatar loup.rimbaud@csiro.au
Browse files

graphLand.R: correction; AgriLand.R, demo_landsepi.R, graphLand.R, logit.R,...

graphLand.R: correction; AgriLand.R, demo_landsepi.R, graphLand.R, logit.R, invlogit.R: documentation
parent 393d8640
Pipeline #60 passed with stage
in 6 minutes and 4 seconds
......@@ -3,7 +3,10 @@
export(AgriLand)
export(HLIRdynamics)
export(demo_landsepi)
export(invlogit)
export(logit)
export(modeleLandsLPI)
export(multiN)
export(plot.land)
import(MASS)
import(RCALI)
......
#' @title AgriLand
#' @description do something
#' @param landscape a saptialpolygon object
#' @description Generate a landscape composed of fields where a susceptible (SC) and one (RC) or two (RC1 and RC2) resistant cultivars
#' are allocated with controlled proportions and spatio-temporal aggregation
#' @param landscape a spatialpolygon object containing field coordinates
#' @param filename output layer name
#' @param propSR num
#' @param isolSR num
#' @param propRR num
#' @param isolRR num
#' @param strat num
#' @param Nhote num
#' @param nYears number of years
#' @param Cmax0 num
#' @param Cmax1 num
#' @param propSR proportion of fields where resistance is deployed: (RC)/(SC+RC) or (RC1+RC2)/(SC+RC1+RC2)
#' @param isolSR spatial aggregation of fields where resistance is deployed (1=highly fragmented, 2=balanced, 3=highly aggregated)
#' @param propRR when applicable (mixtures and mosaics only), relative proportion of the second resistant cultivar: (RC2)/(RC1+RC2)
#' @param isolRR when applicable, spatial (for mosaics) or temporal (for rotations) aggregation of fields cultivated with the second resistant cultivar (1=highly fragmented, 2=balanced, 3=highly aggregated)
#' @param strat deployment strategy ("MO"=mosaic, "MI"=mixture, "RO"=rotations, "PY"=pyramiding)
#' @param Nhote number of cultivars (1, 2 or 3)
#' @param nYears number of simulated years
#' @param Cmax0 carrying capacity of the susceptible cultivar
#' @param Cmax1 carrying capacity of the resistant cultivars
#' @param graphOn if a graph of the landscape must be generated
#' @details An algorithm based on latent Gaussian fields is used to allocate two different crop cultivars across the simulated landscapes
#' (e.g. a susceptible and a resistant cultivar, denoted as SC and RC, respectively). This algorithm allows the control of the proportions
#' of each cultivar in terms of surface coverage, and their level of spatial aggregation. A random vector of values is drawn from a
#' multivariate normal distribution with expectation 0 and a variance-covariance matrix which depends on the pairwise distances between
#' the centroids of the fields. Next, the crop cultivars are allocated to different fields depending on whether the each value drawn from the
#' multivariate normal distribution is above or below a threshold. The proportion of each cultivar in the landscape is controlled by the value
#' of this threshold. The sequential use of this algorithm allows the allocation of more than two crop cultivars (e.g. SC, RC1 and RC2).
#' Therefore, deployment strategies involving two sources of resistance is simulated by:
#' (1) running the allocation algorithm once to segregate the fields where the susceptible cultivar is grown, and
#' (2) applying one of the following deployment strategies to the remaining candidate fields:
#' (i) Mosaics: two resistant cultivars (RC1 and RC2, carrying the first and the second resistance sources, respectively) are assigned to candidate
#' fields by re-running the allocation algorithm;
#' (ii) Mixtures: both RC1 and RC2 are allocated to all candidate fields;
#' (iii) Rotations: RC1 and RC2 are alternatively cultivated in candidate fields, depending on the number of cropping
#' seasons over which a given cultivar is grown before being rotated;
#' (iv) Pyramiding: all candidate fields are cultivated with RC12, a resistant cultivar carrying both resistance sources.
#' @importFrom sf st_as_sf
#' @importFrom sf st_write
#' @importFrom grDevices dev.off graphics.off png tiff
#' @include multiN.R periodic_cov.R graphLand.R
#' @examples
#' ## Generate a landscape consisting in a mosaic of fields cultivated with a susceptible cultivar and two resistant cultivars in balanced proportions and high level of spatial aggregation
#' data(landscapeTEST)
#' AgriLand(landscapeTEST,filename="landscapeTEST",propSR=2/3,isolSR=3,propRR=1/2,isolRR=3,strat="MO",Nhote=3,nYears=30,Cmax0=2,Cmax1=2,graphOn=1)
#' @export
AgriLand <- function(landscape,filename="landscapeTEST",propSR,isolSR,propRR,isolRR,strat,Nhote,nYears,Cmax0,Cmax1){
AgriLand <- function(landscape,filename="landscapeTEST",propSR,isolSR,propRR,isolRR,strat,Nhote,nYears,Cmax0,Cmax1,graphOn){
nPoly.tmp <- length(landscape)
prop <- c(propSR, propRR)
......@@ -82,7 +105,29 @@ habitat <- st_as_sf(landscapeINIT)
results <- st_as_sf(landscape)
st_write(habitat, paste0(filename,".gpkg"), "habitat",layer_options="OVERWRITE=yes")
st_write(results, paste0(filename,".gpkg"), "results", update = TRUE,layer_options="OVERWRITE=yes")
## Graphic representing the landscape
if (graphOn) {
title.hab <- "Simulated landscape"
isol.name <- c("low", "medium", "high", "random")
alloc.name <- c("R/(S+R)","R2/(R1+R2)")
subtitle.hab <- paste("Cropping ratio",alloc.name,"=",round(prop,2), " Aggregation",alloc.name,"=",isol.name[isol])
if (strat!="MO" | Nhote==2)
subtitle.hab <- subtitle.hab[1]
col.hab <- c("white", "gray80","gray45") ## habitat 0 = S ; habitat 1 = R1 ; habitat 2 = R2
dens.hab <- c(0,0,0)
angle.hab <- c(0,0,0)
legend.hab <- c("Susceptible","Resistant")
if ((strat=="MO" | strat=="RO" | strat=="TO") & Nhote>2)
legend.hab <- c("Susceptible","Resistant 1","Resistant 2")
if (strat=="MI")
legend.hab <- c("Susceptible","Resistant 1 + Resistant 2")
if (strat=="PY")
legend.hab <- c("Susceptible","Resistant 1+2")
png(file=paste("landscape.png",sep=""),width=1000,height=1000)
plot.land(landscape, col.hab[habitat1+1], dens.hab[habitat1+1], angle.hab[habitat1+1], col.hab, dens.hab, angle.hab, title.hab, subtitle.hab, legend.hab)
dev.off()
}
return(list(shapefilename=paste0(filename,".gpkg"),layername_hab="habitat",layername_res="results",rotation=as.integer(rotation$habitat),Cmax0=Cmax0,Cmax1=Cmax1,propRR=propRR,strat=strat))
}
......
######################################################
######### analyse.R ########
######################################################
## Computation of the output of the model from the dynamic of H, L, I, R
#source("graphLand.R")
#library("RCALI")
#' @title HLIRRdynamics
#' @description description of the function
#' @title HLIRdynamics
#' @description generate epidemiological and evolutionary outputs from simulations
#' @param pathRES path to resultats
#' @param graphOn calcul graphs data
#' @param graphOn if graphics of the output must be generated
#' @param times times
#' @param landscape a landscape
#' @param host host
......@@ -16,7 +9,7 @@
#' @param evolP evolP
#' @param th_break break
#' @param nMapPY nMapPY
#' @return files bin
#' @return binary files
#' @include RcppExports.R logit.R invlogit.R graphLand.R
#' @importFrom sf st_read
#' @importFrom grDevices colorRampPalette dev.off graphics.off gray png tiff
......@@ -41,28 +34,6 @@ HLIRdynamics <- function(pathRES, graphOn, times, landscape, hostP, epiP, evolP,
Naggr <- evolP$Naggr
adaptation <- evolP$adaptation
# mainPath <- "/home/loup/Documents/model_resistance"
# pathRES <- "/home/loup/Documents/model_resistance/RES/scenarios"
# pathOUT <- "/home/loup/Documents/model_resistance/RES"
# graphOn <- 1
# nYears <- 5
# nTSpY <- 120
# Nhote <- 3
# NpolyTarg <- 150
# idLAN <- 1
# Cmax0 <- 2.0
# Cmax1 <- 2.0
# Naggr <- 6
# th_break <- 50000
# adaptation <- "10000100"
# resistance0 <- "00000000"
# resistance1 <- "10000000"
# resistance2 <- "01000000"
# nMapPY <- 0
# propRR <- 4/6
# strat <- "MO"
#### LANDSCAPE ####
area <- paysage$area
......
#library(RCALI) ## function readpoly2
#library(fields) ## function Exponential
#library(MASS) ## function mvnorm
#library(Matrix) ## function nearPD to convert a covariance matrix in a positive-definite matrix
#library(splancs) ## function polymap
#library(rgdal)
#library(raster)
#library(sp)
......@@ -10,27 +10,31 @@
#library(maptools)
#library(fields)
#library(sf)
#
#' Run demo landepi
#' @title demo_landsepi
#' @description Run a simulated example of mosaic deployment strategy of two resistant cultivars in balanced proportions and high level of spatial aggregation.
#' @include RcppExports.R AgriLand.R graphLand.R multiN.R periodic_cov.R
#' @importFrom utils data
#' @export
demo_landsepi <- function(){
#outputs
pathRES <- getwd()
graphOn <- 1
#landscape parameters
propSR=1/2 ## proportion of cultivars >0
propSR=2/3 ## proportion of cultivars >0
isolSR=3 ## Class of aggregation between S and R (0:3 with 1=low, 2=medium, 3=high, 4=random aggregation)
propRR=4/6 ## relative proportion of cultivars >1
propRR=1/2 ## relative proportion of cultivars >1
isolRR=3 ## Class of aggregation between R1 and R2 (0:3 with 1=low, 2=medium, 3=high, 4=random aggregation)
strat <- "MO"
Nhote <- 3
nYears <- 5
nYears <- 30
Cmax0 <- 2
Cmax1 <- 2
data(landscapeTEST)
landscape <- AgriLand(landscapeTEST,filename="paysageTEST",propSR,isolSR,propRR,isolRR,strat,Nhote,nYears,Cmax0,Cmax1)
landscape <- AgriLand(landscapeTEST,filename="paysageTEST",propSR,isolSR,propRR,isolRR,strat,Nhote,nYears,Cmax0,Cmax1,graphOn)
#time parameters
nTSpY = 120
......@@ -99,9 +103,7 @@ demo_landsepi <- function(){
epiP=paramepi,
evolP=paramevol)
#outputs
pathRES <- getwd()
graphOn <- 1
times <- paramT
hostP <- paramH
epiP <- paramepi
......
#### FUNCTIONS TO ESTIMATE INTERVALS WITH DIFFERENT COLORS ####
#### FUNCTION TO DRAW THE LANDSCAPE ####
#' @title plot
#' @description plot
#' @param landscape a landscape
#' @param COL colors
#' @param DENS density
#' @param ANGLE angle
#' @param COL.LEG leg
#' @param DENS.LEG den
#' @param ANGLE.LEG angle
#' @param TITLE title
#' @param SUBTITLE subtitle
#' @param LEGEND1 legend
#' @param LEGEND2 legend2
#' @param TITLE.LEG2 title leg 2
#' @param XMAX xmax
#' @param YMAX YMAX
#' @title plot.land
#' @description plot a landscape with colors or hatched lines on fields
#' @param landscape a spatialpolygon object containing field coordinates
#' @param COL vector containing the color of each field
#' @param DENS vector containing the density of hatched lines for each field
#' @param ANGLE vector containing the angle of hatched lines for each field
#' @param COL.LEG vector containing the colors in the first legend
#' @param DENS.LEG vector containing the density of hatched lines in the second legend
#' @param ANGLE.LEG vector containing the angle of hatched lines in the second legend
#' @param TITLE title of the graphic
#' @param SUBTITLE subtitle of the graphic
#' @param LEGEND1 labels in the first legend (colors)
#' @param LEGEND2 labels in the second legend (hatched lines)
#' @param TITLE.LEG2 title for the second legend
#' @param XMAX maximal coordinate on horizontal axis
#' @param YMAX maximal coordinate on vertical axis
#' @examples
#' ## Draw a landscape with various colours
#' data(landscapeTEST)
#' plot.land(landscapeTEST, COL=1:length(landscapeTEST), DENS=rep(0,length(landscapeTEST)), ANGLE=rep(30,length(landscapeTEST)))
#' @include RcppExports.R logit.R invlogit.R
# @S3method plot land
#' @export
......@@ -30,13 +31,15 @@ plot.land <- function(landscape, COL=rep(0,length(landscape)), DENS=rep(0,length
if (length(SUBTITLE)>1)
mtext(SUBTITLE[2], side=3, line=0, padj=1.5, las=1, cex=1.4)
for (i in 1:nPoly) {
polymap(landscape[[i]], add=TRUE, col=COL[i], border="black", lwd=2.5)
polymap(landscape[[i]], add=TRUE, col="black", density=DENS[i], angle=ANGLE[i], border=NA)
polymap(landscape@polygons[[i]]@Polygons[[1]]@coords, add=TRUE, col=COL[i], border="black", lwd=2.5)
polymap(landscape@polygons[[i]]@Polygons[[1]]@coords, add=TRUE, col="black", density=DENS[i], angle=ANGLE[i], border=NA)
}
if (TITLE.LEG2=="")
legend(XMAX/2.66, -YMAX/40, legend=LEGEND1, fill=COL.LEG, bty="n")
else {
legend(XMAX/2.66, -YMAX/40, legend=LEGEND1, col="black", density=2*DENS.LEG, angle=ANGLE.LEG, bty="n")
legend(-XMAX/5, YMAX, legend=LEGEND2, fill=COL.LEG, bty="n", title=TITLE.LEG2)
if (LEGEND1[1]!=""){
if (TITLE.LEG2=="")
legend(XMAX/2.66, -YMAX/40, legend=LEGEND1, fill=COL.LEG, bty="n")
else {
legend(XMAX/2.66, -YMAX/40, legend=LEGEND1, col="black", density=2*DENS.LEG, angle=ANGLE.LEG, bty="n")
legend(-XMAX/5, YMAX, legend=LEGEND2, fill=COL.LEG, bty="n", title=TITLE.LEG2)
}
}
}
invlogit <- function(n) exp(n) / (1+exp(n))
#' @title invlogit
#' @description Given a numeric object return the invlogit of the values. Missing values (NAs) are allowed.
#' @param x a numeric object
#' @details The invlogit is defined by exp(x) / (1+exp(x)). Values in x of -Inf or Inf return invlogits of 0 or 1 respectively. Any NAs in the input will also be NAs in the output.
#' @return An object of the same type as x containing the invlogits of the input values.
#' @export
invlogit <- function(x) exp(x) / (1+exp(x))
logit <- function(p) log(p/(1-p))
#' @title logit
#' @description Given a numeric object return the logit of the values. Missing values (NAs) are allowed.
#' @param x a numeric object containing values between 0 and 1
#' @details The logit is defined by log(x/(1-x)). Values in x of 0 or 1 return logits of -Inf or Inf respectively. Any NAs in the input will also be NAs in the output.
#' @return An object of the same type as x containing the logits of the input values.
#' @export
logit <- function(x) log(x/(1-x))
#' @title multiN
#' @description Use a multivariate normal distribution to An algorithm based on latent Gaussian fields is used to allocate two different crop cultivars across the simulated landscapes Given a numeric object return the invlogit of the values. Missing values (NAs) are allowed.
#' @param x a numeric object
#' @details The invlogit is defined by exp(x) / (1+exp(x)). Values in x of -Inf or Inf return invlogits of 0 or 1 respectively. Any NAs in the input will also be NAs in the output.
#' @return An object of the same type as x containing the invlogits of the input values.
#' @export
#'
## take the matrix d to compute the allocation of cultivars
multiN <- function(d, area, aggreg, prop) {
nPoly <- nrow(area)
......
......@@ -5,31 +5,57 @@
\title{AgriLand}
\usage{
AgriLand(landscape, filename = "landscapeTEST", propSR, isolSR, propRR,
isolRR, strat, Nhote, nYears, Cmax0, Cmax1)
isolRR, strat, Nhote, nYears, Cmax0, Cmax1, graphOn)
}
\arguments{
\item{landscape}{a saptialpolygon object}
\item{landscape}{a spatialpolygon object containing field coordinates}
\item{filename}{output layer name}
\item{propSR}{num}
\item{propSR}{proportion of fields where resistance is deployed: (RC)/(SC+RC) or (RC1+RC2)/(SC+RC1+RC2)}
\item{isolSR}{num}
\item{isolSR}{spatial aggregation of fields where resistance is deployed (1=highly fragmented, 2=balanced, 3=highly aggregated)}
\item{propRR}{num}
\item{propRR}{when applicable (mixtures and mosaics only), relative proportion of the second resistant cultivar: (RC2)/(RC1+RC2)}
\item{isolRR}{num}
\item{isolRR}{when applicable, spatial (for mosaics) or temporal (for rotations) aggregation of fields cultivated with the second resistant cultivar (1=highly fragmented, 2=balanced, 3=highly aggregated)}
\item{strat}{num}
\item{strat}{deployment strategy ("MO"=mosaic, "MI"=mixture, "RO"=rotations, "PY"=pyramiding)}
\item{Nhote}{num}
\item{Nhote}{number of cultivars (1, 2 or 3)}
\item{nYears}{number of years}
\item{nYears}{number of simulated years}
\item{Cmax0}{num}
\item{Cmax0}{carrying capacity of the susceptible cultivar}
\item{Cmax1}{num}
\item{Cmax1}{carrying capacity of the resistant cultivars}
\item{graphOn}{if a graph of the landscape must be generated}
}
\description{
do something
Generate a landscape composed of fields where a susceptible (SC) and one (RC) or two (RC1 and RC2) resistant cultivars
are allocated with controlled proportions and spatio-temporal aggregation
}
\details{
An algorithm based on latent Gaussian fields is used to allocate two different crop cultivars across the simulated landscapes
(e.g. a susceptible and a resistant cultivar, denoted as SC and RC, respectively). This algorithm allows the control of the proportions
of each cultivar in terms of surface coverage, and their level of spatial aggregation. A random vector of values is drawn from a
multivariate normal distribution with expectation 0 and a variance-covariance matrix which depends on the pairwise distances between
the centroids of the fields. Next, the crop cultivars are allocated to different fields depending on whether the each value drawn from the
multivariate normal distribution is above or below a threshold. The proportion of each cultivar in the landscape is controlled by the value
of this threshold. The sequential use of this algorithm allows the allocation of more than two crop cultivars (e.g. SC, RC1 and RC2).
Therefore, deployment strategies involving two sources of resistance is simulated by:
(1) running the allocation algorithm once to segregate the fields where the susceptible cultivar is grown, and
(2) applying one of the following deployment strategies to the remaining candidate fields:
(i) Mosaics: two resistant cultivars (RC1 and RC2, carrying the first and the second resistance sources, respectively) are assigned to candidate
fields by re-running the allocation algorithm;
(ii) Mixtures: both RC1 and RC2 are allocated to all candidate fields;
(iii) Rotations: RC1 and RC2 are alternatively cultivated in candidate fields, depending on the number of cropping
seasons over which a given cultivar is grown before being rotated;
(iv) Pyramiding: all candidate fields are cultivated with RC12, a resistant cultivar carrying both resistance sources.
}
\examples{
## Generate a landscape consisting in a mosaic of fields cultivated with a susceptible cultivar and two resistant cultivars in balanced proportions and high level of spatial aggregation
data(landscapeTEST)
AgriLand(landscapeTEST,filename="landscapeTEST",propSR=2/3,isolSR=3,propRR=1/2,isolRR=3,strat="MO",Nhote=3,nYears=30,Cmax0=2,Cmax1=2,graphOn=1)
}
......@@ -2,10 +2,13 @@
% Please edit documentation in R/demo_landsepi.R
\name{demo_landsepi}
\alias{demo_landsepi}
\title{Run demo landepi}
\title{demo_landsepi}
\usage{
demo_landsepi()
}
\description{
Run a simulated example of mosaic deployment strategy of two resistant cultivars in balanced proportions and high level of spatial aggregation.
}
\details{
Run demo landepi
}
......@@ -2,7 +2,7 @@
% Please edit documentation in R/graphLand.R
\name{plot.land}
\alias{plot.land}
\title{plot}
\title{plot.land}
\usage{
plot.land(landscape, COL = rep(0, length(landscape)), DENS = rep(0,
length(landscape)), ANGLE = rep(30, length(landscape)),
......@@ -12,34 +12,39 @@ plot.land(landscape, COL = rep(0, length(landscape)), DENS = rep(0,
TITLE.LEG2 = "", XMAX = 2000, YMAX = 2000)
}
\arguments{
\item{landscape}{a landscape}
\item{landscape}{a spatialpolygon object containing field coordinates}
\item{COL}{colors}
\item{COL}{vector containing the color of each field}
\item{DENS}{density}
\item{DENS}{vector containing the density of hatched lines for each field}
\item{ANGLE}{angle}
\item{ANGLE}{vector containing the angle of hatched lines for each field}
\item{COL.LEG}{leg}
\item{COL.LEG}{vector containing the colors in the first legend}
\item{DENS.LEG}{den}
\item{DENS.LEG}{vector containing the density of hatched lines in the second legend}
\item{ANGLE.LEG}{angle}
\item{ANGLE.LEG}{vector containing the angle of hatched lines in the second legend}
\item{TITLE}{title}
\item{TITLE}{title of the graphic}
\item{SUBTITLE}{subtitle}
\item{SUBTITLE}{subtitle of the graphic}
\item{LEGEND1}{legend}
\item{LEGEND1}{labels in the first legend (colors)}
\item{LEGEND2}{legend2}
\item{LEGEND2}{labels in the second legend (hatched lines)}
\item{TITLE.LEG2}{title leg 2}
\item{TITLE.LEG2}{title for the second legend}
\item{XMAX}{xmax}
\item{XMAX}{maximal coordinate on horizontal axis}
\item{YMAX}{YMAX}
\item{YMAX}{maximal coordinate on vertical axis}
}
\description{
plot
plot a landscape with colors or hatched lines on fields
}
\examples{
## Draw a landscape with various colours
data(landscapeTEST)
plot.land(landscapeTEST, COL=1:length(landscapeTEST), DENS=rep(0,length(landscapeTEST)), ANGLE=rep(30,length(landscapeTEST)))
}
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment