Plotting predicted values in GLM in R I am trying to find a more aesthetic way to present an interaction with a quadratic term in a logistic regression (categorisation of continuous variable is not appropriate).
For a simpler example I use a linear term.
set.seed(1)

df<-data.frame(y=factor(rbinom(50,1,0.5)),var1=rnorm(50),var2=factor(rbinom(50,1,0.5)))
mod<-glm(y ~ var2*var1  , family="binomial" , df)

 #plot of predicted probabilities of two levels

new.df<-with(df,data.frame(expand.grid(var1=seq(-2,3,by=0.01),var2=levels(var2))))
pred<-predict(mod,new.df,se.fit=T,type="r")

with(new.df,plot(var1,pred$fit))

 #plot the difference in predicted probabilities

trans.logit<-function(x) exp(x)/(1+exp(x))

pp<-trans.logit(coef(mod)[1] + seq(-2,3,by=0.01) * coef(mod)[3]) -trans.logit((coef(mod)[1]+coef(mod)[2]) + seq(-2,3,by=0.01) * (coef(mod)[3]+coef(mod)[4]))

plot(seq(-2,3,by=0.01),pp)

Questions


*

*How can I plot the predicted probability difference between the two levels of var2 (rather than the 2 levels separately)  at different values of var1?

*Is there a way to define contrasts so I can use these in the glm so I can then pass this to predict? - I need a CI for the difference in probabilities

 A: The rms package has a general contrast.rms function that also works with the glht function in the multcomp package to give simultaneous confidence intervals.  This works for log odds ratios (and hence odds ratios).  I will investigate whether it is easy to add an option to get bootstrap confidence intervals for differences in probabilities.  We generally use the odds ratio scale because odds ratios can be independent of the settings of other variables in the model.  Risk differences must be considered for a whole array of settings of the other covariates in the model.
A: Here are a few options:
You can use the glht function in the multcomp package for R and specify your own contrasts/comparisons.  The result can be used with the confint function to compute the confidence intervals.
You can use the Predict.Plot function in the TeachingDemos package for R (and the related TkPredict function) to create plots that will demonstrate how the predictions change with the variables.
You fit the model using Bayesian methods and MCMC, then you just do the calculation that you want to get the posterior distribution of the combination of interest and plot that, or the intervals based on them.
