目录：
- 一个简单的演练示例
- 给用户的建议
- 我们可以从拟合模型对象中获得的有用信息
- 好的，我知道现在是什么问题，但是如何让
predict 工作？
- 有没有更好的方法可以完全避免此类问题？
一个简单的演练示例
这里有一个足够简单的可重复的例子来提示你发生了什么。
train <- data.frame(y = runif(4), x = c(runif(3), NA), f = factor(letters[1:4]))
test <- data.frame(y = runif(4), x = runif(4), f = factor(letters[1:4]))
fit <- lm(y ~ x + f, data = train)
predict(fit, newdata = test)
#Error in model.frame.default(Terms, newdata, na.action = na.action, xlev = object$xlevels) :
# factor f has new levels d
我正在拟合一个参数多于数据的模型，因此该模型排名不足（将在最后解释）。但是，这不会影响 lm 和 predict 的工作方式。
如果您只检查table(train$f) 和table(test$f)，它没有用，因为问题不是由变量f 引起的，而是由x 中的NA 引起的。 lm 和 glm 删除不完整的案例，即至少有一个 NA 的行（参见 ?complete.cases）用于模型拟合。他们必须这样做，否则 QR 分解的底层 FORTRAN 例程将失败，因为它无法处理 NA。如果你查看?lm 的文档，你会看到这个函数有一个参数na.action，默认为na.omit。您也可以将其设置为na.exclude，但na.pass 保留NA 会导致FORTRAN 错误：
fit <- lm(y ~ x + f, data = train, na.action = na.pass)
#Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) :
# NA/NaN/Inf in 'x'
让我们从训练数据集中删除NA。
train <- na.omit(train)
train$f
#[1] a b c
#Levels: a b c d
f 现在有一个未使用的级别 "d"。 lm 和 glm 将在构建模型框架（以及后来的模型矩阵）时丢弃未使用的级别：
## source code of lm; don't run
mf$drop.unused.levels <- TRUE
mf[[1L]] <- quote(stats::model.frame)
mf <- eval(mf, parent.frame())
这是用户无法控制的。原因是如果包含一个未使用的级别，它将在模型矩阵中生成一列零。
mf <- model.frame(y ~ x + f, data = train, drop.unused.levels = FALSE)
model.matrix(y ~ x + f, data = mf)
# (Intercept) x fb fc fd
#1 1 0.90021178 0 0 0
#2 1 0.10188534 1 0 0
#3 1 0.05881954 0 1 0
#attr(,"assign")
#[1] 0 1 2 2 2
#attr(,"contrasts")
#attr(,"contrasts")$f
#[1] "contr.treatment"
这是不希望的，因为它会为虚拟变量 fd 生成 NA 系数。由drop.unused.levels = TRUE 强制由lm 和glm：
mf <- model.frame(y ~ x + f, data = train, drop.unused.levels = TRUE)
model.matrix(y ~ x + f, data = mf)
# (Intercept) x fb fc
#1 1 0.90021178 0 0
#2 1 0.10188534 1 0
#3 1 0.05881954 0 1
#attr(,"assign")
#[1] 0 1 2 2
#attr(,"contrasts")
#attr(,"contrasts")$f
#[1] "contr.treatment"
fd 不见了，并且
mf$f
#[1] a b c
#Levels: a b c
现在不存在的"d" 级别将导致predict 中出现“新因子级别”错误。
给用户的建议
强烈建议所有用户在拟合模型时手动执行以下操作：
-
[没有。 1] 删除不完整的案例；
-
[没有。 2] 降低未使用的因子水平。
这正是这里推荐的过程：How to debug "contrasts can be applied only to factors with 2 or more levels" error? 这让用户了解lm 和glm 在后台做了什么，并使他们的调试工作更加轻松。
注意，列表中应该还有其他推荐：
用户可能偶尔会使用subset 参数。但是有一个潜在的陷阱：不是所有的因子水平都可能出现在子集数据集中，因此您以后使用predict 时可能会得到“新的因子水平”。
当您编写包装lm 或glm 的函数时，上述建议尤其重要。你希望你的函数是健壮的。要求您的函数返回一个信息性错误，而不是等待 lm 和 glm 抱怨。
我们可以从拟合模型对象中获得的有用信息
lm 和 glm 在拟合对象中返回 xlevels 值。它包含实际用于模型拟合的因子水平。
fit$xlevels
#$f
#[1] "a" "b" "c"
因此，如果您没有遵循上面列出的建议并且在因子水平方面遇到问题，那么应该首先检查这个xlevels。
如果您想使用 table 之类的东西来计算每个因子水平有多少个案例，这里有一种方法：Get number of data in each factor level (as well as interaction) from a fitted lm or glm [R]，尽管制作模型矩阵可能会占用大量 RAM。
好的，我知道现在是什么问题，但是如何让predict 工作？
如果您不能选择使用不同的 train 和 test 数据集（请参阅下一节），您需要在 test 而不是 xlevels 中将这些因子水平设置为 @ 987654383@。然后predict 将只预测NA 用于此类不完整的情况。
有没有更好的方法来避免此类问题？
人们将数据拆分为train 和test，因为他们想要进行交叉验证。第一步是在您的完整数据集上应用na.omit 以消除NA 噪声。然后我们可以对剩下的内容进行随机分区，但是这种天真的方式可能会以
-
test 中的某些因子水平，但 train 中没有
-
train 中的某些因子变量在移除未使用的级别后只有 1 个级别（糟糕，使用 lm 和 glm 时出现“对比”错误）；
因此，强烈建议您进行一些更复杂的分区，例如分层抽样。
其实还有一个危险，但不会导致编程错误：
-
train 的模型矩阵秩不足（糟糕，使用predict 时，我们收到“秩不足模型的预测可能具有误导性”警告）。
关于模型拟合中的rank-deficiency，请参阅lme4::lmer reports "fixed-effect model matrix is rank deficient", do I need a fix and how to?rank-deficiency不会对模型估计和检查造成问题，但可能会对预测造成危害：R lm, Could anyone give me an example of the misleading case on “prediction from a rank-deficient”?但是，这样的问题更难避免，尤其是当你有很多因素并且可能与交互时。