lme4::lmer报告“固定效应模型矩阵秩不足”，我需要修复吗？如何修复？

Question

我正尝试使用混合效应模型，将其余的列作为预测变量来预测F2_difference，但是我收到了一个错误信息：

fixed-effect model matrix is rank deficient so dropping 7 columns / coefficients.

根据这个链接Fixed-effects model is rank deficient，我认为我应该在R包caret中使用findLinearCombos。然而，当我尝试使用findLinearCombos(data.df)时，它会给我一个错误信息：

Error in qr.default(object) : NA/NaN/Inf in foreign function call (arg 1) In addition: Warning message: In qr.default(object) : NAs introduced by coercion

我的数据中没有任何缺失值，可能是什么原因导致这种情况？（如果答案显而易见，抱歉 - 我是R新手。）我的所有数据都是因子，除了我要预测的数值。这是我的数据的一个小样本。

sex <- c("f", "m", "f", "m")
nasal <- c("TRUE", "TRUE", "FALSE", "FALSE")
vowelLabel <- c("a", "e", "i", "o")
speaker <- c("Jim", "John", "Ben", "Sally")
word_1 <- c("going", "back", "bag", "back")
type <- c("coronal", "coronal", "labial", "velar")
F2_difference <- c(345.6, -765.8, 800, 900.5)
data.df <- data.frame(sex, nasal, vowelLabel, speaker,
                      word_1, type, F2_difference
                      stringsAsFactors = TRUE)

编辑： 这里有更多的代码，如果有帮助的话。

formula <- F2_difference ~ sex + nasal + type + vowelLabel + 
           type * vowelLabel + nasal * type +
           (1|speaker) + (1|word_1)

lmer(formula, REML = FALSE, data = data.df)

编辑编辑：

原帖没有提供足够数量的测试数据，不能让读者实际运行lmer模型。但这不是太大的问题。这仍然是一篇非常好的帖子！

- Lisa

Answer 1

您对警告信息有些过于关注：

固定效应模型矩阵秩不足，因此删除7列/系数。

这只是一个警告而不是错误。既没有误用lmer也没有模型公式的不当规定，因此您将获得一个估计的模型。但为了回答您的问题，我会努力解释它。

---

在执行lmer期间，您的模型公式被分成一个固定效应公式和一个随机效应公式，并且为每个公式构建了一个模型矩阵。 固定效应公式的构建是通过标准模型矩阵构造器model.matrix进行的; 随机效应公式的构建很复杂，但与您的问题无关，因此我将跳过它。

对于您的模型，您可以通过以下方式检查固定效应模型矩阵的外观：

fix.formula <- F2_difference ~ sex + nasal + type + vowelLabel + 
               type * vowelLabel + nasal * type

X <- model.matrix (fix.formula, data.df)

你的所有变量都是因子，所以X将是二元的。虽然model.matrix应用contrasts来处理每个因子和它们的交互作用, 但仍有可能X最终没有完整的列秩，因为一列可能是其他某些列的线性组合（这可能是精确的或数值上接近的）。在你的情况下，一个因子的某些水平可能嵌套在另一个因子的某些水平中。

列秩不足可能以许多不同的方式出现。另一个回答分享了一个CrossValidated答案，提供了大量的讨论，我会做一些评论。

• 对于情况1，人们可以通过LASSO等方法实际进行特征选择模型。
• 情况2和3与数据收集过程有关。良好的实验设计是防止秩缺陷的最佳方法，但对于许多构建模型的人来说，数据已经存在，无法进行改进（如获取更多数据）。然而，我想强调的是，即使没有秩缺陷的数据集，如果我们不小心使用它，仍然可能出现这个问题。例如，交叉验证是一种比较模型的好方法。为此，我们需要将完整数据集分成一个训练集和一个测试集，但如果不小心，我们可能会从训练数据集中得到一个秩缺陷模型。
• 情况4是一个非常棘手的问题，可能完全超出我们的控制范围。也许减少模型复杂度是一个自然的选择，但另一种选择是尝试惩罚回归。
• 情况5是导致数值秩缺陷的数字问题，this是一个很好的例子。
• 情况6和7说明数值计算是以有限精度进行的。通常，如果情况5得到妥善处理，这些问题就不会成为问题。

因此，有时我们可以解决这种缺陷，但并不总是可能实现这一点。因此，任何良好编写的模型拟合例程，如lm，glm，mgcv::gam，都将对X应用QR分解，仅使用其全秩子空间进行估计，即X的列的最大子集，该子集给出一个全秩空间，将与其余列相关的系数修正为0或NA。你收到的警告只是意味着这一点。原本有ncol(X)个系数需要估计，但由于缺陷，只会估计ncol(X) - 7个，其余的为0或NA。这种数值解决方案确保以最稳定的方式获得最小二乘解。

---

为了更好地理解这个问题，您可以使用 lm 来拟合一个带有 fix.formula 的线性模型。

fix.fit <- lm(fix.formula, data.df, method = "qr", singular.ok = TRUE)

method = "qr" 和 singular.ok = TRUE 是默认设置，因此我们实际上不需要进行设置。但如果我们指定 singular.ok = FALSE，lm 将停止并抱怨秩缺乏。

lm(fix.formula, data.df, method = "qr", singular.ok = FALSE)
#Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) : 
#  singular fit encountered

您可以在fix.fit中检查返回的值。

p <- length(coef)
coef <- fix.fit$coef
no.NA <- sum(is.na(coef))
rank <- fix.fit$rank

保证p = ncol(X)，但您应该查看no.NA = 7和rank + no.NA = p。

lmer内部发生的事情完全相同。 lm不会报告缺陷，而lmer会。实际上，这是有益的，因为我经常看到人们问为什么lm返回某些系数的NA。

---

更新1（2016-05-07）：

让我看看我是否理解正确：简而言之，我的一个预测变量与另一个变量相关，但我不用担心。使用因子是合适的吗？我仍然可以通过使用anova或查看BIC来比较模型吗？

不用担心使用summary或anova。方法被编写成使用正确数量的参数（自由度）以生成有效的摘要统计信息。

更新2（2016-11-06）：

让我们听听lme4软件包作者会说什么：rank deficiency warning mixed model lmer。Ben Bolker也提到了caret::findLinearCombos，特别是因为那里的OP想自己解决缺陷问题。

更新3（2018-07-27）：

排名缺陷对于有效模型估计和比较并不是问题，但在预测中可能会有危险。我最近在CrossValidated上用模拟例子写了一个详细的答案: R lm, Could anyone give me an example of the misleading case on “prediction from a rank-deficient”? 因此，在理论上，我们应该避免使用排名缺陷估计。但是，在现实中，没有所谓的“真实模型”：我们试图从数据中学习它。我们永远无法将估计的模型与“真相”进行比较；最好的选择是从我们建立的多个模型中选择最佳模型。因此，如果“最佳”模型最终出现排名缺陷，我们可以对其持怀疑态度，但可能立即无法做任何事情。

- Zheyuan Li

Answer 2

这个回答很好地解释了什么是秩亏，以及可能的原因。

具体如下：

1. 数据量过少：小于n个数据点无法精确估计n个参数。
2. 有太多重复的点。
3. 信息放错位置。
4. 模型过于复杂（变量过多）。
5. 单位和比例尺度。
6. 数字间的差异：12.001与12.005、44566与44555等。
7. 数据精度：即使双精度变量也有限制。

- Mox