[R-sig-ME] Plotting MCMCglmm regression with credible intervals
jonnations
jonnations at gmail.com
Fri Feb 9 21:26:02 CET 2018
Hi all,
I am working on a logistic regression model using MCMCglmm. I would like to
plot the mean with 95% credible intervals on my regression plot, however
this is not so straightforward as the MCMCglmm output is in a very
different format than most software. I think I can manually bang it out
with a page of code, but I was curious if anyone knows of a simple, elegant
way to do this.
Below is a sample script. Imagine this plot with 2 additional lines (95%
CI's) on either side of the mean, shaded in between.
If r-sig-mixed-models is not the correct venue for this kind of question, I
apologize.
Best,
Jon
library(MCMCglmm)
#generate binary data
intercept1 = -6.0
b = 2.75
x = rnorm(100, 1, 3)
linpred = intercept1 + x * b
prob1 = exp(linpred) / (1 + exp(linpred))
run1 = runif(100, 0, 1)
y = ifelse(run1 < prob1, 1, 0)
df <- data.frame(y = y, x = x)
model1<-MCMCglmm(y ~ x, data = df, family = "categorical",
verbose = F, nitt =2000, burnin = 500, thin = 1)
x <- df$x
y <- df$y
mod.x <- mean(model1$Sol[,1])
mod.y <- mean(model1$Sol[,2])
plot(x,y)
curve(plogis(mod.x+mod.y*x),col="red",add=TRUE)
--
Jonathan A. Nations
PhD Candidate
Esselstyn Lab <http://www.museum.lsu.edu/esselstyn>
Museum of Natural Sciences <http://sites01.lsu.edu/wp/mns>
Louisiana State University
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list