jSDM is an R package for fitting joint species distribution models (JSDM) in a hierarchical Bayesian framework.

The Gibbs sampler is written in C++. It uses Rcpp, Armadillo and GSL to maximize computation efficiency.



The package includes the following functions to fit various species distribution models :

  • jSDM_binomial_probit :

    Ecological process: $$y_{ij} \sim \mathcal{B}ernoulli(\theta_{ij})$$ where

    if n_latent=0 and site_effect="none"probit\((\theta_{ij}) = X_i \beta_j\)
    if n_latent>0 and site_effect="none"probit\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j\)
    if n_latent=0 and site_effect="fixed"probit\((\theta_{ij}) = X_i \beta_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
    if n_latent>0 and site_effect="fixed"probit\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j + \alpha_i\)
    if n_latent=0 and site_effect="random"probit\((\theta_{ij}) = X_i \beta_j + \alpha_i\)
    if n_latent>0 and site_effect="random"probit\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
  • jSDM_binomial_probit_sp_constrained :

    This function allows to fit the same models than the function jSDM_binomial_probit except for models not including latent variables, indeed n_latent must be greater than zero in this function. At first, the function fit a JSDM with the constrained species arbitrarily chosen as the first ones in the presence-absence data-set. Then, the function evaluates the convergence of MCMC \(\lambda\) chains using the Gelman-Rubin convergence diagnostic (\(\hat{R}\)). It identifies the species (\(\hat{j}_l\)) having the higher \(\hat{R}\) for \(\lambda_{\hat{j}_l}\). These species drive the structure of the latent axis \(l\). The \(\lambda\) corresponding to this species are constrained to be positive and placed on the diagonal of the \(\Lambda\) matrix for fitting a second model. This usually improves the convergence of the latent variables and factor loadings. The function returns the parameter posterior distributions for this second model.

  • jSDM_binomial_logit :

    Ecological process : $$y_{ij} \sim \mathcal{B}inomial(\theta_{ij},t_i)$$ where

    if n_latent=0 and site_effect="none"logit\((\theta_{ij}) = X_i \beta_j\)
    if n_latent>0 and site_effect="none"logit\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j\)
    if n_latent=0 and site_effect="fixed"logit\((\theta_{ij}) = X_i \beta_j + \alpha_i\)
    if n_latent>0 and site_effect="fixed"logit\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j + \alpha_i\)
    if n_latent=0 and site_effect="random"logit\((\theta_{ij}) = X_i \beta_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
    if n_latent>0 and site_effect="random"logit\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
  • jSDM_poisson_log :

    Ecological process : $$y_{ij} \sim \mathcal{P}oisson(\theta_{ij})$$ where

    if n_latent=0 and site_effect="none"log\((\theta_{ij}) = X_i \beta_j\)
    if n_latent>0 and site_effect="none"log\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j\)
    if n_latent=0 and site_effect="fixed"log\((\theta_{ij}) = X_i \beta_j + \alpha_i\)
    if n_latent>0 and site_effect="fixed"log\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j + \alpha_i\)
    if n_latent=0 and site_effect="random"log\((\theta_{ij}) = X_i \beta_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
    if n_latent>0 and site_effect="random"log\((\theta_{ij}) = X_i \beta_j + W_i \lambda_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
  • jSDM_binomial_probit_long_format :

    Ecological process: $$y_n \sim \mathcal{B}ernoulli(\theta_n)$$ such as \(species_n=j\) and \(site_n=i\), where

    if n_latent=0 and site_effect="none"probit\((\theta_n) = D_n \gamma + X_n \beta_j\)
    if n_latent>0 and site_effect="none"probit\((\theta_n) = D_n \gamma + X_n \beta_j + W_i \lambda_j\)
    if n_latent=0 and site_effect="fixed"probit\((\theta_n) = D_n \gamma + X_n \beta_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)
    if n_latent>0 and site_effect="fixed"probit\((\theta_n) = D_n \gamma + X_n \beta_j + W_i \lambda_j + \alpha_i\)
    if n_latent=0 and site_effect="random"probit\((\theta_n) = D_n \gamma + X_n \beta_j + \alpha_i\)
    if n_latent>0 and site_effect="random"probit\((\theta_n) = D_n \gamma + X_n \beta_j + W_i \lambda_j + \alpha_i\) and \(\alpha_i \sim \mathcal{N}(0,V_\alpha)\)


Chib, S. and Greenberg, E. (1998) Analysis of multivariate probit models. Biometrika, 85, 347-361.

Warton, D. I.; Blanchet, F. G.; O'Hara, R. B.; O'Hara, R. B.; Ovaskainen, O.; Taskinen, S.; Walker, S. C. and Hui, F. K. C. (2015) So Many Variables: Joint Modeling in Community Ecology. Trends in Ecology & Evolution, 30, 766-779.

Ovaskainen, O., Tikhonov, G., Norberg, A., Blanchet, F. G., Duan, L., Dunson, D., Roslin, T. and Abrego, N. (2017) How to make more out of community data? A conceptual framework and its implementation as models and software. Ecology Letters, 20, 561-576.


Ghislain Vieilledent <ghislain.vieilledent@cirad.fr>

Jeanne Clément <jeanne.clement16@laposte.net>

Frédéric Gosselin <frederic.gosselin@inrae.fr>