This function fits the isocape as a mixed model. The fitting procedures are done by the package spaMM which we use to jointly fit the mean isotopic values and their associated residual dispersion variance in a spatially explicit manner.
isofit(iso.data, mean.model.fix = list(elev = TRUE, lat.abs = TRUE, lat.2 =
FALSE, long = FALSE, long.2 = FALSE), disp.model.fix = list(elev = FALSE,
lat.abs = FALSE, lat.2 = FALSE, long = FALSE, long.2 = FALSE),
mean.model.rand = list(uncorr = TRUE, spatial = TRUE),
disp.model.rand = list(uncorr = TRUE, spatial = TRUE),
uncorr.terms = list(mean.model = "lambda", disp.model = "lambda"),
spaMM.method = list(mean.model = "fitme", disp.model = "fitme"),
dist.method = "Earth", control.mean = list(), control.disp = list(),
verbose = interactive())
The dataframe containing the data used for fitting the isoscape model
A list of logical indicating which fixed effects to consider in mean.fit
A list of logical indicating which fixed effects to consider in disp.fit
A list of logical indicating which random effects to consider in mean.fit
A list of logical indicating which random effects to consider in disp.fit
A list of two strings defining the parametrization used to model the uncorrelated random effects for mean.fit and disp.fit
A list of two strings defining the spaMM functions used for mean.fit and disp.fit
A string indicating the distance method
A list of additional arguments to be passed to the call of mean.fit
A list of additional arguments to be passed to the call of disp.fit
A logical indicating whether information about the
progress of the procedure should be displayed or not while the function is
running. By default verbose is TRUE
if users use an interactive R
session and FALSE
otherwise.
This function returns a list of class isofit
containing
two inter-related fits: mean.fit
and disp.fit
. The returned
list also contains the object info.fit
that contains all the
call arguments.
The detailed statistical definition of the isoscape model will be soon
available in another document. Briefly, the fitting procedure of the
isoscape model is divided into two fits: mean.fit
and
disp.fit
. mean.fit
corresponds to the fit of the "mean model",
which we will use to predict the mean isotopic values at any location in
other functions of the package. disp.fit
corresponds to the fit of
the "residual dispersion model", which we will use to predict the residual
dispersion variance associated to the mean predictions. mean.fit
is a
linear mixed-effects model (LMM) with fixed effects, an optional spatial
random effect with a Matern correlation structure and an optional
uncorrelated random effect accounting for variation between weather station
unrelated to their location. disp.fit
is a Gamma Generalized LMM
(Gamma GLMM) that also has fixed effects, an optional spatial random effect
with a Matern correlation structure and an optional uncorrelated random
effect. For the GLMM the residual variance is fixed to its theoretical
expectation.
The dataframe iso.data
must contain a single row per source
location with the following columns: isoscape.value
(the isotopic
value), var.isoscape.value
(the unbiased variance estimate of the
isotopic value at the location), n.isoscape.value
(the number of
measurements performed at the location, could be 1) and stationID
(a
factor defining the identity of the sources at a given location).
The arguments mean.model.fix
and disp.model.fix
allow the user
to choose among different fixed-effect structures for each model. These
arguments are lists of booleans (TRUE
or FALSE
), which define
which of the following fixed effects must be considered: the elevation
(elev
), the absolute value of the latitude (lat.abs
), the
squared latitude (lat.2
), the longitude (long
) and the squared
longitude (long.2
). An intercept is always considered in both models.
By default, only intercept are being fitted.
In the models, the mean (for the mean model) or the log residual variance
(for the residual dispersion model) follow a Gaussian distribution around a
constant value. The arguments mean.model.rand
and
disp.model.rand
allow to choose among different random effects for
each model influencing the realizations of these Gaussian random processes.
For each model one can choose not to include random effects or to include an
uncorrelated random effect, a spatial random effect, or both (default).
Setting "uncorr" = TRUE
implies that the different realizations are
identical for a given weather station (e.g. some micro-climate or some
measurement errors trigger a shift in all measurements (mean model) or a
shift in the variance between measurements (residual dispersion model)
performed on a given weather station by the same amount). Setting
"spatial" = TRUE
(default) implies that the random realizations of the
Gaussian process follow a Matern correlation structure. Put simply, this
implies that the closer two locations are, the more similar the means (or
the log residual variance) in isotopic values are (e.g. because they are
likely to be traversed by the same air masses).
The arguments uncorr.terms
allow the choice between two alternative
parameterizations for the uncorrelated random effect in the fits:
"lambda"
or "nugget"
for each model. When using
"lambda"
, the variance of the uncorrelated random terms is
classically modeled by a variance. When a spatial random effect is
considered, one can alternatively choose "nugget"
, which modifies the
Matern correlation value when distance between location tends to zero. If no
random effect is considered, one should stick to the default setting and it
will be ignored by the function. The choice of the parametrization is a
matter of personal preferences and it does not change the underlying models,
so the estimations for all the other parameters of the models will not be
impacted by whether one chooses lambda
or nugget
. Depending on
the data one parametrization may lead to faster fit than the other.
The argument spaMM.method
is also a list of two strings where
the first element defines the spaMM functions used for fitting the mean
model and the second element defines the spaMM method used for fitting the
residual dispersion model. The possible options are "HLfit", "corrHLfit" and
"fitme". Note that "HLfit" shall only be used in the absence of a Matern
correlation structure and "corrHLfit" shall only be used in the presence of
it. In contrast, "fitme" should work in all situations. Which method is best
remains to be determined and it is good practice to try different methods
(if applicable) to check for the robustness of the results. If all is well
one should obtain very similar results with the different methods. If this
is not the case, carefully check the model output to see if one model fit
did not get stuck at a local minimum during optimization (which would
translate in a lower likelihood).
The argument dist.method
allows modifying how the distance between
locations is computed to estimate the spatial correlation structure. By
default, we consider the so-called "Earth" distances which are technically
called orthodromic distances. They account for earth curvature. The
alternative "Euclidean" distances do not. For studies performed on a small
geographic scale, both distance methods should lead to similar results.
The arguments control.mean
and control.dist
are lists that are
transmitted to the spaMM fitting functions (defined by
spaMM.method
). These lists can be used to finely control the fitting
procedure, so advanced knowledge of the package spaMM is
required before messing around with these inputs.
We highly recommend users to examine the output produced by isofit
.
Sometimes, poor fit may occur and such models should therefore not be used
for building isoscapes or performing assignments.
Rousset, F., Ferdy, J. B. (2014). Testing environmental and genetic effects in the presence of spatial autocorrelation. Ecography, 37(8):781-790.
Bowen, G. J., Wassenaar, L. I., Hobson, K. A. (2005). Global application of stable hydrogen and oxygen isotopes to wildlife forensics. Oecologia, 143(3):337-348.
spaMM for an overview of the spaMM package
fitme
and corrHLfit
for
information about the two possible fitting procedures that can be used here
Matern.corr
for information about the Matern
correlation structure
IsoriX
for the complete work-flow of our package
# NOT RUN {
## The examples below will only be run if sufficient time is allowed
## You can change that by typing e.g. IsoriX.options(example_maxtime = XX)
## if you want to allow for examples taking up to ca. XX seconds to run
## (so don't write XX but put a number instead!)
if(IsoriX.getOption("example_maxtime") > 10) {
## Fitting the models for Germany
GNIPDataDEagg <- prepdata(data = GNIPDataDE)
GermanFit <- isofit(iso.data = GNIPDataDEagg)
GermanFit
## Diagnostics for the fits
plot(GermanFit)
}
# }
Run the code above in your browser using DataLab