Skip to content
Prev 16063 / 20628 Next

Plotting partial residuals from a glmmADMB model

Here are some examples of partial residual plots ... the only fanciness
available in the base-R machinery [e.g. for partial residuals for lm()
objects] that can't be implemented with the simple machinery here is
collapsing entire terms (due to e.g. factors with more than two levels,
orthogonal polynomials, splines ...) into a single effect

library(glmmADMB)
## ran without FoodTreatment, had trouble with convergence ...
owlfit <- glmmadmb(SiblingNegotiation~ArrivalTime+SexParent+
                       (1|Nest),family="nbinom1",data=Owls)

X <- model.matrix(~ArrivalTime+SexParent,data=Owls)

beta <- fixef(owlfit)

## <https://en.wikipedia.org/wiki/Partial_residual_plot> gives the
## definition of a partial residual as (residuals + \hat beta_i X_i)

## multiply each column of X by the corresponding beta
beta_X <- sweep(X,MARGIN=2,STATS=beta,FUN="*")

## add residuals (columnwise addition works automatically,
## although we could use sweep with MARGIN=1)
p_resid <- sweep(beta_X,MARGIN=1,STATS=residuals(owlfit),FUN="+")

par(mfrow=c(1,2))
for (i in 2:3) {
    plot(X[,i],p_resid[,i],xlab=colnames(X)[i],ylab="partial residuals")
}

(Would be better to plot partial residuals for the binary predictor
as a boxplot ...)

If you were doing this in glmmTMB you can take a couple of shortcuts:

library(glmmTMB)
owlfit2 <- glmmTMB(SiblingNegotiation~ArrivalTime+SexParent+FoodTreatment+
                       (1|Nest),family="nbinom1",data=Owls)

X <- getME(owlfit2,"X")
beta <- fixef(owlfit2)$cond
beta_X <- sweep(X,MARGIN=2,STATS=beta,FUN="*")
p_resid <- sweep(beta_X,MARGIN=1,STATS=residuals(owlfit2),FUN="+")
par(mfrow=c(1,3))
for (i in 2:4) {
    plot(X[,i],p_resid[,i],xlab=colnames(X)[i],ylab="partial residuals")
}