如何使用公式在 R 中可靠地使用矩阵（多变量）响应？

Question

我试图在R 中使用多变量响应，但遇到了该死的公式，并且出现了各种意外行为，主要是在函数和包中使用它们时。这个问题是双重的

我可以输入一个多变量响应，然后能够使用模型进行预测吗？

此 MVE 使用 rpart 包作为示例。这里y 是一个两列矩阵（响应），我想使用x 预测y，即x 中的每一列（此MVE 中的两列）。请注意method 本身和y 中每一列的含义无关，这只是重现问题的MVE：

library(rpart)

set.seed(1)
y <- rpois(10, lambda = 1.25)
y <- cbind(c(1,4,10,11,12, 14,16,17,20, 21), y)
print(y)
x <- matrix(1:20, ncol = 2) # just two dummy predictors
print(x)

mymodel <- rpart(y ~ x, method = "poisson", minbucket = 1)

newx <- matrix(11:20, ncol = 2) # just some dummy test predictors, note that we have less rows
predict(mymodel, newdata = data.frame(newx))
# output:
        1          2          3          4          5          6          7          8          
 9         10 
 0.12500000 0.12500000 0.12500000 0.20000000 0.04761905 0.17948718 0.17948718 0.11538462 
 0.04000000 0.04000000 
 Warning message:
 'newdata' had 5 rows but variables found have 10 rows 

如您所见，我无法预测新数据集。我一直在弄乱列名和行名，但无法使其正常工作。

另外，

如何制作“安全”的包装器？

例如，在这个 MVE 中：

mywrapper <- function(y, x){
  mymodel <- rpart(y ~ x, method = "poisson", minbucket = 1)
  
  return(mymodel)
}

并提供了 R 文档中提供的帮助：

公式对象具有关联的环境，model.frame 使用此环境（而不是父环境）来评估在提供的数据参数中找不到的变量。 使用 ~ 运算符创建的公式使用创建它们的环境。使用 as.formula 创建的公式将为其环境使用 env 参数。

我真的不明白这是什么意思。据我了解，不输入也不输入y 或x 到mywrapper() 将导致错误（这是预期的行为）。我问是因为我正在处理 r 包和包内的函数，并且我想确保公式没有意外行为。

Answer 1

一般来说，R 中的公式适用于数据框。 rpart 适用于矩阵，虽然数据框可以保存矩阵，但它们往往会转换为单独的列。为避免这种情况，请将矩阵包装在I()：

# Same as your code to start...then this:

predict(mymodel, newdata = data.frame(x = I(newx)))
#>    1    2    3    4    5 
#> 0.04 0.04 0.04 0.04 0.04

在问题的第二部分中，您将在 mywrapper 函数中创建一个公式，因此如果变量不包含在 newdata 数据框中，它将在那里查找变量。 R 中的“环境”类似于其他语言中的“栈帧”；主要区别在于环境只有一个父对象，如果在原始环境中找不到对象，则搜索会在那里进行。

一般来说，父级不是调用者的框架，它是环境被创建的框架，或者被特别列为父级的东西。

如果您对来自mywrapper 的返回值运行predict 会发生什么，它会查看公式以找到它需要的变量。预测只需要右侧的变量，所以这只是x。如果您在newdata 参数中为predict 提供x，一切都会好起来并像以前一样继续进行，但如果您不这样做，事情就会不同。

由于在newdata 中没有找到x，所以它进入公式的环境。那是mywrapper 的评估框架，它会在那里看到x，因为它是该函数的参数。

如果它正在寻找z，它不会在那里找到它。下一个要查看的地方是父环境，它是在创建mywrapper 时生效的环境，即全局环境。如果那里没有z，它将搜索search() 列出的环境链，这些环境通常是包导出。 search() 列表链接在一起，因此每个条目都是前一个条目的父条目。

我希望这不是太多的信息.​​.....

Answer 2

我没有使用过rpart:predict，但是根据文档，对于这个函数，您需要一个与原始数据集具有相同变量的新数据集。

因此，您应该使用正确的列名启动 newx：

newx = matrix(11:20, ncol = 2,dimnames=list(NULL,c("x1","x2")))

现在，这些列被标记为 x1 和 x2，就像模型和预测中的变量知道如何处理这些列一样。