В R, как можно извлечь матрицу или значения шляпы / проекции / влияния из объекта модели `nls`?

За lm или же glm тип объектов или даже lmer тип объектов, вы можете извлечь значения шляпы из модели с помощью функции R hatvalues(), Тем не менее, это не работает с nls объекты, по-видимому. Я гуглил в разные стороны, но не могу найти способ получить эти значения. Есть ли nls просто не создавайте матрицу шляп, или значения шляп, полученные из нелинейной модели наименьших квадратов, просто ненадежны?

Воспроизводимый пример:

xs = rep(1:10, times = 10)
ys = 3 + 2*exp(-0.5*xs)
for (i in 1:100) {
xs[i] = rnorm(1, xs[i], 2)
}
df1 = data.frame(xs, ys)
nls1 = nls(ys ~ a + b*exp(d*xs), data=df1, start=c(a=3, b=2, d=-0.5))

1 ответ

Есть хорошая статья (Об обнаружении выбросов в нелинейной регрессии), в которой матрица шляп аппроксимируется градиентной матрицей, вычисленной в расчетной точке.

В твоем случае:

# gradient of the model function at the current parameter values
V <- nls1$m$gradient()

# tangent plane leverage matrix (it plays a similar role as the Hat matrix)
H <- V %*% solve(t(V) %*% V) %*% t(V)

# 'hat' values for nls
nls1.hat_values <- diag(H)

И если вы будете следовать этой статье, вы можете рассчитать H немного быстрее:

Q1 <- qr.Q(qr(V)) # V is the same matrix as above
H <- Q1 %*% t(Q1)

поскольку H может быть довольно большим, и если вам нужны только значения хет, вы можете вообще пропустить умножение матриц. Нам нужна только диагональ H матрица.

###
#' Approximation of hat values for nls.
#'
#' @param model An 'nls' object
#' @param ... Additional parameters (ignored)
#' @return Vector of approximated hat values
###
hatvalues.nls <- function(model, ...) {
  stopifnot(is(model, 'nls'))
  list(...) # ignore additional parameters
  V <- model$m$gradient()
  Q1 <- qr.Q(qr(V))
  rowSums(Q1*Q1)
}
Другие вопросы по тегам