R语言lm回归模型predict()返回行数不符预期的问题求助

我按如下步骤拟合了lm回归模型：

y <- hiedata_short$income
x1 <- hiedata_short$education
x2 <- hiedata_short$age
x3 <- hiedata_short$ghi

lm.more <- lm(formula = y ~ x1 + x2 + x3,
             data = hiedata_short
)

summary(lm.more)

随后希望固定age和ghi，对0到25的education值进行预测，创建了新数据框：

education.pts <- seq(from = 0, # 起始值
                     to = 25, # 结束值
                     by = 1 # 步长
)

fitted.values.data <- data.frame(age = 30, # 固定age为30
                                 education = education.pts,
                                 ghi = 70
)

接着使用模型进行预测并转换为数据框：

fitted.values <- predict(object = lm.more,
                         newdata = fitted.values.data, 
                         interval = "confidence", 
                         level = 0.95
)

fitted.values <- data.frame(fitted.values) # 将矩阵转为数据框

但绘图时失败，出现警告：

Warning message: 'newdata' had 26 rows but variables found have 1932 rows

fitted.values$fit的行数远多于education.pts的26行，predict函数返回的行数与模型原数据集一致，而非newdata的行数，这是为什么？

绘图代码：

plot(x = education.pts, 
     y = fitted.values$fit,
     type = "l",
     col = "black",
     xlab = "受教育年限",
     ylab = "预测家庭收入（单位：千美元）", 
     main = "受教育程度对应的预测收入\nRAND健康保险实验\n",
     sub = "年龄固定为30；GHI固定为70", 
     ylim = c(0, 20) 
)

原因分析 #

问题出在模型公式的变量名与newdata的列名不匹配：

• 拟合模型时，你用的是全局环境中定义的x1、x2、x3作为公式变量，而非原数据框hiedata_short中的列名（education、age、ghi）。
• 当调用predict()时，R会优先在全局环境中寻找与模型公式里同名的变量（也就是x1、x2、x3，对应原数据集的1932行数据），而非newdata中的education、age、ghi列。这就导致predict直接使用了原数据集的变量长度，忽略了newdata的26行数据。

解决方法 #

有两种可行的修正方式：

方法1：拟合模型时直接使用数据框列名写公式

修改模型拟合代码，直接用原数据框的列名构建公式，这样模型会关联数据框的列名，predict时会匹配newdata的列名：

# 无需单独定义y、x1、x2、x3
lm.more <- lm(formula = income ~ education + age + ghi,
             data = hiedata_short
)

summary(lm.more)

之后再运行原有的newdata创建和predict代码即可，此时predict会正确识别newdata中的列名，返回26行结果。

方法2：让newdata的列名与模型公式变量名一致

如果不想修改模型，就把newdata的列名改成x1、x2、x3，对应模型里的变量：

fitted.values.data <- data.frame(x2 = 30, # x2对应age
                                 x1 = education.pts, # x1对应education
                                 x3 = 70 # x3对应ghi
)

再运行predict代码，就能得到与newdata行数一致的结果。

内容的提问来源于stack exchange，提问作者yuskox