在循环中使用predict.lm()。警告：从秩缺失的拟合中进行预测可能会导致误导。

Question

这段 R 代码会抛出一个警告

# Fit regression model to each cluster
y <- list() 
length(y) <- k
vars <- list() 
length(vars) <- k
f <- list()
length(f) <- k

for (i in 1:k) {
  vars[[i]] <- names(corc[[i]][corc[[i]]!= "1"])
  f[[i]]  <- as.formula(paste("Death ~", paste(vars[[i]], collapse= "+")))
  y[[i]]  <- lm(f[[i]], data=C1[[i]]) #training set
  C1[[i]] <- cbind(C1[[i]], fitted(y[[i]]))
  C2[[i]] <- cbind(C2[[i]], predict(y[[i]], C2[[i]])) #test set
}

我有一个训练数据集（C1）和一个测试数据集（C2），每个数据集都有129个变量。我对C1进行了k均值聚类分析，然后基于簇成员身份拆分了我的数据集，并创建了不同簇的列表（C1 [ [1]]，C1 [ [2]]，...，C1 [ [k]]）。我还为C2中的每个案例分配了一个簇成员身份，并创建了C2 [ [1]]，...，C2 [ [k]]。然后我对C1中的每个簇进行线性回归拟合。我的因变量是“死亡”。在每个簇中，我的预测变量都不同，并且vars [ [i]]（i = 1，...，k）显示了预测变量的名称列表。我想预测测试数据集（C2 [ [1]]，...，C2 [ [k]]）中每个案例的死亡情况。当我运行以下代码时，对于其中一些簇。

我收到了这个警告：

In predict.lm(y[[i]], C2[[i]]) :
prediction from a rank-deficient fit may be misleading

我读了很多关于这个警告的内容，但是我还是没能弄清楚问题出在哪里。

- Mahsa

---

问题在于你得到了秩不足的拟合结果。你需要找出哪些拟合结果出现了警告并仔细检查它们。 - Roland

Answer 1

你可以使用 body(predict.lm) 来检查预测函数。你会看到这一行：

if (p < ncol(X) && !(missing(newdata) || is.null(newdata))) 
    warning("prediction from a rank-deficient fit may be misleading")

这个警告检查您的数据矩阵的秩是否至少等于您想要拟合的参数数量。触发它的一种方式是具有一些共线性协变量：

data <- data.frame(y=c(1,2,3,4), x1=c(1,1,2,3), x2=c(3,4,5,2), x3=c(4,2,6,0), x4=c(2,1,3,0))
data2 <- data.frame(x1=c(3,2,1,3), x2=c(3,2,1,4), x3=c(3,4,5,1), x4=c(0,0,2,3))
fit <- lm(y ~ ., data=data)

predict(fit, data2)
       1        2        3        4 
4.076087 2.826087 1.576087 4.065217 
Warning message:
In predict.lm(fit, data2) :
  prediction from a rank-deficient fit may be misleading

注意在data中x3和x4具有相同的方向。其中一个是另一个的倍数。可以使用length(fit$coefficients) > fit$rank进行检查。

另一种方法是参数数量多于可用变量：

fit2 <- lm(y ~ x1*x2*x3*x4, data=data)
predict(fit2, data2)
Warning message:
In predict.lm(fit2, data2) :
  prediction from a rank-deficient fit may be misleading

- Karolis Koncevičius

---

您IP地址为143.198.54.68，由于运营成本限制，当前对于免费用户的使用频率限制为每个IP每72小时10次对话，如需解除限制，请点击左下角设置图标按钮（手机用户先点击左上角菜单按钮）。 - Mahsa

在您的模型中，拥有比变量更多的案例非常重要。您可以尝试在案例不足的情况下进行操作，但是请记住，在这种情况下，您的预测可能不可靠。这就是为什么R会给出“警告”而不是错误的原因。只是为了引起您的注意。即使收到警告（它们不是错误），您也应该能够获得答案并继续工作，但最好尝试简化您的模型。 - Karolis Koncevičius

Answer 2

这是一条警告：

In predict.lm(model, test) :
  prediction from a rank-deficient fit may be misleading

