使用装有 poly(..., 2) 的 nlme::gls 二次模型进行预测的问题

Question

我已经拟合了一个具有方差结构的二次模型，该结构允许每个因子水平有不同的方差水平，但我在预测只有 2 个条目的新数据集时遇到了问题。这是一个可重现的例子：

library(nlme)
set.seed(101)
mtcars$amf <- factor(mtcars$am)
modGLS <- gls(mpg ~ amf*poly(hp, 2),
           weights = varIdent(form = ~ 1|amf), data = mtcars)
minhp <- min(mtcars$hp); maxhp <- max(mtcars$hp)
newdata <- data.frame(amf = as.factor(c(0, 1)),
                        hp = round(runif(2, min = minhp, max = maxhp)))
newdata2 <- data.frame(amf = as.factor(c(0, 0, 1)),
                        hp = round(runif(3, min = minhp, max = maxhp)))
predict(modGLS, newdata = newdata)
# Error in poly(hp, 2) : 'degree' must be less than number of unique points

predict(modGLS, newdata = newdata2)
## [1]  5.973306 13.758955 44.037921
## attr(,"label")
## [1] "Predicted values"

但是，预测在 lm 框架上运行良好：

modLM <- lm(mpg ~ amf*poly(hp, 2), mtcars)
predict(modLM, newdata = newdata)
##        1        2 
## 25.22253 16.83943

为什么会这样？ emmeans 的包维护者之一似乎认为这可能与 attr(, “predvars”) 上缺少信息有关（请在此处查看我们的讨论 https://github.com/rvlenth/emmeans/issues/133）

我已将此事报告给 Bates 博士（nlme 联系人），但我想我也会联系更广泛的社区。提前致谢

Answer 1

感谢@BenBolker 和@russ-lenth 确认问题与 GLS 对象中缺少的 terms 属性 "predvars" 有关，它提供了 [=13= 的拟合系数].请注意这是如何在 LM 框架（原始 post）中工作的，并且属性在那里（另请参见 ?makepredictcall）。请注意，这可能对预测有潜在影响。

## e.g. this fails
poly(newdata$hp, 2)
## but this is okay, because the polynomial has been estimated
polyFit <- poly(mtcars$hp, 2)
predict(polyFit, newdata = newdata$hp)

attr(terms(modLM), "predvars")
## list(mpg, amf, poly(hp, 2, coefs = list(alpha = c(146.6875, 198.071514938048), norm2 = c(1, 32, 145726.875, 977168457.086594))))

attr(terms(modGLS), "predvars")
## NULL

使用装有 poly(..., 2) 的 nlme::gls 二次模型进行预测的问题

Issues predicting with `nlme::gls`quadratic model fitted with `poly(..., 2)`

r

nlme

poly