0

I'm using nimble package to run MCMC chains using parallel computation. Nimble normally has nice progress bars for each chain run in sequence with either the nimbleMCMC integrated function or the runMCMC component function.

With nimble you need to create a custom function containing all model code including the runMCMC function, and then using parLapply to run the custom function for separate chains in parallel on a cluster.

https://r-nimble.org/nimbleExamples/parallelizing_NIMBLE.html

This computes fine but it loses the original progress bar and does not show progress.

There are solutions here but I don't think they work with the Nimble framework: How to show the progress of code in parallel computation in R? It is possible to create a progress bar with pblapply, but this calculates progress as a % based on completed chains/total chains. If running 8 chains in parallel over 8 cores for example, this is essentially uninformative because it stays at 0% until all chains finish simultaneously in parallel and jumps to 100%.

Is there a way to print the original progress bar from the runMCMC from within the custom function? Or alternatively, show the progress in some other way?

This is a toy CJS hidden markov model to illustrate:

library(nimble)
library(parallel)
library(pbapply)

ymat <- matrix(sample(c(0,1),120,replace = TRUE),ncol = 6, nrow = 20)
y <- ymat

detectCores() # Number of cores on the machine = 8
this_cluster <- makeCluster(4) 

index <- c(1:4)

paraChain_code <- function(index, y){

library(nimble)  
first <- apply(y, 1, function(x) min(which(x!= 0),na.rm = TRUE))

hmm.model<- nimbleCode({
 
  # Priors
  phi ~ dunif(0,1)    # prior for survival φ
  p ~ dunif(0,1)      # prior for detection p
  
  # Vector of initial state probabilities δ
  delta[1] <- 1     # Alive
  delta[2] <- 0     # Dead
  
  # Transition matrix Γ
  gamma[1,1] <- phi      # transition Pr(alive t --> alive t+1)
  gamma[1,2] <- 1-phi    # transition Pr(alive t --> dead t+1)
  gamma[2,1] <- 0        # transition Pr(dead t --> alive t+1)
  gamma[2,2] <- 1        # transition Pr(dead t --> dead t+1)
  
  # Observation matrix Ω
  omega[1,1] <- 1-p   # observation Pr(alive t --> not detected t)
  omega[1,2] <- p     # observation Pr(alive t --> detected t)
  omega[2,1] <- 1     # observation Pr(dead t --> not detected t)
  omega[2,2] <- 0     # observation Pr(dead t --> detected t)
  
  # Likelihood
  for(i in 1:N){    # Loops over each individual (row)
    z[i,first[i]] ~ dcat(delta[1:2])  # Draws an initial state   
    for(j in (first[i]+1):T){    # Loops through each time point (column)
      z[i,j] ~ dcat(gamma[z[i,j-1], 1:2]) # Draws current state based on previous state
      y[i,j] ~ dcat(omega[z[i,j], 1:2])  # Draws an observation based on current state.
    }
  }
})

# Defines constants in the model
my.constants <- list(N = nrow(y),   # N individual
                     T = ncol(y),   # N of time points
                     first = first) # occasions of first capture

# Defines data in the model
# To remove 0s. 1 is non-detected, 2 is detected
my.data <- list(y = y+1)

# Defines initial values for states and probabilities
zinits <- y + 1
zinits[zinits == 2] <- 1
zinits[is.na(zinits) == TRUE] <- 1

initial.values <-  list(phi = runif(1,0,1),
                        p = runif(1,0,1),
                        z = zinits)

# Specfies the parameters of interest to record in output
parameters.to.save <- c("phi","p","z")

#MCMC model details
n.iter <- 1000
n.burnin <- 200

myModel <- nimbleModel(code = hmm.model,        
                               constants = my.constants,      
                               data = my.data,
                               inits = initial.values)
                               #monitors = parameters.to.save)
                               #niter = n.iter,
                               #nburnin = n.burnin)
                               #nchains = n.chains,
                               
CmyModel <- compileNimble(myModel)

myMCMC <- buildMCMC(CmyModel, monitors = parameters.to.save, enableWAIC = TRUE)

CmyMCMC <- compileNimble(myMCMC)

### This is the Nimble function that normally produces the progress bar. ###
results <- runMCMC(CmyMCMC, niter = n.iter, nburnin = n.burnin, WAIC = TRUE)
############################################################################

return(results)

}

# Function works correctly for multiple chains with parLapply
system.time({
chain_output <- parLapply(cl = this_cluster, index,
                          fun = paraChain_code,
                          y = ymat)
})

# The pblapply function includes a progress bar, but only for completed chains/total chains.
system.time({
  chain_output <- pblapply(X = index, cl = this_cluster, 
                            FUN = paraChain_code,
                            y = ymat)
 })

stopCluster(this_cluster)
Roasty247
  • 455
  • 3
  • 14

0 Answers0