1
$\begingroup$

I have data on 3 species abundances as fractional compositions. Within each sample, the sum of all species abundances is 1. Essentially, each observation is a proportional abundance. Here is an example of this type of data using R:

y1 <- rnorm(100, 100, 5)
y2 <- y1*2 + rnorm(100,0,10)
y3 <- rnorm(100,50,5)
spp.y <- data.frame(y1,y2,y3)
tot <- rowSums(spp.y)
spp.y <- spp.y / tot
spp.y <- as.matrix(spp.y)

In this example, the abundances of y1 and y2 are correlated. y1 and y2 are uncorrelated with y3. I would like to specify a model that can capture and quantify the magnitude and sign of this correlation, while accounting for the fact that these observations are constrained by summing to 1. To do this I specified a dirlichet model in JAGS, the multivariate extension of the beta model:

Setup JAGS data as an R object, jags.data:

jags.data <- list(y = spp.y,
                  N = nrow(spp.y),
                  N.spp = ncol(spp.y))

Specify dirlichet model as an R object, jags.model:

jags.model = "
model {
    #setup priors for each species
    for(j in 1:N.spp){
      m0[j] ~ dgamma(1.0E-3, 1.0E-3) #intercept prior
    }

    #implement dirlichet
    for(i in 1:N){
    for(j in 1:N.spp){
         log(a0[i,j]) <- m0[j] ## eventually add linear model here
      }
     y[i,1:N.spp] ~ ddirch(a0[i,1:N.spp]) 
    }

} #close model loop.
"

Fit the model in R using the runjags package:

library(runjags)
jags.out <- run.jags(jags.model,
                     data=jags.data,
                     adapt = 100,
                     burnin = 100,
                     sample = 200,
                     n.chains=3,
                     monitor=c('m0'))

This works well. The intercepts capture the within species relative abundances.

colMeans(spp.y)
a.m <- exp(summary(jags.out)[,"Median"])
a.m/sum(a.m)

However, this model assumes the species abundances are independent, which they are not. Furthermore, I want to capture and quantify the correlation between y1 and y2. I believe this is possible by adding covariance among species abundances, but I do not know how to implement this in JAGS. Any advice (and especially code!) would be super appreciated.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.