Search code examples
rjagsr2jags

Assigning a prior distribution for multiple parameters with for loop


I've read the JAGS manual but did not find a way to assign the same prior distribution to multiple parameters in a JAGS / R2JAGS model.

For example, currently I have to repeat a lot of code like this:

reg.model <- function() {
  # Model structure
  for(i in 1:N){
    Y[i] ~ dnorm(mu[i], phi)
    mu[i] <- beta0 + beta1*x1[i] + beta2*x2[i]
  }
  sigma2 <- 1 / phi

  # Priors (Lots of re-typing here for beta0, beta1, beta2)
  phi   ~ dgamma(1,1)
  beta0 ~ dnorm(0, 0.01*phi)
  beta1 ~ dnorm(0, 0.01*phi)
  beta2 ~ dnorm(0, 0.01*phi)
}

How to DRY this code?


Solution

  • You can do this if you pass your independent variables as a matrix, put your coefficients into a vector, and use matrix multiplication (with inprod or %*%) to calculate your linear predictor.

    For example:

    M <- function() {
      for (i in 1:n) {
        y[i] ~ dnorm(mu[i], sd^-2)
        mu[i] <- X[i, ] %*% beta
      }
    
      sd ~ dunif(0, 100)
      for (j in 1:J) {
        beta[j] ~ dnorm(0, 0.0001)
      }
    }
    

    We generate some dummy data, below:

    set.seed(1)
    n <- 1000
    J <- 3
    X <- cbind(1, matrix(runif(n * J), ncol=J))
    
    head(X)
    #      [,1]      [,2]       [,3]      [,4]
    # [1,]    1 0.2655087 0.53080879 0.8718050
    # [2,]    1 0.3721239 0.68486090 0.9671970
    # [3,]    1 0.5728534 0.38328339 0.8669163
    # [4,]    1 0.9082078 0.95498800 0.4377153
    # [5,]    1 0.2016819 0.11835658 0.1919378
    # [6,]    1 0.8983897 0.03910006 0.0822944
    
    b <- rnorm(J + 1, 0, 10)
    sigma <- runif(1)
    y <- c(X %*% b) + rnorm(n, 0, sigma)
    

    And fit the model:

    library(R2jags)
    out <- jags(list(y=y, X=X, n=n, J=ncol(X)), NULL, c('beta', 'sd'), M)
    

    Now compare the estimated coefficients and standard deviation to their true values:

    out$BUGSoutput$summary[, '50%']
    #     beta[1]      beta[2]      beta[3]      beta[4]     deviance           sd 
    #   8.5783315   -9.3849277    8.9632598   -9.4242272 2196.1535645    0.7264754 
    
    b # True betas
    # [1]  8.500435 -9.253130  8.935812 -9.410097
    
    sigma # True sd
    # [1] 0.70504