Modeling community-level data controlling compositional autocorrelation.

The function uses fitme (package spaMM) to model compositional autocorrelation when analyzing community-level data. Phylogenetic relationship between species can be provided or optimized (assuming a Brownian motion in the latter case)

Usage

EcoCoMix(formula,
          data,
          VCV_sp,
          comm,
          force.PD=FALSE,
          optim.lambda=TRUE,
          original.VCV=TRUE,
          AIC_threshold=-4,
          method.spaMM="REML",
          true_VCV=NULL,
          control.optim=NULL,
          comm_kronecker=NULL,
          int_model=TRUE,
          ...)

Arguments

formula: a formula specified in spaMM format (see fitme). Use "comp_id"" to represent community-level random effect controlling compositional autocorrelation (see Example).
data: a data frame containing the variables in the model.
VCV_sp: a covariance (not correlation!) matrix illsutrating phylogenetic relationship between species
comm: a site(row)-species(column) matrix. Pleasure ensure species names are consistent between VCV_sp and comm.
force.PD: Forcing Positive Definite using nearPD. Defaulted as FALSE.
optim.lambda: TRUE/FALSE. Whether lambda should be optimized (through finding value with the lowest cAIC).
original.VCV: Return model results based on the provided VCV_sp. Defaulted as TRUE.
AIC_threshold: cAIC threshold to include compositional random effect. Defaulted as -4 (i.e., including compositional random effect needs to reduce cAIC by 4 to be considered as improving the model).
method.spaMM: "ML" (Default), "REML", or PGL/L. Will be passed to the method argument of fitme. See details for model comparisons when using REML.
true_VCV: True phylogenetic relationship between species. Useful in simulations but not empirical analyses.
control.optim: A list of control parameters for optim (See the documentation of optim when using this argument).
comm_kronecker: Arguments used to speed up internal calculations. Defaulted as NULL.
int_model: Obtain results for intercept-only models (both with and without the compositional random effects). Defaulted as TRUE.
...: Optional arguments passed to fitme.

Details

Users can provide extra specifications to the fitme function to speed up calculations/obtain more reliable results. These technical details, however, are non-trivial. Please see the documentation of fitme, inits, and optim. If ignored, the model will be run using the defaults of these functions.

Note that when models are fitted using REML, AICs are only comparable when the fixed effects are the same. Use ML if you want to compare models with different fixed effects!

Value

a list of 19 objects.

best_model: Results of the best model (i.e., the model with optimized lambda if its dAIC is lower than the intercept-only model by 4 (or other self-defined threshold.)
best_model_satt: Results of the best model based on the Satterthwaite approximation, which provides better estimations of p-values.
optimized_lambda_model,optimized_lambda_model_satt: Same as best_model and best_model_satt, but this is for the model with optimized phylogenetic structure.
original_VCV_model,original_VCV_m_satt: Same as above, but for models based on the provided phylogenetic structure.
true_model,true_model_satt: Same as above, but for models based on the true phylogenetic structure (only useful for simulations).
without_comp_model: Model results without considering compositional autocorrelation.
without_comp_anova: Model results without considering compositional autocorrelation in ANOVA format. Easier to compare with results from mixed models based on Satterthwaite approximation.
intercept_only_model: Intercept-only model results (with compositional correlation).
star_model: Model results without considering phylogenetic dependence between species (i.e., all species are phylogenetically independent).
AIC: cAIC of all models above.
optimized_lambda,optimized_lambda_int: The optimized lambda based on the full model (former) and intercept-only model(latter).
conv: whether the following models have reached convergence: the best model, the full model based on optimized phylogenetic structure, and the model based on the provided phylogenetic relationship. Yes=1, No=0.

Author

Toby P.N. Tsang

Examples

library(phytools)
library(tidyverse)
library(EcoCoMix)
library(ape)
library(DHARMa)

data(KSR)
data(KSR_MLtree)
data(KSR_EF)

m_LAI <- EcoCoMix(LAI~Real.rich+corrMatrix(1|comp_id),
          data=KSR_EF,
          comm=KSR,
          VCV_sp = vcv(KSR_MLtree),
          method.spaMM="REML")

m_LAI$best_model
m_LAI$best_model_satt
plot(simulateResiduals(m_LAI$optimized_lambda_model,re.form=NULL))