I have a rather large and complicated Poisson GLM with a mix of categorical (many binary, yes-no) and continuous values. To calculate the incident rate ratios (IRR) for the binary variables, I simply exponentiate the fitted coefficient and for the confidence intervals I also exponentiate the lower and upper bounds (I'm doing this in R so I've used the
confint() function, as well as simply
+/- 1.96*SE when I want a quick calculation. They don't align exactly but are very close).
One variable in our model of particular interest now is a continuous variable with value between 0 and 100 where we've essentially manually created a spline within our model. For example, let's call this variable
grade. For this value, we've created new columns with the following logical rules (all values are 0 if logic=
grade0_70_slope=grade, else if
grade70_95_slope=grade, else if
grade95_slope=log(grade/(1-grade)). We are particularly interested in
grade=c(90,95,99). We want to create a plot of the IRRs (with CIs) comparing the IRRs of
grade=c(90,95,99) against a couple of our binary variables. To calculate the IRR point estimate for each
grade of interest, I used
predict.glm on a new dataset where each covariate value is exactly the same except the
grade=c(0,90,95,99). I calculate each fitted value, as well as the
se.fit for each prediction. Using
grade=0 as a baseline, I calculate the IRRs for each
grade of interest by
exp(fitted grade=90|95|99)/exp(fitted grade=0).
After this is where I'm unsure of how to proceed. I'm not sure I know how to accurately calculate the CIs for the
grade IRRs. With my approach I get a
se.fit value for every prediction point, including
grade=0. The "logical" approach is to use my same approach as above but for each lower and upper bound (so I calculate the lower bound for each grade, exponentiate it and then divide it by the exponentiated lower bound for
grade=0)...but my instinct tells me this is not a valid approach.