r - plot the results glm with multiple explanatories with 95% CIs -


data:

df <- structure(list(x = c(9.5638945103927, 13.7767187698566, 6.0019477258207,  10.1897072092089, 15.4019854273531, 10.9746646056535, 12.9429073949468,  20.7513493525379, 18.5764146937149, 2.91302077116471, 13.6523222711501,  10.0920467755108), y = c(83.949498880077, 18.066881289085, 71.3052196358606,  39.8975644317452, 57.2933166677927, 87.8484256883889, 92.6818329896141,  49.8297961197214, 56.3650103496898, 14.7950650020996, 37.9271392096266,  50.4357237591891), z = c("a", "c", "e", "f", "b", "a", "b", "a",  "b", "a", "c", "d")), .names = c("x", "y", "z"), row.names = c(na,  -12l), class = "data.frame") 

my model:

mod <- glm(y ~ x + i(x^2) + z, family=quasipoisson, data = df) summary(mod) 

i want plot this:

ggplot(df, aes(x=x,y=y)) +    geom_point() +   stat_smooth(method="lm",se=false,               formula=y~x+i(x^2),fill="transparent",               colour="black") +   stat_smooth(method="lm",geom="ribbon",               formula=y~x+i(x^2),fill="transparent",               colour="red",linetype="dashed",fullrange=true) +           scale_x_continuous(limits=c(-2,35)) +           coord_cartesian(xlim=c(2,25),                           ylim=range(pretty(df$y)))  

enter image description here

however, model have plotted not same model mod, there no z , not quasiposson.

how can plot ggplot using actual model? have looked @ predict however, lost when there more 1 explanatory, in case. don't care doing in ggplot2

it seems trivially adapt example new model using stat_smooth(method='glm', family=quasipoisson, ...), adding z formula leads errors. looking ggplot2 docs, can see predictdf used generate limits intervals. looking @ code function, looks meant work predictions across x dimension. can write our own version works in multiple dimensions , plot limits separate layers.

mypredictdf <- function (model, newdata, level=0.95){   pred <- stats::predict(model, newdata = newdata, se =true, type = "link")   std <- qnorm(level/2 + 0.5)   data.frame(newdata,              y = model$family$linkinv(as.vector(pred$fit)),              ymin = model$family$linkinv(as.vector(pred$fit - std * pred$se)),              ymax = model$family$linkinv(as.vector(pred$fit + std * pred$se)),               se = as.vector(pred$se)) } px <- with(df, seq(from=min(x), to=max(x), length=100)) pdf <- expand.grid(x=px, z=unique(df$z)) pdf <- mypredictdf(mod, newdata=pdf) g <- ggplot(data=pdf, aes(group=z)) g <- g + geom_point(data=df, aes(x=x, y=y, color=z)) g <- g + geom_ribbon(aes(x=x, ymin=ymin, ymax=ymax),                      alpha=0.2) g <- g + geom_line(aes(x=x, y=y, color=z)) 

one panel

it seems faceting idea:

g <- g + facet_wrap(~z) 

facetted version


Comments

Popular posts from this blog

c# - How to get the current UAC mode -

postgresql - Lazarus + Postgres: incomplete startup packet -

javascript - Ajax jqXHR.status==0 fix error -