How to extract posterior data from output_posteriors object?

[naalipalo]

I have similar data to the alligator example with a continuous variable (mine is weight) but my data is significantly different by location (6 sites). The output settings work to produce the graphs I want, but I need to modify the graphs because they change the color of each prey item by location and the titles are wrong etc.... Anyway, I am not certain what data I need to extract and where.
Ive been following the code for modifying MixSIAR plots, http://brianstock.github.io/MixSIAR/articles/modify_output.html#set-output-options-to-return-ggplot-objects-1

but it still doesn't have the flexibility I know how to manipulate, so would like to know what to extract from the output_posteriors object created so I can make my own ggplot of the data.

[hthalmann]

I am having this same issue. In the output options for my model with one fixed categorical variable (Year) and one continuous variable (Age), I can get plots of proportion by age for each year, but I would like to be able to modify these in ggplot. However, when I extract the code using a similar set of functions as the alligator code, I get the age and proportion plot for all years combined. How would I extract posterior data for each level of my categorical variable?

Here is my code to replicate extracting the posteriors for the continuous variable with all levels of my categorical variable combined:

R2jags::attach.jags(jags.2)
n.sources <- source$n.sources
source_names <- source$source_names

calc_eps <- function(f){
n.sources <- length(f)
gam <- rep(1/n.sources,n.sources)
phi <- rep(0,n.sources)
phi[1] <- 1
sqrt(sum((f-gam)^2))/sqrt(sum((phi-gam)^2))
}

ce=1
label <- mix$cont_effects[ce]
cont <- mix$CE[[ce]]
ilr.cont <- get(paste("ilr.cont",ce,sep=""))

n.plot = 200
chain.len = dim(p.global)[1]
Cont1.plot <- seq(from=round(min(cont),1), to=round(max(cont),1), length.out=n.plot)
ilr.plot <- array(NA,dim=c(n.plot, n.sources-1, chain.len))
for(src in 1:n.sources-1){
for(i in 1:n.plot){
ilr.plot[i,src,] <- ilr.global[,src] + ilr.cont[,src]*Cont1.plot[i]
}
}

e <- matrix(rep(0,n.sources*(n.sources-1)),nrow=n.sources,ncol=(n.sources-1))
for(i in 1:(n.sources-1)){
e[,i] <- exp(c(rep(sqrt(1/(i*(i+1))),i),-sqrt(i/(i+1)),rep(0,n.sources-i-1)))
e[,i] <- e[,i]/sum(e[,i])
}

cross <- array(data=NA,dim=c(n.plot, chain.len, n.sources, n.sources-1))
tmp <- array(data=NA,dim=c(n.plot, chain.len, n.sources))
p.plot <- array(data=NA,dim=c(n.plot, chain.len, n.sources))
for(i in 1:n.plot){
for(d in 1:chain.len){
for(j in 1:(n.sources-1)){
cross[i,d,,j] <- (e[,j]^ilr.plot[i,j,d])/sum(e[,j]^ilr.plot[i,j,d]);
}
for(src in 1:n.sources){
tmp[i,d,src] <- prod(cross[i,d,src,]);
}
for(src in 1:n.sources){
p.plot[i,d,src] <- tmp[i,d,src]/sum(tmp[i,d,]);
}
}
}

get_high <- function(x){return(quantile(x,.95))} # 90% CI
get_low <- function(x){return(quantile(x,.05))}
p.low <- apply(p.plot, c(1,3), get_low)
p.high <- apply(p.plot, c(1,3), get_high)
p.median <- apply(p.plot, c(1,3), median)
colnames(p.median) <- source_names

Cont1.plot <- Cont1.plot*mix$CE_scale + mix$CE_center # transform Cont1.plot (x-axis) back to the original scale
df <- data.frame(reshape2::melt(p.median)[,2:3], rep(Cont1.plot,n.sources), reshape2::melt(p.low)[,3], reshape2::melt(p.high)[,3])
colnames(df) <- c("source","median","ifbf.per","low","high")