# Plotting confidence bands around fitted values from a binomial GLMM

I have some parameter estimates and confidence intervals estimated from a set of model-averaged binomial GLMMs: two main effects and their interaction. I would like to plot [population level] fitted values for the main effects and their interaction, as well as 95% confidence bands. I can plot the population level fitted values, but I'm not sure how to combine the parameter estimates with their CI's to plot confidence bands.

Code for plotting the fitted values, and resulting graph, is below. I'm plotting rat caputre probability within forest patches (Rat) against distance from forest edge (Dist) and whether the forest patch was grazed or not by livestock (Grazed). Distance was square-root transformed in the model (Dist_T).

Could someone please tell me how to combine the parameter estimates with their CI's to plot confidence bands?

Parameter estimates from the model were:

Intercept: -1.557; CI=-3.034, -0.08

Dist_T: 0.097; CI=-0.017, 0.212

GrazedY: -1.898; CI=-4.17, 0.375

Dist_T:GrazedY: 0.182; CI=-0.038, 0.401

# estimated abund where grazing=N:
PlotData <- data.frame(Dist_T=seq(from=min(Rat$Dist_T), to=max(Rat$Dist_T), by = 0.1))
g <- -1.557 + 0.097*PlotData$Dist_T Rat.fitted <- exp(g)/(1+exp(g)) Dist.bt <- PlotData$Dist_T^2
plot(Rat$Dist, jitter(Rat$Rat, amount=0.02), ylim=c(0,1), cex=0.7, cex.lab=0.7, cex.main=0.4, ylab="Rat capture probability",
xlab="Distance from forest edge (m)", pch=ifelse(Rat$Grazed=="Y", 1,2), col=ifelse(Rat$Grazed=="Y", "grey","black"))
lines(Dist.bt,Rat.fitted, lty=1, col="black")

# estimated abund where grazing=Y:
par(new=T)
PlotData <- data.frame(Dist_T=seq(from=min(Rat$Dist_T), to=max(Rat$Dist_T[Rat$Grazed=="Y"]), by = 0.1))#only plotting within range; where Grazed=Y g <- -1.557 + (0.097+0.182)*PlotData$Dist_T -1.898*1 #additional grazing term and interaction term
Rat.fitted <- exp(g)/(1+exp(g))
Dist.bt <- PlotData$Dist_T^2 plot(Rat$Dist, jitter(Rat\$Rat, amount=0.05), ylim=c(0,1), xlab="", ylab="", type="n")
lines(Dist.bt,Rat.fitted, lty=1, col="grey")