### ------------------------------------------------------------------- ###
###                Load INLA library and data file                      ###
### ------------------------------------------------------------------- ###

library(INLA)
setwd("C:\\Users\\Shawn\\Dropbox\\current_projects\\wolf_ddhs_applied_1\\Manuscript")
dat = read.csv("data.csv")

### ------------------------------------------------------------------- ###
###   Function to run hierarchical model in INLA with default priors    ###
### ------------------------------------------------------------------- ###

run.inla = function(form, data, link) {
  mod = inla(form, data=data, family="binomial", control.compute=list(waic=TRUE, cpo=TRUE),
             control.predictor=list(compute=TRUE, link=link))
  return(mod)
}

### ------------------------------------------------------------------- ###
###                      Model formula                                  ###
### ------------------------------------------------------------------- ###

form.all.inla = Used ~ 
  # Main effects
  (Buck + PPC1 + Stream + Elevation + I(Elevation^2) + 
     Slope + I(Slope^2) + Edge + Road + PLand + Impervious)*Density +
  # Random effects: Year
  f(BioYear, model="iid") + 
  # Random effects: Pack
  f(Pack, model="iid") + 
  f(Pack1, PPC1, model="iid") +
  f(Pack2, Buck, model="iid") + 
  f(Pack3, Stream, model="iid") + 
  f(Pack4, Elevation, model="iid") +   
  f(Pack5, Slope, model="iid") + 
  f(Pack6, Edge, model="iid") +
  f(Pack7, Road, model="iid") +  
  f(Pack8, PLand, model="iid") + 
  f(Pack9, Impervious, model="iid") + 
  f(Pack10, I(Elevation^2), model="iid") + 
  f(Pack11, I(Slope^2), model="iid")

### ------------------------------------------------------------------- ###
###                      Run the Model                                  ###
### ------------------------------------------------------------------- ###

link = rep(1, nrow(dat))
Model.Out = run.inla(form.all.inla, dat, link)