drnickisaac / sMon-Amphibia

Analysis of trends in German amphibia
3 stars 2 forks source link

Convergence code #2

Closed drnickisaac closed 6 years ago

drnickisaac commented 6 years ago

share code to check convergence with @bowlerbear

bowlerbear commented 6 years ago

I have done this now. It was pretty straightforward to pull out each chain.

The function called 'plotModels' in now at the bottom of the sparta_wrapper_function.r script in the R folder. I paste it below too:

plotModels<-function(models,param="psi.fs\["){

create a directory to put the traceplots

newdir <- file.path("model-outputs", "traceplots") dir.create(newdir,showWarnings = FALSE) if (!grepl(newdir,getwd())){ setwd(newdir) }

remove any previous plotting files, if there are any

if(length(list.files())>0){ file.remove(list.files()) }

for each species, do the following

lapply(models,function(x){ require(reshape2)

get the chains

df <- melt(x$BUGSoutput$sims.array)
#give sensible names
names(df)<-c("Index","Chain","Parameter","Estimate")
#subset to occupancy parameter
df <- subset(df,grepl(param,df$Parameter))

#if there are many parameters (i.e. years), just plot every other year
if(length(unique(df$Parameter))>20){
  df$ParamNu <- as.numeric(sub(".*\\[([^][]+)].*", "\\1", df$Parameter))
  df <- df[df$ParamNu%%2==0,]
}

#plot it
require(ggplot2)
ggplot(df)+
  geom_path(aes(x=Index,y=Estimate,colour=factor(Chain)))+
  facet_wrap(~Parameter,scales="free",ncol=4)+
  theme_bw()+
  theme(legend.position="none")
ggsave(paste0(x$SPP_NAME,"_traceplot.png"))

})

}