0
$\begingroup$

So I've built a (long computational time) jags mixed model that includes an interaction term in the fixed effects. After the fact I want to multiply two of the posterior means for the regression coefficients - with credible intervals.

Example:

library(lme4)
library(runjags)

model <- template.jags(total.fruits ~ reg + amd*gen + (1 | reg), data=Arabidopsis, family = "poisson")

# can ignore warning for this demo

ml.jags <- run.jags(model)
summary(ml.jags)

So I want to work out the value of amd*gen for level 2 of amd afterwards. # For the Mean or Median I can multiply the terms

> results[row.names(results) == "amd_effect[2]", "Mean"] * results[row.names(results) == "gen_coefficient_amd_level[2]", "Mean"]
[1] 0.0001745244

But how do I put credible intervals around this since multiplying the credible intervals directly is wrong? I know in a frequentist framework I can do this using the standard errors but I was not sure if that is true also in Bayesian stats?

$\endgroup$
  • $\begingroup$ Normally when I see interaction models they would be of the form total.fruits ~ reg + amd + gen + amd*gen.... Is there a reason you didn't go this route? Also can you provide a bit more detail about what you're trying to do? I'm not sure I follow, but normally you work with the entire set of posterior samples and then calculate your quantiles of interest directly (e.g. 0.025, 0.5, 97.5). $\endgroup$ – Tdisher Jan 18 '18 at 16:25
  • $\begingroup$ The notation amd*gen in this model implies: amd + gen + amd:gen allowing for main effects of both terms and interaction term. I no longer have the posterior samples only the summary data. I didn't anticipate wanting to do this you see and they were taking up mountains of disk space. $\endgroup$ – user2498193 Jan 18 '18 at 16:36
  • $\begingroup$ Also just as an aside, I would highly recommend you switch to rstanarm. Uses lme4 syntax but is being kept to date by some premier applied Bayesians. I only use JAGS when I have to adapt someone's JAGS code by hand. $\endgroup$ – Tdisher Jan 18 '18 at 16:53
  • 1
    $\begingroup$ Ah OK. Well if you have the posterior standard deviation still and are willing to assuming your posterior was normally distributed than you could use that in place of standard error in whatever calculations you would normally do. It's not strictly correct though. $\endgroup$ – Tdisher Jan 18 '18 at 16:55
  • $\begingroup$ Yeah learning STAN is in my list of things to do. If they release within thread parallelisation via MPI I'll be all over it :) Ok thanks I might do that as a temporary solution.....it will take weeks to rerun everything and I'll have gone insane before then! $\endgroup$ – user2498193 Jan 18 '18 at 17:10
1
$\begingroup$

Error propagation of this sort is easy in MCMC-based Bayesian statistics. Instead of multiplying summary statistics of the posterior (like the means or medians), multiply the ACTUAL MCMC samples from the posterior chains. Since each iteration of the chain is a valid sample from the joint posterior distribution, you can do the multiplication at each iteration of the chain to get a new chain that represents samples from the posterior distribution of the derived parameter. Then, simply conduct your inference on this new chain (e.g. find the median, mode, credible interval, or whatever feature of the posterior you are interested in).

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.