Почему использование "mgcv::s" в "gam(y ~ mgcv::s...)" приводит к ошибке?

Я хотел быть ясным и использовать :: обозначения в строках для подгонки mgcv::gam, Я наткнулся на одну вещь при использовании обозначения в вызове модели для mgcv::s, Код с воспроизводимым примером / ошибкой показан ниже.

Причина, вероятно, в том, что я использую это обозначение в формуле модели, но я не мог понять, почему это не работает / не допускается. Вероятно, это что-то довольно специфическое в отношении синтаксиса (возможно, не в отношении mgcv), но, возможно, кто-то может помочь мне понять это и мое понимание R. Заранее спасибо.

library(mgcv)
dat <- data.frame(x = 1:10, y = 101:110)
# this results in an error: invalid type (list)...
mgcv::gam(y ~ mgcv::s(x, bs = "cs", k = -1), data = dat)
# after removing the mgcv:: in front of s everything works fine
mgcv::gam(y ~ s(x, bs = "cs", k = -1), data = dat)

# outside of the model call, both calls return the desired function
class(s)
# [1] "function"
class(mgcv::s)
# [1] "function"
Источник

2 ответа

Решение

объяснение

library(mgcv)
#Loading required package: nlme
#This is mgcv 1.8-24. For overview type 'help("mgcv-package")'.

f1 <- ~ s(x, bs = 'cr', k = -1)
f2 <- ~ mgcv::s(x, bs = 'cr', k = -1)

OK <- mgcv:::interpret.gam0(f1)$smooth.spec
FAIL <- mgcv:::interpret.gam0(f2)$smooth.spec

str(OK)
# $ :List of 10
#  ..$ term   : chr "x"
#  ..$ bs.dim : num -1
#  ..$ fixed  : logi FALSE
#  ..$ dim    : int 1
#  ..$ p.order: logi NA
#  ..$ by     : chr "NA"
#  ..$ label  : chr "s(x)"
#  ..$ xt     : NULL
#  ..$ id     : NULL
#  ..$ sp     : NULL
#  ..- attr(*, "class")= chr "cr.smooth.spec"

str(FAIL)
# list()

4-я строка исходного кода interpret.gam0 раскрывает проблему:

head(mgcv:::interpret.gam0)

1 function (gf, textra = NULL, extra.special = NULL)              
2 {                                                               
3     p.env <- environment(gf)                                    
4     tf <- terms.formula(gf, specials = c("s", "te", "ti", "t2", 
5         extra.special))                                         
6     terms <- attr(tf, "term.labels")

поскольку "mgcv::s" не должно совпадать, вы получите проблему. Но mgcv позволяет вам комнату обойти это, передавая "mgcv::s" через аргумент extra.special:

FIX <- mgcv:::interpret.gam0(f, extra.special = "mgcv::s")$smooth.spec
all.equal(FIX, OK)
# [1] TRUE

Просто это не контролируется пользователем при выполнении процедуры высокого уровня:

head(mgcv::gam, n = 10)

#1  function (formula, family = gaussian(), data = list(), weights = NULL, 
#2      subset = NULL, na.action, offset = NULL, method = "GCV.Cp",        
#3      optimizer = c("outer", "newton"), control = list(), scale = 0,     
#4      select = FALSE, knots = NULL, sp = NULL, min.sp = NULL, H = NULL,  
#5      gamma = 1, fit = TRUE, paraPen = NULL, G = NULL, in.out = NULL,    
#6      drop.unused.levels = TRUE, drop.intercept = NULL, ...)             
#7  {                                                                      
#8      control <- do.call("gam.control", control)                         
#9      if (is.null(G)) {                                                  
#10         gp <- interpret.gam(formula)  ## <- default to extra.special = NULL

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

Больше понимания:

s, teи т. д. в mgcv не работает в той же логике с stats::poly а также splines::bs,

  • Когда вы делаете, например, X <- splines::bs(x, df = 10, degree = 3)оценивает x и создать матрицу дизайна X непосредственно.
  • Когда вы делаете s(x, bs = 'cr', k = 10)оценка не производится; это анализируется

Гладкая конструкция в mgcv проходит несколько этапов:

  1. парсинг / интерпретация mgcv::interpret.gam, который генерирует профиль для сглаживания;
  2. начальная конструкция mgcv::smooth.construct, которая устанавливает матрицу базиса / дизайна и матрицу штрафов (в основном это делается на уровне C);
  3. вторичное строительство mgcv::smoothCon, который выбирает переменную "by" (например, дублирующее сглаживание для фактора "by"), линейные функциональные члены, штраф за пустое пространство (если вы используете select = TRUE), изменение масштаба штрафа, ограничение центрирования и т. д.;
  4. окончательная интеграция mgcv:::gam.setup, который объединяет все сглаживатели вместе, возвращая матрицу модели и т. д.

Так что это гораздо более сложный процесс.

Источник

Это выглядит как mgcv вопрос. Например, lm() функция принимает poly() или же stats::poly() и дает те же результаты (кроме названий вещей):

> x <- 1:100
> y <- rnorm(100)
> lm(y ~ poly(x, 3))

Call:
lm(formula = y ~ poly(x, 3))

Coefficients:
(Intercept)  poly(x, 3)1  poly(x, 3)2  poly(x, 3)3  
    0.07074      0.13631     -1.52845     -0.93285  

> lm(y ~ stats::poly(x, 3))

Call:
lm(formula = y ~ stats::poly(x, 3))

Coefficients:
       (Intercept)  stats::poly(x, 3)1  stats::poly(x, 3)2  stats::poly(x, 3)3  
           0.07074             0.13631            -1.52845            -0.93285

Это также работает с splines::bs функция, так что это не относится к poly(),

Вы должны связаться с mgcv сопровождающий и указать на эту ошибку в этом пакете. Я предполагаю, что это ищет специально для s, а не для выражения как mgcv::s это оценивает одно и то же.

Источник
Другие вопросы по тегам