CODEWITHSUNDEEP
rglm
Fomite
Fomite

Reputation: 2273

How to get an estimate and confidence interval for a contrast in R with offset

I've got a Poisson GLM model fitted in R, looking something like this:

glm(Outcome~Exposure + Var1 + offset(log(persontime)),family=poisson,data=G))

Where Outcome will end up being a rate, the Exposure is a continuous variable, and Var1 is a factor with three levels.

It's easy enough from the output of that:

    Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)   -5.6998     0.1963 -29.029  < 2e-16
Exposure       4.7482     1.0793   4.399 1.09e-05
Var1Thing1    -0.2930     0.2008  -1.459 0.144524    
Var1Thin       1.0395     0.2037   5.103 3.34e-07
Var1Thing3     0.7722     0.2201   3.508 0.000451

To get the estimate of a one-unit increase in Exposure. But a one-unit increase isn't actually particularly meaningful. An increase of 0.025 is actually far more likely. Getting an estimate for that isn't particularly difficult either, but I'd like a confidence interval along with the estimate. My intuition is that I need to use the contrast package, but the following generated an error:

diff <- contrast(Fit,list(Exposure=0.030,Var1="Thing1"),list(Exposure=0.005,Type="Thing1"))

"Error in offset(log(persontime)) : object 'persontime' not found"

Any idea what I'm doing wrong?

Upvotes: 1

Views: 2444

Answers (1)

Jthorpe
Jthorpe

Reputation: 10167

you want to use the confint function (which in this case will call the MASS:::confint.glm method), as in:

confint(Fit)

Since the standard errors is the model scale linearly with the linear changes in the scale of the variable 'Exposure' in your model, you can simply multiply the confidence interval by the difference in scale to get the confidence for a smaller 'unit' change.

Dumb example:

Lets say you want to test the hypothesis that people fall down more often when they've had more alcohol. You test this by randomly serving individuals varying amounts of alcohol (which you measure in ml) and counting the number of times each person falls down. Your model is:

Fit <- glm(falls ~ alcohol_ml,data=myData, family=poisson)

and the coef table is 

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept)   -5.6998     0.1963 -29.029  < 2e-16
Alcohol_ml     4.7482     1.0793   4.399 1.09e-05

and the Confidence interval for alcohol is 4-6 (just to keep is simple). Now a colegue asks you to give the the confidence interval in ounces. All you have to do is to scale by the confidence interval by the conversion factor (29.5735 ounces per ml), as in:

c(4,6) * 29.5735 # effect per ounce alcohol [notice multiplication is used to rescale here]

alternatively you could re-scale your data and re-fit the model:

mydata$alcohol_oz <- mydata$alcohol_ml / 29.5735 #[notice division is used to rescale here]
Fit <- glm(falls ~ alcohol_oz,data=myData, family=poisson)

or you could re-scale your data right in the model:

#[again notice that division is used here]
Fit <- glm(falls ~ I(alcohol_ml/29.5735),data=myData, family=poisson)

Either way, you will get the same confidence intervals on the new scale.

Back to your example: if you're units of Exposure are so large that you are unlikely to observe such a change within an individual and a smaller change is more easily interpreted, just re-scale your variable 'Exposure' (as in myData$Exposure_newScale = myData$Exposure / 0.030 so Exposure_newScale is in multiples of 0.030) or rescale the confidence intervals using either of these methods.

Upvotes: 1

Related Questions