Set the offset value in lsmeans

Question

I'm using a negative binomial glmm for analysis.

The y values of my samples are summed over 12 individual counts of insect mines on leaves of a tree. You can get the data (and session info) from here!.

Treat has 2 levels; Year has 2 levels; Habitat has 3 levels; leaves is the same number of 12 leaves for each sample; location has 4 levels and is almost similar to location. Two of the locations fall in the same "habitat" category. I use sublocation as random effect, since I want to get results for habitat not for all sublocations separately while being aware of the more similar measures at a single location.

This is the model:

library(MASS)
library(lme4)
library(optimx)

glmm1 <- glmer.nb(response ~ year*(Habitat+Treat) + offset(log(leaves)) + (1|location), data = df1, control = glmerControl(optimizer = "optimx", calc.derivs = FALSE, optCtrl = list(method = "nlminb", starttests = FALSE, kkt = FALSE)))

Now I use the following code to create a prediction plot:

# Prediction plot
predframe <- expand.grid(year= levels(df1$year),
                         Habitat= levels(df1$Habitat),
                         Treat= levels(KMM2$Treat),
                         leaves = 10) # Set the offset!!

Here I am able to alter the offset to get results for the more commonly used reference of 10 leaves !

mm <- model.matrix(~ year*(Habitat + Treat) + offset(log(leaves)), data = predframe)

predframe$fit <- predict(glmm1,newdata=predframe, re.form = NA) # Don't use type="response"!!!
# Alternatively use: predframe$fit <- mm %*% fixef(glmm1)

pvar1 <- diag(mm %*% tcrossprod(vcov(glmm1),mm))
# alternatively use: pvar1 <- diag(mm %*%vcov(glmm1) %*% t(mm)) 
tvar1 <- pvar1+VarCorr(glmm1)$Standort[1] # must be adapted for more complex models
predframe <- data.frame(predframe
                        , p_lwr = predframe$fit-1.96*sqrt(pvar1) # for confidence interval
                        , p_upr = predframe$fit+1.96*sqrt(pvar1)
                        , t_lwr = predframe$fit-1.96*sqrt(tvar1) # for prediciton interval
                        , t_upr = predframe$fit+1.96*sqrt(tvar1)
)

predframe2 <- aggregate(cbind(fit, p_lwr, p_upr) ~ year + Treat,predframe, mean)

ggplot(aes(y=interaction(year, Treat), x=exp(fit)), data=predframe2) + geom_point() +
        geom_errorbarh(aes(xmax=exp(p_upr), xmin=exp(p_lwr)))

Here I use the lsmeans package to do the same:

library(lsmeans)
lsm.options(save.ref.grid = TRUE)
lsm <- lsmeans(glmm1, ~ Netz*Jahr, offset=10)
.Last.ref.grid

Offset was not adjusted to 10!

plot(cld(lsm, type="response"))

Besides slightly different confidence intervals (why?) these results are the same as if I used leaves=12 for the predframe above

How can I adjust the offset in lsmeans?

Thanks in advance

Answer 1

This is how it works, thanks to rvl!!

lsm.options(save.ref.grid = TRUE)
lsm <- lsmeans(glmm1, ~ Treat*year, at=list(leaves=10))
.Last.ref.grid
plot(cld(lsm, type="response", sort=FALSE))

Now the reference grid is similar to the prediction frame above.

However: same same but different! Here are the two plots from both different procedures. I PREFER THE ONE FROM lsmeans!!!!!!