从 R 的 predict.lm 抛出。请参阅：http://stat.ethz.ch/R-manual/R-devel/library/stats/html/predict.lm.html

理解秩缺失：请求 R 告诉您矩阵的秩：

train <- data.frame(y=c(1234, 325, 152, 403), 
                   x1=c(3538, 324, 382, 335), 
                   x2=c(2985, 323, 223, 288), 
                   x3=c(8750, 322, 123, 935))
test <- data.frame(x1=c(3538, 324, 382, 335), 
                   x2=c(2985, 323, 223, 288), 
                   x3=c(8750, 322, 123, 935))
library(Matrix)
cat(rankMatrix(train), "\n")   #prints 4
cat(rankMatrix(test), "\n")    #prints 3

如果矩阵没有“满秩”，则称其为“秩不足”。如果矩阵的秩等于其列数或行数（或两者都相等），则称该矩阵具有满秩。

问题在于，即使你的矩阵具有完全秩（而非秩不足），predict.lm 也会抛出警告，因为它在幕后进行了一些操作，丢弃了它认为无用的特征，修改了您的完全秩输入以成为秩不足。然后，它通过警告来抱怨。

此外，这个警告似乎是一个通用警告，适用于其他情况，例如你的输入特征太多，数据密度太低，它会提出自己的意见，认为预测结果是脆弱的。

传递完全秩矩阵的示例，但predict.lm仍抱怨秩不足

train <- data.frame(y=c(1,2,3,4),
                        x1=c(1,1,2,3),
                        x2=c(3,4,5,2),
                        x3=c(4,2,6,0),
                        x4=c(2,1,3,0)
                   )
test <- data.frame(x1=c(1, 2,  3,  9),
                   x2=c(3, 5,  1, 15),
                   x3=c(5, 9,  5, 22),
                   x4=c(9, 13, 2, 99))
library(Matrix)
cat(rankMatrix(train), "\n")    #prints 4, is full rank, good to go
cat(rankMatrix(test), "\n")     #prints 4, is full rank, good to go
myformula = as.formula("y ~ x1+x2+x3+x4")
model <- lm(myformula, train)
predict(model, test) 
    #Warning: prediction from a rank-deficient fit may be misleading

解决方法:

如果 predict 返回的预测结果良好，您可以忽略警告。predict.lm 在缺乏足够角度的情况下提供其意见，而您正在这里。

因此，请在预测步骤上禁用警告，如下所示：

options(warn=-1)      #turn off warnings
predict(model, test)
options(warn=1)      #turn warnings back on

- Eric Leschinski

---

我认为这个答案是错误的。你需要考虑设计矩阵的秩，而不是包含响应变量列的数据框的秩！第二个例子中的“train”数据（通过x4列的y）是秩缺失的。尝试X <- model.matrix(y ~ ., train); rankMatrix(X) < ncol(X)，你会发现设计矩阵是秩缺失的。我认为你所说的R会修改“你的完整秩输入为秩缺失”的说法是错误的。 - Cameron Bieganek

Answer 3

这是因为，在lm（..）函数输出的系数中，你的一个依赖变量具有NA值。由于多重共线性问题，这样的变量对模型没有任何影响，即该预测变量在其他预测变量上线性相关，或者该预测变量对所有记录（行）都是常数。最好的做法是从lm（..）函数的公式中删除该变量，并重新进行回归。这不会降低模型的准确性。 在我的情况下，

model <- lm(Happiness.Score ~ Economy..GDP.per.Capita.+year+Health..Life.Expectancy., data=dfTrain)

> model
Call:
lm(formula = Happiness.Score ~ Economy..GDP.per.Capita. + year + 
    Health..Life.Expectancy., data = dfTrain)

Coefficients:
             (Intercept)  Economy..GDP.per.Capita.                      year  
                   3.036                     1.569                        NA  
Health..Life.Expectancy.  
                   1.559

变量“year”对于所有记录具有相同的值。在删除“year”变量后

model <- lm(Happiness.Score ~ Economy..GDP.per.Capita.+Health..Life.Expectancy., data=dfTrain)

preds <- predict.lm(model, dfTest[, c(1:nrow(dfTest)-1]))

这不会给出任何警告消息

- Manjunath Bhadrannavar