model_gamm
accounts for temporal autocorrelation (TAC) in the time series
by fitting Generalized Additive Mixed Models (GAMMs) that include AR or ARMA
correlation structures using the gamm
function. The GAMMs
are applied to all IND~pressure combinations provided as input or only those
with significant TAC in the GAM residuals (using filter argument).
model_gamm( init_tbl, k = 5, family = stats::gaussian(), excl_outlier = NULL, filter = NULL, lme_control = NULL )
init_tbl  The output tibble of the 

k  Choice of knots (for the smoothing function 
family  A description of the error distribution and link to be used in the GAM.
This needs to be defined as a family function (see also 
excl_outlier  A list of values identified as outliers in specific
IND~pressure GAMMs, which should be excluded in this modeling step
(the output tibble of this function includes the variable
`pres_outlier`, which is a columnlist containing
all indices of values with cook`s distance > 1 (see below). The function
can be rerun again, then excluding all these outliers provided in

filter  logical; a filter used to select specific rows in init_tbl
(row gets selected if value TRUE). That could be the 
lme_control  A specification list of the control values for the lme fit.
The default is NULL so the control settings in 
Returns a model output tibble that contains for each filtered IND~pressure pair
6 rows with the individual GAMM outputs. The structure remains the same as in
model_gam
except for the explained deviance, which is not computed
by the gamm function.
The selection of the final correlation structure for each IND~pressure model can be
done manually on this tibble or with an automatized routine using
select_model
.
Modeling firstdifferenced indicator time series can be an alternative solution to avoid temporal dependence between observations. However, this approach does often not help reducing the significant autocorrelation while GAMMs do as found in Otto et al. (2018). Such an extension implies that the single elements of the response variable are not independent anymore and that the correlation between the residuals at time t1 and t2 only depends on their time difference t1 – t2 (Wood, 2006).
In model_gamm
six GAMMs are computed for each filtered IND~pressure pair, i.e.
no correlation structure (for AIC comparison)
autoregressive error structure of order p=1 (AR1)
autoregressive error structure of order p=2 (AR2)
autoregressive moving average of order p=1 and q=1 (ARMA11)
autoregressive moving average of order p=1 and q=2 (ARMA12)
autoregressive moving average of order p=2 and q=1 (ARMA21)
Otto, S.A., Kadin, M., Casini, M., Torres, M.A., Blenckner, T. (2018) A quantitative framework for selecting and validating food web indicators. Ecological Indicators, 84: 619631, doi: https://doi.org/10.1016/j.ecolind.2017.05.045
Wood, S.N. (2006) Generalized Additive Models: An Introduction with R. Chapman and Hall/CRC Press
gamm
for more information on GAMMs and
plot_diagnostics
for assessing the model diagnostics
Other IND~pressure modeling functions:
find_id()
,
ind_init()
,
model_gam()
,
plot_diagnostics()
,
plot_model()
,
scoring()
,
select_model()
,
test_interaction()
# Using the Baltic Sea demo data in this package dat_init < ind_init( ind_tbl = data.frame(Cod = ind_ex$Sprat), press_tbl = press_ex[, c("Fsprat", "Fher")], time = ind_ex[ ,1]) gam_tbl < model_gam(dat_init) # Any temporal autocorrelation gam_tbl$tac#> [1] TRUE TRUE# Applying model_gamm function and passing the $tac variable as filter gamm_tbl < model_gamm(dat_init, filter = gam_tbl$tac) # Use different control values for lme fit: gamm_tbl < model_gamm( dat_init, filter = gam_tbl$tac, lme_control = nlme::lmeControl(niterEM = 100, msMaxIter = 200) ) # \donttest{ # run analysis will full demo dataset gam_tbl < model_gam(ind_init_ex) gamm_tbl < model_gamm(ind_init_ex, filter = gam_tbl$tac) # now run the same analysis but exclude all outliers in the GAMMs gamm_tbl2 < model_gamm(ind_init_ex, filter = gam_tbl$tac, excl_outlier = gamm_tbl$pres_outlier) # }