BIOMOD_PresenceOnly {biomod2} | R Documentation |
Evaluate models with presence-only metrics
Description
This function computes presence-only evaluation metrics (Boyce index and Minimal
Predicted Area) for BIOMOD.models.out
or BIOMOD.ensemble.models.out
objects that can be obtained with the BIOMOD_Modeling
or
BIOMOD_EnsembleModeling
functions.
Usage
BIOMOD_PresenceOnly(bm.mod = NULL, bm.em = NULL, bg.env = NULL, perc = 0.9)
Arguments
bm.mod |
a |
bm.em |
a |
bg.env |
(optional, default |
perc |
a |
Details
em.by
parameter of BIOMOD_EnsembleModeling
must have been set to
PA+run
in order to have an ensemble for each RUN of the
NbRunEval
parameter of the BIOMOD_Modeling
function for evaluation.
The Boyce index returns NA
values for SRE
models because it can not be
calculated with binary predictions.
This is also the reason why some NA
values
might appear for GLM
models if they do not converge.
Value
A data.frame
containing evaluation scores both for the evaluation metrics used in the
BIOMOD_Modeling
function and additional Boyce index and Minimal Predicted Area.
Note
In order to break dependency loop between packages biomod2 and ecospat,
code of ecospat.boyce()
and ecospat.mpa()
in ecospat)
functions have been copied within this file from version 3.2.2 (august 2022).
Author(s)
Frank Breiner, Maya Gueguen
References
Engler, R., Guisan, A., and Rechsteiner L. 2004. An improved approach for predicting the distribution of rare and endangered species from occurrence and pseudo-absence data. Journal of Applied Ecology, 41(2), 263-274.
Hirzel, A. H., Le Lay, G., Helfer, V., Randin, C., and Guisan, A. 2006. Evaluating the ability of habitat suitability models to predict species presences. Ecological Modelling, 199(2), 142-152.
See Also
ecospat.boyce()
and ecospat.mpa()
in ecospat,
BIOMOD.models.out
, BIOMOD_Modeling
,
BIOMOD.ensemble.models.out
, BIOMOD_EnsembleModeling
Other Main functions:
BIOMOD_EnsembleForecasting()
,
BIOMOD_EnsembleModeling()
,
BIOMOD_FormatingData()
,
BIOMOD_LoadModels()
,
BIOMOD_ModelingOptions()
,
BIOMOD_Modeling()
,
BIOMOD_Projection()
,
BIOMOD_RangeSize()
,
BIOMOD_Tuning()
Examples
library(terra)
# Load species occurrences (6 species available)
data(DataSpecies)
head(DataSpecies)
# Select the name of the studied species
myRespName <- 'GuloGulo'
# Get corresponding presence/absence data
myResp <- as.numeric(DataSpecies[, myRespName])
# Get corresponding XY coordinates
myRespXY <- DataSpecies[, c('X_WGS84', 'Y_WGS84')]
# Load environmental variables extracted from BIOCLIM (bio_3, bio_4, bio_7, bio_11 & bio_12)
data(bioclim_current)
myExpl <- terra::rast(bioclim_current)
# --------------------------------------------------------------- #
file.out <- paste0(myRespName, "/", myRespName, ".AllModels.models.out")
if (file.exists(file.out)) {
myBiomodModelOut <- get(load(file.out))
} else {
# Format Data with true absences
myBiomodData <- BIOMOD_FormatingData(resp.var = myResp,
expl.var = myExpl,
resp.xy = myRespXY,
resp.name = myRespName)
# Create default modeling options
myBiomodOptions <- BIOMOD_ModelingOptions()
# Model single models
myBiomodModelOut <- BIOMOD_Modeling(bm.format = myBiomodData,
modeling.id = 'AllModels',
models = c('RF', 'GLM'),
bm.options = myBiomodOptions,
CV.strategy = 'random',
CV.nb.rep = 2,
CV.perc = 0.8,
metric.eval = c('TSS','ROC'),
var.import = 3,
seed.val = 42)
}
file.EM <- paste0(myRespName, "/", myRespName, ".AllModels.ensemble.models.out")
if (file.exists(file.EM)) {
myBiomodEM <- get(load(file.EM))
} else {
# Model ensemble models
myBiomodEM <- BIOMOD_EnsembleModeling(bm.mod = myBiomodModelOut,
models.chosen = 'all',
em.by = 'all',
em.algo = c('EMmean', 'EMca'),
metric.select = c('TSS'),
metric.select.thresh = c(0.7),
metric.eval = c('TSS', 'ROC'),
var.import = 3,
seed.val = 42)
}
# --------------------------------------------------------------- #
# # Evaluate models with Boyce index and MPA
# myBiomodPO <- BIOMOD_PresenceOnly(bm.mod = myBiomodModelOut,
# bm.em = myBiomodEM)
# myBiomodPO
#
# # Evaluate models with Boyce index and MPA (using background data)
# myBiomodPO <- BIOMOD_PresenceOnly(bm.mod = myBiomodModelOut,
# bm.em = myBiomodEM,
# bg.env = myExpl)
# myBiomodPO