2
$\begingroup$

I'm attempting to use latent variable modelling, as described by Hui (2016) and using the boral package in R, to explore relationships between plant species composition and environmental variables while allowing for possible interactions between species. The species data consist of presence-absence of 165 taxa over 66 survey sites.

I'd like to estimate the overall proportion of co-variation between species accounted for by the environmental variables. Hui (2016) and Warton et al. (2015) describe using the trace of the residual covariance matrix for this purpose. They compare the trace value of a model with latent variables but no covariates to that from a model with latent variables plus covariates. One expects that the trace value will be lower for the second model, with the proportional reduction being an estimate of how much of the species co-occurrence is explained by covariates. The boral package provides the function get.residual.cor which returns, among other things, the trace value.

Attempting to apply the suggested method to my data and models gives surprising results. Alternative models with environmental covariates included are returning trace values substantially higher (e.g. 30%) than that for a base model with only latent variables. This leads me to worry that the models are somehow badly mis-posed but, if that is the case, it isn't obvious. Fitted coefficients for species against environmental variables appear to make ecological sense and posterior predictions (e.g. for species richness of field survey sites) seem reasonable.

A boral model is, behind the scenes, a Bayesian model fitted using MCMC via JAGS. I understand (in a rudimentary way) that this means the residual covariance matrix should be treated with the same caveats as any other parameter sampled by MCMC. The package documentation for function get.residual.cor states:

Of course, the trace itself is random due to the MCMC sampling, and so it is not always guaranteed to produce sensible answers!

So my question is how does one know when to trust the trace value? Are there any recommended post-hoc checks? My models appear to have converged (as judged by the Geweke's diagnostic returned by the boral function) and auto-correlation of parameter samples seems acceptable. Perhaps there is a posterior sample size requirement for the trace to be reliable? I'm presently running 100,000 MCMC iterations and thinning to every 25th sample.

Improve this question
$\endgroup$
2
  • $\begingroup$ The word "trace" has multiple meanings. You are interested in the trace of a matrix (the sum of the main diagonal elements). By contrast, a "trace plot" is a standard MCMC graphical diagnostic tool. I mention this because your question could involve both meanings. In fact, you could make a trace plot of the trace! $\endgroup$
    – mef
    Commented Mar 22, 2018 at 9:01
  • $\begingroup$ @mef Thanks. I meant the trace of a matrix (the residual covariance matrix estimated from the MCMC samples). $\endgroup$
    – michael
    Commented Mar 24, 2018 at 7:24

3 Answers 3

Reset to default
1
$\begingroup$

Have you checked out the get.enviro.cor function implemented in Boral? Sounds exactly what you're looking for. This retrieves the covariance(or correlation) between species due to shared environmental response.

Improve this answer
$\endgroup$
1
  • $\begingroup$ Yes, the get.enviro.cor function is really useful and a nice complement to the get.residual.cor function to look for patterns at the species level. However, I was hoping to also report on an overall measure of species covariation explained by alternative models, hence my interest in the covariance matrix trace values. $\endgroup$
    – michael
    Commented Mar 28, 2018 at 0:59
1
$\begingroup$

I have had same problem with the residual correlation trace in BORAL. I have not found a solution yet.
An alternative approach you can use is the var.part function, this partitions the variance explained by the latent variables and the covariates in your model.

Improve this answer
$\endgroup$
2
  • $\begingroup$ If I understand correctly, the calc.varpart function decomposes variance at the species level. Ideally I would like to report on something similar at the overall level. $\endgroup$
    – michael
    Commented Jul 30, 2018 at 6:23
  • $\begingroup$ I see your point. I know this is not exactly what you want, but you could take the mean/median of the variance decomposition across your set of species. This could act as a proxy to show the average variance explain by your covariates and the residual variance. One more thing,how does the ordination of your correlated response model compare to the pure latent variable model? are the points less dispersed? $\endgroup$
    – Ch16S
    Commented Aug 8, 2018 at 0:47
0
$\begingroup$

Ovaskainen et al. (2017), used Tjur’s R2 (Tjur 2009) averaged over species, as an overall community effect.

Ovaskainen et al. (2017) Ecology Letters, 20: 561–576 doi: 10.1111/ele.12757

Improve this answer
$\endgroup$
2
  • $\begingroup$ Say a bit more as to why you think this answers the question, please. $\endgroup$
    – Carl
    Commented Sep 5, 2018 at 23:33
  • $\begingroup$ Many thanks for the reply. I'd also be interested in more detail please - if that's possible. More generally, I wonder if my main problem is seeking a single-value summary statistic for the contribution of covariates in the first place. Following @Ch16S 's suggestion to use the calc.varpart function, perhaps comparing density curves for the proportion of variance attributed to the latent variables for each response (species in my case) under alternative models would be more instructive and a little less prone to artefacts? $\endgroup$
    – michael
    Commented Sep 7, 2018 at 1:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.