我想清楚一点,并在各行中使用::
表示法来适合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"
最佳答案
说明
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()
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
中的平滑构造分几个阶段:mgcv::interpret.gam
解析/解释,生成更平滑的配置文件;mgcv::smooth.construct
的初始构造,建立基础/设计矩阵和惩罚矩阵(主要在C级完成);mgcv::smoothCon
的二次构造,它拾取“ by”变量(例如,对因子“ by”进行平滑复制),线性函数项,空空间罚分(如果使用select = TRUE
),罚分重定比例,居中约束等;mgcv:::gam.setup
的最终集成,它将所有平滑器组合在一起,返回模型矩阵,等等。因此,这是一个复杂得多的过程。