R/POST_FATE.graphic_evolutionStability.R
In leca-dev/RFate: FATE with R
Defines functions
POST_FATE.graphic_evolutionStability
Documented in
POST_FATE.graphic_evolutionStability
### HEADER #####################################################################
##' @title Create a graphical representation of the evolution of habitat
##' composition through time for a \code{FATE} simulation
##'
##' @name POST_FATE.graphic_evolutionStability
##'
##' @author Ceres Barros, Maya Guéguen
##'
##' @description This script is designed to produce one graphical representation
##' for a \code{FATE} simulation : the evolution through time of the total
##' abundance and evenness of each habitat.
##'
##'
##' @param name.simulation a \code{string} corresponding to the main directory
##' or simulation name of the \code{FATE} simulation
##' @param file.simulParam default \code{NULL}. \cr A \code{string}
##' corresponding to the name of a parameter file that will be contained into
##' the \code{PARAM_SIMUL} folder of the \code{FATE} simulation
##' @param movingWindow_size default \code{3}. \cr An \code{integer}
##' corresponding to the size (\emph{in years}) of the moving window that will
##' be used to calculate metrics of habitat stability
##' @param movingWindow_step default \code{1}. \cr An \code{integer}
##' corresponding to the step (\emph{in years}) between the years of the moving
##' window
##' @param opt.doPlot (\emph{optional}) default \code{TRUE}. \cr If \code{TRUE},
##' plot(s) will be processed, otherwise only the calculation and reorganization
##' of outputs will occur, be saved and returned
##'
##' @details
##'
##' This function allows to obtain, for a specific \code{FATE} simulation and
##' a specific parameter file within this simulation, one preanalytical graphic :
##'
##' \itemize{
##' \item{the evolution of \strong{total abundance} (\code{FATE}
##' \emph{arbitrary unit}) and \strong{evenness} (\emph{between \code{0} and
##' \code{1}}) of each habitat through simulation time, with \emph{evenness}
##' representing the uniformity of the species composition of the habitat
##' (similar to \strong{Shannon entropy}) :
##' \deqn{
##' \text{evenness} = - \frac{\Sigma(\text{proportion}_{\text{ PFG}_i} *
##' log(\text{proportion}_{\text{ PFG}_i}))}{log(\text{no.PFG})}
##' }
##' with \deqn{
##' \text{proportion}_{\text{ PFG}_i} = \frac{abund_{\text{ PFG}_i
##' \text{, }\text{Habitat}_j}}{abund_{\text{ PFG}_{all}\text{, }
##' \text{Habitat}_j}}
##' }
##' }
##' }
##'
##' If the information has been provided (see
##' \code{\link{POST_FATE.temporalEvolution}}), the graphics will be also done
##' per habitat. \cr \cr
##'
##' \strong{It requires} that the \code{\link{POST_FATE.temporalEvolution}}
##' function has been run and that the file
##' \code{POST_FATE_TABLE_PIXEL_evolution_abundance.csv} exists.
##'
##'
##'
##' @return A \code{list} containing two \code{data.frame} objects with the
##' following columns, and one \code{ggplot2} object :
##'
##' \describe{
##' \item{tab.hab}{
##' \describe{
##' \item{\code{HAB}}{concerned habitat}
##' \item{\code{year}}{concerned simulation year}
##' \item{\code{totalAbundance}}{total abundance over all the pixels
##' within the concerned habitat}
##' \item{\code{no.PFG}}{number of PFG over all the pixels within the
##' concerned habitat}
##' \item{\code{evenness}}{evenness over all the pixels within the
##' concerned habitat}
##' }
##' }
##' \item{tab.stab}{
##' \describe{
##' \item{\code{HAB}}{concerned habitat}
##' \item{\code{no.years}}{number of simulation years used (\emph{moving
##' window size})}
##' \item{\code{yearStep}}{step between each simulation year of the moving
##' window}
##' \item{\code{yearStart}}{first simulation year of the moving window}
##' \item{\code{yearEnd}}{last simulation year of the moving window}
##' \item{\code{metric}}{concerned metric (either \code{totalAbundance} or
##' \code{evenness})}
##' \item{\code{mean}}{mean value of the concerned metric over the years
##' of the concerned moving window}
##' \item{\code{sd}}{value of standard deviation of the concerned metric
##' over the years of the concerned moving window}
##' \item{\code{cv}}{value of coefficient of variation of the concerned
##' metric over the years of the concerned moving window}
##' }
##' }
##' \item{plot.stab}{\code{ggplot2} object, representing the evolution of
##' total abundance and evenness of each habitat \cr \cr}
##' }
##'
##' Two \code{POST_FATE_TABLE_HAB_evolution_[...].csv} files are created :
##' \describe{
##' \item{\file{stability1}}{always, containing \code{tab.hab}}
##' \item{\file{stability2}}{\emph{if successive years available}, containing
##' \code{tab.stab}}
##' }
##'
##'
##' One \code{POST_FATE_[...].pdf} files is created :
##' \describe{
##' \item{\file{GRAPHIC_A \cr stability}}{to visualize for each habitat the
##' evolution of its total abundance and its evenness through simulation time}
##' }
##'
##'
##' @keywords FATE, outputs, abundance through time, evenness
##'
##' @seealso \code{\link{POST_FATE.temporalEvolution}}
##'
##'
##' @export
##'
##' @importFrom utils write.csv
##' @importFrom data.table fread
##' @importFrom foreach foreach %do%
##' @importFrom reshape2 melt
##' @importFrom grDevices colorRampPalette
##'
##' @importFrom ggplot2 ggplot ggsave aes_string
##' geom_line geom_point geom_rect
##' scale_color_manual scale_fill_manual
##' facet_grid labs theme element_text element_blank
##' @importFrom ggthemes theme_fivethirtyeight
##'
## END OF HEADER ###############################################################
POST_FATE.graphic_evolutionStability = function(
name.simulation
, file.simulParam = NULL
, movingWindow_size = 3
, movingWindow_step = 1
, opt.doPlot = TRUE
){
#############################################################################
## CHECK parameter name.simulation
.testParam_existFolder(name.simulation, "PARAM_SIMUL/")
.testParam_existFolder(name.simulation, "RESULTS/")
.testParam_existFolder(name.simulation, "DATA/")
name.simulation = sub("/", "", name.simulation)
## CHECK parameter file.simulParam
abs.simulParams = .getParam_abs.simulParams(file.simulParam, name.simulation)
## CHECK parameter movingWindow_size
.testParam_notInteger.m("movingWindow_size", movingWindow_size)
## CHECK parameter movingWindow_step
.testParam_notInteger.m("movingWindow_step", movingWindow_step)
cat("\n\n #------------------------------------------------------------#")
cat("\n # POST_FATE.graphic_evolutionStability")
cat("\n #------------------------------------------------------------# \n")
#############################################################################
res = foreach (abs.simulParam = abs.simulParams) %do%
{
cat("\n+++++++\n")
cat("\n Simulation name : ", name.simulation)
cat("\n Simulation file : ", abs.simulParam)
cat("\n")
## Get results directories ------------------------------------------------
GLOB_DIR = .getGraphics_results(name.simulation = name.simulation
, abs.simulParam = abs.simulParam)
## Get number of PFGs -----------------------------------------------------
## Get PFG names ----------------------------------------------------------
GLOB_SIM = .getGraphics_PFG(name.simulation = name.simulation
, abs.simulParam = abs.simulParam)
## Get raster mask --------------------------------------------------------
GLOB_MASK = .getGraphics_mask(name.simulation = name.simulation
, abs.simulParam = abs.simulParam)
## Get the abundance table ------------------------------------------------
file.abundance = paste0(name.simulation
, "/RESULTS/POST_FATE_TABLE_PIXEL_evolution_abundance_"
, basename(GLOB_DIR$dir.save)
, ".csv")
.testParam_existFile(file.abundance)
tab.totalAbundance = fread(file.abundance)
tab.totalAbundance = as.data.frame(tab.totalAbundance, stringsAsFactors = FALSE)
years = colnames(tab.totalAbundance)
years = years[which(!(years %in% c("PFG", "ID.pixel", "X", "Y", "HAB")))]
years = as.numeric(years)
no_years = length(years)
hab_names = sort(unique(tab.totalAbundance$HAB))
no_hab = length(hab_names)
cat("\n Number of years : ", no_years)
cat("\n Number of habitat : ", no_hab)
cat("\n")
###########################################################################
## Transform the data inside the table ------------------------------------
cat("\n ---------- GETTING TOTAL ABUNDANCE and EVENNESS within each habitat... \n")
## Getting total abundance within each habitat
cat("\n> Habitat total abundance...")
tab.split = split(tab.totalAbundance, list(tab.totalAbundance$HAB))
tab.hab = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
hab = names(tab.split)[i]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
tab = tab[, as.character(years), drop = FALSE]
return(data.frame(HAB = hab, year = years
, totalAbundance = colSums(tab, na.rm = TRUE)
, stringsAsFactors = FALSE))
}
}
## Getting PFG relative abundance within each habitat
cat("\n> PFG relative abundance...")
tab.split = split(tab.totalAbundance, list(tab.totalAbundance$PFG
, tab.totalAbundance$HAB))
tab.prop = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
pfg = strsplit(names(tab.split)[i], "[.]")[[1]][1]
hab = strsplit(names(tab.split)[i], "[.]")[[1]][2]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
tab = tab[, as.character(years), drop = FALSE]
tot.pfg = colSums(tab, na.rm = TRUE)
tot.hab = tab.hab[which(tab.hab$HAB == hab), ]
rownames(tot.hab) = tot.hab$year
prop = tot.pfg / tot.hab[as.character(years), "totalAbundance"]
return(data.frame(PFG = pfg, HAB = hab, year = years
, PROP = prop, stringsAsFactors = FALSE))
}
}
## Calculating evenness within each habitat
cat("\n> Habitat evenness...")
tab.split = split(tab.prop, list(tab.prop$HAB, tab.prop$year))
tab.even = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
hab = strsplit(names(tab.split)[i], "[.]")[[1]][1]
ye = strsplit(names(tab.split)[i], "[.]")[[1]][2]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
even = (- sum(tab$PROP * log(tab$PROP)) / log(nrow(tab)))
return(data.frame(HAB = hab, year = ye, no.PFG = nrow(tab)
, evenness = even, stringsAsFactors = FALSE))
}
}
## Merging total abundance and evenness for each habitat
if (nrow(tab.even) > 0) {
tab.HAB = merge(tab.hab, tab.even, by = c("HAB", "year"))
} else {
tab.HAB = tab.hab
}
cat("\n")
write.csv(tab.HAB
, file = paste0(name.simulation
, "/RESULTS/POST_FATE_TABLE_HAB_evolution_stability1_"
, basename(GLOB_DIR$dir.save)
, ".csv")
, row.names = FALSE)
message(paste0("\n The output file \n"
, " > POST_FATE_TABLE_HAB_evolution_stability1_"
, basename(GLOB_DIR$dir.save)
, ".csv \n"
, "has been successfully created !\n"))
###########################################################################
cat("\n ---------- EVALUATE, if possible, ABUNDANCE and EVENNESS STABILITY within each habitat... \n")
.getSucc = function(vec.years)
{
if (length(vec.years) > 1)
{
## gaps between years
gaps = vec.years[2:length(vec.years)] - vec.years[1:(length(vec.years) - 1)]
## periods with successive years
gaps = ifelse(gaps == movingWindow_step, 1, 0)
succ = c(gaps[1], rep(0, length(gaps) - 1))
for (i in 2:length(gaps))
{
if (gaps[i] == 1){
succ[i] = succ[i - 1] + 1
} else {
succ[i] = 0
}
}
sel = vec.years[which(succ > 0) + 1]
if (succ[1] == 1) sel = c(vec.years[1], sel)
return(list(selected = sel, max = max(succ, na.rm = TRUE) + 1))
} else
{
return(list(selected = vec.years, max = length(vec.years)))
}
}
## Calculate number of successive years -----------------------------------
max_years.succ = .getSucc(years)$max
years.succ = .getSucc(years)$selected
no_years.succ = length(years.succ)
no_movingWindow = no_years.succ - movingWindow_size + 1
## Evaluate if stability can be done --------------------------------------
if (length(years.succ) < movingWindow_size ||
movingWindow_size > max_years.succ ||
no_movingWindow < 1)
{
tab.STAB = pp = NULL
warning(paste0("Missing data!\n No moving window can be defined "
, "with the current parameters and available data :\n"
, " > mw_size = ", movingWindow_size, "\n"
, " > number of successive years = ", no_years.succ, "\n"
, " > maximum number of successive years = ", max_years.succ
, "\n Please check."))
} else
{
tab.split = split(tab.HAB, list(tab.HAB$HAB))
tab.STAB = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
hab = names(tab.split)[i]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
combi = expand.grid(mw = 1:no_movingWindow, metric = c("totalAbundance", "evenness")
, stringsAsFactors = FALSE)
tab.mw = foreach(mw = combi$mw
, metric = combi$metric
, .combine = "rbind"
) %do%
{
mw_years = years[mw:(mw + movingWindow_size - 1)]
mw_years = .getSucc(mw_years)$selected
tab.mw = tab[which(tab$year %in% c(mw_years)), ]
if (nrow(tab.mw) > 0)
{
val = sum(tab.mw[, metric])
val.sq = sum(tab.mw[, metric] * tab.mw[, metric])
val.mean = val / movingWindow_size
val.sd = sqrt(val.sq / movingWindow_size - val.mean * val.mean)
val.cv = val.sd / val.mean
return(data.frame(HAB = hab
, no.years = length(mw_years)
, yearStep = movingWindow_step
, yearStart = min(mw_years)
, yearEnd = max(mw_years)
, metric = metric
, mean = val.mean
, sd = val.sd
, cv = val.cv
, stringsAsFactors = FALSE))
}
}
if (nrow(tab.mw) > 0)
{
return(tab.mw)
}
}
}
write.csv(tab.STAB
, file = paste0(name.simulation
, "/RESULTS/POST_FATE_TABLE_HAB_evolution_stability2_"
, basename(GLOB_DIR$dir.save)
, ".csv")
, row.names = FALSE)
message(paste0("\n The output file \n"
, " > POST_FATE_TABLE_HAB_evolution_stability2_"
, basename(GLOB_DIR$dir.save)
, ".csv \n"
, "has been successfully created !\n"))
## Compare quantiles ----------------------------------------------------
## TO BE DONE ??
## produce the plot -----------------------------------------------------
if (opt.doPlot && !is.null(tab.HAB))
{
cat("\n ---------- PRODUCING PLOT \n")
col_vec = c("#1B9E77", "#D95F02", "#7570B3", "#E7298A", "#66A61E", "#E6AB02", "#A6761D", "#666666")
col_fun = colorRampPalette(col_vec)
## Evolution of abundance and evenness through time -------------------
tab.plot1 = melt(tab.HAB, id.vars = c("HAB", "year"))
colnames(tab.plot1) = c("HAB", "year", "metric", "value")
tab.plot1$metric = factor(tab.plot1$metric, c("totalAbundance", "evenness", "no.PFG"))
## plot
pp = ggplot(tab.plot1, aes_string(x = "year", y = "value", color = "factor(HAB, hab_names)"))
if (!is.null(tab.STAB))
{
## Evolution of stability through time --------------------------------
tab.plot2 = tab.STAB
tab.plot2$year_median = sapply(1:nrow(tab.plot2), function(x) {
median(as.numeric(as.character(tab.plot2[x, c("yearStart", "yearEnd")]))) })
tab.plot2 = tab.plot2[which(tab.plot2$sd > 0.00001), ]
## plot
pp = pp +
geom_rect(data = tab.plot2, inherit.aes = FALSE, alpha = 0.5
, aes_string(xmin = "yearStart", xmax = "yearEnd"
, ymin = "mean - sd"
, ymax = "mean + sd"
, fill = "factor(HAB, hab_names)"))
}
## plot
pp = pp +
geom_line(lwd = 1) +
geom_point() +
facet_grid("metric ~ .", scales = "free_y"
, labeller = as_labeller(c("totalAbundance" = "Total abundance"
, "evenness" = "Evenness"
, "no.PFG" = "Number of PFG"))) +
scale_color_manual("Habitat", values = col_fun(no_hab)) +
scale_fill_manual(guide = "none", values = col_fun(no_hab)) +
labs(x = "", y = ""
, title = paste0("GRAPH A : evolution of habitat composition")) +
.getGraphics_theme()
ggsave(filename = paste0(name.simulation
, "/RESULTS/POST_FATE_GRAPHIC_A_evolution_stability_"
, basename(GLOB_DIR$dir.save), ".pdf")
, plot = pp, width = 10, height = 8)
} else
{
pp = NULL
} ## END opt.doPlot
}
## ------------------------------------------------------------------------
cat("\n> Done!\n")
return(list(tab.hab = tab.HAB
, tab.stab = tab.STAB
, plot.stab = pp))
} ## END loop on abs.simulParams
names(res) = abs.simulParams
return(res)
}
leca-dev/RFate documentation built on Sept. 19, 2024, 6:09 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions POST_FATE.graphic_evolutionStability
Documented in POST_FATE.graphic_evolutionStability
### HEADER #####################################################################
##' @title Create a graphical representation of the evolution of habitat
##' composition through time for a \code{FATE} simulation
##'
##' @name POST_FATE.graphic_evolutionStability
##'
##' @author Ceres Barros, Maya Guéguen
##'
##' @description This script is designed to produce one graphical representation
##' for a \code{FATE} simulation : the evolution through time of the total
##' abundance and evenness of each habitat.
##'
##'
##' @param name.simulation a \code{string} corresponding to the main directory
##' or simulation name of the \code{FATE} simulation
##' @param file.simulParam default \code{NULL}. \cr A \code{string}
##' corresponding to the name of a parameter file that will be contained into
##' the \code{PARAM_SIMUL} folder of the \code{FATE} simulation
##' @param movingWindow_size default \code{3}. \cr An \code{integer}
##' corresponding to the size (\emph{in years}) of the moving window that will
##' be used to calculate metrics of habitat stability
##' @param movingWindow_step default \code{1}. \cr An \code{integer}
##' corresponding to the step (\emph{in years}) between the years of the moving
##' window
##' @param opt.doPlot (\emph{optional}) default \code{TRUE}. \cr If \code{TRUE},
##' plot(s) will be processed, otherwise only the calculation and reorganization
##' of outputs will occur, be saved and returned
##'
##' @details
##'
##' This function allows to obtain, for a specific \code{FATE} simulation and
##' a specific parameter file within this simulation, one preanalytical graphic :
##'
##' \itemize{
##' \item{the evolution of \strong{total abundance} (\code{FATE}
##' \emph{arbitrary unit}) and \strong{evenness} (\emph{between \code{0} and
##' \code{1}}) of each habitat through simulation time, with \emph{evenness}
##' representing the uniformity of the species composition of the habitat
##' (similar to \strong{Shannon entropy}) :
##' \deqn{
##' \text{evenness} = - \frac{\Sigma(\text{proportion}_{\text{ PFG}_i} *
##' log(\text{proportion}_{\text{ PFG}_i}))}{log(\text{no.PFG})}
##' }
##' with \deqn{
##' \text{proportion}_{\text{ PFG}_i} = \frac{abund_{\text{ PFG}_i
##' \text{, }\text{Habitat}_j}}{abund_{\text{ PFG}_{all}\text{, }
##' \text{Habitat}_j}}
##' }
##' }
##' }
##'
##' If the information has been provided (see
##' \code{\link{POST_FATE.temporalEvolution}}), the graphics will be also done
##' per habitat. \cr \cr
##'
##' \strong{It requires} that the \code{\link{POST_FATE.temporalEvolution}}
##' function has been run and that the file
##' \code{POST_FATE_TABLE_PIXEL_evolution_abundance.csv} exists.
##'
##'
##'
##' @return A \code{list} containing two \code{data.frame} objects with the
##' following columns, and one \code{ggplot2} object :
##'
##' \describe{
##' \item{tab.hab}{
##' \describe{
##' \item{\code{HAB}}{concerned habitat}
##' \item{\code{year}}{concerned simulation year}
##' \item{\code{totalAbundance}}{total abundance over all the pixels
##' within the concerned habitat}
##' \item{\code{no.PFG}}{number of PFG over all the pixels within the
##' concerned habitat}
##' \item{\code{evenness}}{evenness over all the pixels within the
##' concerned habitat}
##' }
##' }
##' \item{tab.stab}{
##' \describe{
##' \item{\code{HAB}}{concerned habitat}
##' \item{\code{no.years}}{number of simulation years used (\emph{moving
##' window size})}
##' \item{\code{yearStep}}{step between each simulation year of the moving
##' window}
##' \item{\code{yearStart}}{first simulation year of the moving window}
##' \item{\code{yearEnd}}{last simulation year of the moving window}
##' \item{\code{metric}}{concerned metric (either \code{totalAbundance} or
##' \code{evenness})}
##' \item{\code{mean}}{mean value of the concerned metric over the years
##' of the concerned moving window}
##' \item{\code{sd}}{value of standard deviation of the concerned metric
##' over the years of the concerned moving window}
##' \item{\code{cv}}{value of coefficient of variation of the concerned
##' metric over the years of the concerned moving window}
##' }
##' }
##' \item{plot.stab}{\code{ggplot2} object, representing the evolution of
##' total abundance and evenness of each habitat \cr \cr}
##' }
##'
##' Two \code{POST_FATE_TABLE_HAB_evolution_[...].csv} files are created :
##' \describe{
##' \item{\file{stability1}}{always, containing \code{tab.hab}}
##' \item{\file{stability2}}{\emph{if successive years available}, containing
##' \code{tab.stab}}
##' }
##'
##'
##' One \code{POST_FATE_[...].pdf} files is created :
##' \describe{
##' \item{\file{GRAPHIC_A \cr stability}}{to visualize for each habitat the
##' evolution of its total abundance and its evenness through simulation time}
##' }
##'
##'
##' @keywords FATE, outputs, abundance through time, evenness
##'
##' @seealso \code{\link{POST_FATE.temporalEvolution}}
##'
##'
##' @export
##'
##' @importFrom utils write.csv
##' @importFrom data.table fread
##' @importFrom foreach foreach %do%
##' @importFrom reshape2 melt
##' @importFrom grDevices colorRampPalette
##'
##' @importFrom ggplot2 ggplot ggsave aes_string
##' geom_line geom_point geom_rect
##' scale_color_manual scale_fill_manual
##' facet_grid labs theme element_text element_blank
##' @importFrom ggthemes theme_fivethirtyeight
##'
## END OF HEADER ###############################################################
POST_FATE.graphic_evolutionStability = function(
name.simulation
, file.simulParam = NULL
, movingWindow_size = 3
, movingWindow_step = 1
, opt.doPlot = TRUE
){
#############################################################################
## CHECK parameter name.simulation
.testParam_existFolder(name.simulation, "PARAM_SIMUL/")
.testParam_existFolder(name.simulation, "RESULTS/")
.testParam_existFolder(name.simulation, "DATA/")
name.simulation = sub("/", "", name.simulation)
## CHECK parameter file.simulParam
abs.simulParams = .getParam_abs.simulParams(file.simulParam, name.simulation)
## CHECK parameter movingWindow_size
.testParam_notInteger.m("movingWindow_size", movingWindow_size)
## CHECK parameter movingWindow_step
.testParam_notInteger.m("movingWindow_step", movingWindow_step)
cat("\n\n #------------------------------------------------------------#")
cat("\n # POST_FATE.graphic_evolutionStability")
cat("\n #------------------------------------------------------------# \n")
#############################################################################
res = foreach (abs.simulParam = abs.simulParams) %do%
{
cat("\n+++++++\n")
cat("\n Simulation name : ", name.simulation)
cat("\n Simulation file : ", abs.simulParam)
cat("\n")
## Get results directories ------------------------------------------------
GLOB_DIR = .getGraphics_results(name.simulation = name.simulation
, abs.simulParam = abs.simulParam)
## Get number of PFGs -----------------------------------------------------
## Get PFG names ----------------------------------------------------------
GLOB_SIM = .getGraphics_PFG(name.simulation = name.simulation
, abs.simulParam = abs.simulParam)
## Get raster mask --------------------------------------------------------
GLOB_MASK = .getGraphics_mask(name.simulation = name.simulation
, abs.simulParam = abs.simulParam)
## Get the abundance table ------------------------------------------------
file.abundance = paste0(name.simulation
, "/RESULTS/POST_FATE_TABLE_PIXEL_evolution_abundance_"
, basename(GLOB_DIR$dir.save)
, ".csv")
.testParam_existFile(file.abundance)
tab.totalAbundance = fread(file.abundance)
tab.totalAbundance = as.data.frame(tab.totalAbundance, stringsAsFactors = FALSE)
years = colnames(tab.totalAbundance)
years = years[which(!(years %in% c("PFG", "ID.pixel", "X", "Y", "HAB")))]
years = as.numeric(years)
no_years = length(years)
hab_names = sort(unique(tab.totalAbundance$HAB))
no_hab = length(hab_names)
cat("\n Number of years : ", no_years)
cat("\n Number of habitat : ", no_hab)
cat("\n")
###########################################################################
## Transform the data inside the table ------------------------------------
cat("\n ---------- GETTING TOTAL ABUNDANCE and EVENNESS within each habitat... \n")
## Getting total abundance within each habitat
cat("\n> Habitat total abundance...")
tab.split = split(tab.totalAbundance, list(tab.totalAbundance$HAB))
tab.hab = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
hab = names(tab.split)[i]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
tab = tab[, as.character(years), drop = FALSE]
return(data.frame(HAB = hab, year = years
, totalAbundance = colSums(tab, na.rm = TRUE)
, stringsAsFactors = FALSE))
}
}
## Getting PFG relative abundance within each habitat
cat("\n> PFG relative abundance...")
tab.split = split(tab.totalAbundance, list(tab.totalAbundance$PFG
, tab.totalAbundance$HAB))
tab.prop = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
pfg = strsplit(names(tab.split)[i], "[.]")[[1]][1]
hab = strsplit(names(tab.split)[i], "[.]")[[1]][2]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
tab = tab[, as.character(years), drop = FALSE]
tot.pfg = colSums(tab, na.rm = TRUE)
tot.hab = tab.hab[which(tab.hab$HAB == hab), ]
rownames(tot.hab) = tot.hab$year
prop = tot.pfg / tot.hab[as.character(years), "totalAbundance"]
return(data.frame(PFG = pfg, HAB = hab, year = years
, PROP = prop, stringsAsFactors = FALSE))
}
}
## Calculating evenness within each habitat
cat("\n> Habitat evenness...")
tab.split = split(tab.prop, list(tab.prop$HAB, tab.prop$year))
tab.even = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
hab = strsplit(names(tab.split)[i], "[.]")[[1]][1]
ye = strsplit(names(tab.split)[i], "[.]")[[1]][2]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
even = (- sum(tab$PROP * log(tab$PROP)) / log(nrow(tab)))
return(data.frame(HAB = hab, year = ye, no.PFG = nrow(tab)
, evenness = even, stringsAsFactors = FALSE))
}
}
## Merging total abundance and evenness for each habitat
if (nrow(tab.even) > 0) {
tab.HAB = merge(tab.hab, tab.even, by = c("HAB", "year"))
} else {
tab.HAB = tab.hab
}
cat("\n")
write.csv(tab.HAB
, file = paste0(name.simulation
, "/RESULTS/POST_FATE_TABLE_HAB_evolution_stability1_"
, basename(GLOB_DIR$dir.save)
, ".csv")
, row.names = FALSE)
message(paste0("\n The output file \n"
, " > POST_FATE_TABLE_HAB_evolution_stability1_"
, basename(GLOB_DIR$dir.save)
, ".csv \n"
, "has been successfully created !\n"))
###########################################################################
cat("\n ---------- EVALUATE, if possible, ABUNDANCE and EVENNESS STABILITY within each habitat... \n")
.getSucc = function(vec.years)
{
if (length(vec.years) > 1)
{
## gaps between years
gaps = vec.years[2:length(vec.years)] - vec.years[1:(length(vec.years) - 1)]
## periods with successive years
gaps = ifelse(gaps == movingWindow_step, 1, 0)
succ = c(gaps[1], rep(0, length(gaps) - 1))
for (i in 2:length(gaps))
{
if (gaps[i] == 1){
succ[i] = succ[i - 1] + 1
} else {
succ[i] = 0
}
}
sel = vec.years[which(succ > 0) + 1]
if (succ[1] == 1) sel = c(vec.years[1], sel)
return(list(selected = sel, max = max(succ, na.rm = TRUE) + 1))
} else
{
return(list(selected = vec.years, max = length(vec.years)))
}
}
## Calculate number of successive years -----------------------------------
max_years.succ = .getSucc(years)$max
years.succ = .getSucc(years)$selected
no_years.succ = length(years.succ)
no_movingWindow = no_years.succ - movingWindow_size + 1
## Evaluate if stability can be done --------------------------------------
if (length(years.succ) < movingWindow_size ||
movingWindow_size > max_years.succ ||
no_movingWindow < 1)
{
tab.STAB = pp = NULL
warning(paste0("Missing data!\n No moving window can be defined "
, "with the current parameters and available data :\n"
, " > mw_size = ", movingWindow_size, "\n"
, " > number of successive years = ", no_years.succ, "\n"
, " > maximum number of successive years = ", max_years.succ
, "\n Please check."))
} else
{
tab.split = split(tab.HAB, list(tab.HAB$HAB))
tab.STAB = foreach(i = 1:length(tab.split), .combine = "rbind") %do%
{
hab = names(tab.split)[i]
tab = tab.split[[i]]
if (nrow(tab) > 0)
{
combi = expand.grid(mw = 1:no_movingWindow, metric = c("totalAbundance", "evenness")
, stringsAsFactors = FALSE)
tab.mw = foreach(mw = combi$mw
, metric = combi$metric
, .combine = "rbind"
) %do%
{
mw_years = years[mw:(mw + movingWindow_size - 1)]
mw_years = .getSucc(mw_years)$selected
tab.mw = tab[which(tab$year %in% c(mw_years)), ]
if (nrow(tab.mw) > 0)
{
val = sum(tab.mw[, metric])
val.sq = sum(tab.mw[, metric] * tab.mw[, metric])
val.mean = val / movingWindow_size
val.sd = sqrt(val.sq / movingWindow_size - val.mean * val.mean)
val.cv = val.sd / val.mean
return(data.frame(HAB = hab
, no.years = length(mw_years)
, yearStep = movingWindow_step
, yearStart = min(mw_years)
, yearEnd = max(mw_years)
, metric = metric
, mean = val.mean
, sd = val.sd
, cv = val.cv
, stringsAsFactors = FALSE))
}
}
if (nrow(tab.mw) > 0)
{
return(tab.mw)
}
}
}
write.csv(tab.STAB
, file = paste0(name.simulation
, "/RESULTS/POST_FATE_TABLE_HAB_evolution_stability2_"
, basename(GLOB_DIR$dir.save)
, ".csv")
, row.names = FALSE)
message(paste0("\n The output file \n"
, " > POST_FATE_TABLE_HAB_evolution_stability2_"
, basename(GLOB_DIR$dir.save)
, ".csv \n"
, "has been successfully created !\n"))
## Compare quantiles ----------------------------------------------------
## TO BE DONE ??
## produce the plot -----------------------------------------------------
if (opt.doPlot && !is.null(tab.HAB))
{
cat("\n ---------- PRODUCING PLOT \n")
col_vec = c("#1B9E77", "#D95F02", "#7570B3", "#E7298A", "#66A61E", "#E6AB02", "#A6761D", "#666666")
col_fun = colorRampPalette(col_vec)
## Evolution of abundance and evenness through time -------------------
tab.plot1 = melt(tab.HAB, id.vars = c("HAB", "year"))
colnames(tab.plot1) = c("HAB", "year", "metric", "value")
tab.plot1$metric = factor(tab.plot1$metric, c("totalAbundance", "evenness", "no.PFG"))
## plot
pp = ggplot(tab.plot1, aes_string(x = "year", y = "value", color = "factor(HAB, hab_names)"))
if (!is.null(tab.STAB))
{
## Evolution of stability through time --------------------------------
tab.plot2 = tab.STAB
tab.plot2$year_median = sapply(1:nrow(tab.plot2), function(x) {
median(as.numeric(as.character(tab.plot2[x, c("yearStart", "yearEnd")]))) })
tab.plot2 = tab.plot2[which(tab.plot2$sd > 0.00001), ]
## plot
pp = pp +
geom_rect(data = tab.plot2, inherit.aes = FALSE, alpha = 0.5
, aes_string(xmin = "yearStart", xmax = "yearEnd"
, ymin = "mean - sd"
, ymax = "mean + sd"
, fill = "factor(HAB, hab_names)"))
}
## plot
pp = pp +
geom_line(lwd = 1) +
geom_point() +
facet_grid("metric ~ .", scales = "free_y"
, labeller = as_labeller(c("totalAbundance" = "Total abundance"
, "evenness" = "Evenness"
, "no.PFG" = "Number of PFG"))) +
scale_color_manual("Habitat", values = col_fun(no_hab)) +
scale_fill_manual(guide = "none", values = col_fun(no_hab)) +
labs(x = "", y = ""
, title = paste0("GRAPH A : evolution of habitat composition")) +
.getGraphics_theme()
ggsave(filename = paste0(name.simulation
, "/RESULTS/POST_FATE_GRAPHIC_A_evolution_stability_"
, basename(GLOB_DIR$dir.save), ".pdf")
, plot = pp, width = 10, height = 8)
} else
{
pp = NULL
} ## END opt.doPlot
}
## ------------------------------------------------------------------------
cat("\n> Done!\n")
return(list(tab.hab = tab.HAB
, tab.stab = tab.STAB
, plot.stab = pp))
} ## END loop on abs.simulParams
names(res) = abs.simulParams
return(res)
}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.