Florian Hartig

← Publications

A new joint species distribution model for faster and more accurate inference of species associations from big community data

Maximilian Pichler & Florian Hartig

Methods in Ecology and Evolution, 12(11), 2159–2173 (2021)
Cite this
@article{pichler2021new,
  author = {Maximilian Pichler and Florian Hartig},
  title = {A new joint species distribution model for faster and more accurate inference of species associations from big community data},
  journal = {Methods in Ecology and Evolution},
  volume = {12},
  number = {11},
  pages = {2159--2173},
  year = {2021},
  doi = {10.1111/2041-210X.13687},
}

DOI: 10.1111/2041-210X.13687
Citation statistics: 140 (Google Scholar) · 99 (OpenAlex), as of 29 October 2026

Google Scholar citations per year, as of 20 October 2026 (yearly counts are fetched separately from the total)

View article (DOI) Preprint Software: sjSDM

Abstract

Joint species distribution models (JSDMs) explain spatial variation in community composition by contributions of the environment, biotic associations and possibly spatially structured residual covariance. They show great promise as a general analytical framework for community ecology and macroecology, but current JSDMs, even when approximated by latent variables, scale poorly on large datasets, limiting their usefulness for currently emerging big (e.g. metabarcoding and metagenomics) community datasets. Here, we present a novel, more scalable JSDM (sjSDM) that circumvents the need to use latent variables by using a Monte Carlo integration of the joint JSDM likelihood together with flexible elastic net regularization on all model components. We implemented sjSDM in PyTorch, a modern machine learning framework, which allows making use of both CPU and GPU calculations. Using simulated communities with known species–species associations and different number of species and sites, we compare sjSDM with state‐of‐the‐art JSDM implementations to determine computational runtimes and accuracy of the inferred species–species and species–environment associations. We find that sjSDM is orders of magnitude faster than existing JSDM algorithms (even when run on the CPU) and can be scaled to very large datasets. Despite the dramatically improved speed, sjSDM produces more accurate estimates of species association structures than alternative JSDM implementations. We demonstrate the applicability of sjSDM to big community data using eDNA case study with thousands of fungi operational taxonomic units (OTU). Our sjSDM approach makes the analysis of JSDMs to large community datasets with hundreds or thousands of species possible, substantially extending the applicability of JSDMs in ecology. We provide our method in an R package to facilitate its applicability for practical data analysis.

What the paper shows and why it matters

(AI-generated)

Joint species distribution models (JSDMs) describe species communities by combining each species’ response to the environment with associations between species (pairs that occur together more or less often than the environment predicts). They suit large datasets from environmental DNA (eDNA) and metabarcoding but scale poorly: the full model’s covariance matrix grows quadratically with the number of species. The usual remedy replaces it by a few latent variables, hidden gradients summarising the associations (Warton et al. 2015). This is reasonable, as it should in theory reduce the variance of the estimates. Yet fitting can still take hours to days for hundreds of species, and this low-rank structure may bias the estimated associations. The paper questions the assumption that scalable JSDMs need latent variables.

The paper shows that joint species distribution models (JSDMs) with a full association matrix can be fitted fast and accurately on big community data without latent variables. Its method, sjSDM, approximates the likelihood by random sampling (Monte Carlo integration), which runs in parallel on a processor or a graphics card (GPU), and limits overfitting with elastic-net regularisation, a penalty shrinking weak associations and environmental effects towards zero. In simulations with known associations, sjSDM scaled far better than Hmsc, gllvm and BayesComm (gllvm was about as fast for fewer than 50 species), with a gap that grew with species number and reached orders of magnitude even without a GPU; on a GPU it fitted 30,000 sites with up to 1,000 species in about 4.5 minutes. It also recovered the signs of associations, and which are non-zero, better than the latent-variable packages Hmsc and gllvm, though accuracy stayed moderate; all predicted occurrences about equally well. The paper treats the latent-variable structure as the likely but unproven cause of the bias. As an example, it fits 3,649 fungal operational taxonomic units (OTUs) from 125 Danish sites.

Figure 5. Inferred operational taxonomic unit (OTU) associations and environmental preferences for the eDNA community data. From Pichler & Hartig (2021), Methods in Ecology and Evolution, licensed under CC BY 4.0. Click the image to enlarge.

The results make JSDMs usable for hundreds to thousands of species, the largest runs using a GPU, and suggest that shrinking a full covariance matrix can estimate associations better than a low-rank structure. The paper cautions that associations are not necessarily causal and that their uncertainty has no analytical estimate yet but can be obtained by resampling (bootstrapping), and offers the R package sjSDM. Rahman et al. (2024), when porting the JSDM package Hmsc to GPUs, state that, to their knowledge, sjSDM is the only existing JSDM package designed for GPUs.

References

  • Warton et al. (2015). So Many Variables: Joint Modeling in Community Ecology. Trends in Ecology & Evolution 30, 766–779. doi:10.1016/j.tree.2015.09.007
  • Rahman et al. (2024). Accelerating joint species distribution modelling with Hmsc-HPC by GPU porting. PLOS Computational Biology 20, e1011914. doi:10.1371/journal.pcbi.1011914