Calculqting leverage and Cook's distance in non-linear mixed effects model

Question

I am trying to calculate leverage and cooks distance for the model that I developed using nlme, but I'm getting the error saying

Error in cooks.distance.lme(model) : not implemented for "nlme" objects"

What are the ways to calculate Cook's distance in nlme?

formula = log_Lobs ~ log(150*((1 + ((150/Lt_1)^(1/exp(p))-1)*exp(-exp(k)*td/365))^(-exp(p))))
model <- do.call(nlme,
                    list(formula,
                         fixed = c(p ~ 1, k ~ 1 + season2),
                         random = k ~ 1 | id,
                         data = data_select,
                         start = list(fixed = c(p, k)),
                         na.action = na.exclude,
                         control=list(maxIter=1e6, msMaxIter = 1e6, msVerbose = TRUE)
                    ))

#Calculate leverage
lev<-hat(model.matrix(model))

Error in x$terms %||% attr(x, "terms") %||% stop("no terms component nor attribute") : no terms component nor attribute

#Calculate Cook's Distance
cd<-cooks.distance(model)

Error in cooks.distance.lme(model) : not implemented for "nlme" objects

Answer 1

You could do it by brute force: from Wikipedia,

Cook's distance of observation i is defined as the sum of all the changes in the regression model when observation i is removed from it

Using the formula from Nieuwenhuis 2012

library(nlme)
mod1 <- lme(distance ~ age*Sex, random = ~age|Subject, data = Orthodont)
coef1 <- fixef(mod1)
cooksd <- rep(NA, nrow(Orthodont))
Vinv <- solve(vcov(mod1))
np <- length(coef1)
for (i in seq(nrow(Orthodont))) {
   newmod <- try(update(mod1, data = Orthodont[-i, ]), silent = TRUE)
   if (inherits(newmod, "try-error")) next
   newcoef <- fixef(newmod)
   cdiff <- coef1-newcoef
   cooksd[i] <- (cdiff %*% Vinv %*% cdiff)/np
}

Nieuwenhuis, Rense. 2012. “Influence.ME: Tools for Detecting Influential Data in Mixed Effects Models” R Journal 4/2 https://journal.r-project.org/archive/2012/RJ-2012-011/RJ-2012-011.pdf