在 R 中使用 LongituRF 包实现纵向随机森林

Question

我有一些高维重复测量数据，我有兴趣拟合随机森林模型来研究此类模型的适用性和预测效用。具体来说，我正在尝试实现 LongituRF 包中的方法。此处详细介绍了此包背后的方法：

Capitaine, L., et al. Random forests for high-dimensional longitudinal data. Stat Methods Med Res (2020) doi:10.1177/0962280220946080.

作者方便地提供了一些有用的数据生成函数用于测试。所以我们有

install.packages("LongituRF")
library(LongituRF)

让我们使用 DataLongGenerator() 生成一些数据，这些数据将 n=样本大小、p=预测变量的数量和 G=具有时间行为的预测变量的数量作为参数。

my_data <- DataLongGenerator(n=50,p=6,G=6)

my_data 是您期望 Y（响应向量）的列表， X（固定效应预测矩阵），Z（随机效应预测矩阵）， id（样本标识符的向量）和时间（时间测量的向量）。简单地拟合随机森林模型

model <- REEMforest(X=my_data$X,Y=my_data$Y,Z=my_data$Z,time=my_data$time,
                    id=my_data$id,sto="BM",mtry=2)

这里大约需要 50 秒，请耐心等待

到目前为止一切顺利。现在我清楚这里除了 Z 之外的所有参数。 什么是 Z 什么时候我要在我的实际数据上拟合这个模型？

正在查看 my_data$Z。

dim(my_data$Z)
[1] 471   2
head(my_data$Z)
      [,1]      [,2]
 [1,]    1 1.1128914
 [2,]    1 1.0349287
 [3,]    1 0.7308948
 [4,]    1 1.0976203
 [5,]    1 1.3739856
 [6,]    1 0.6840415

每行看起来像一个截距项（即 1）和从均匀分布中提取的值 runif()。

REEMforest() 的文档指出“Z [矩阵]：包含随机效应的 q 个预测变量的 Nxq 矩阵。” 使用实际数据时，这个矩阵如何指定？

我的理解是，传统上 Z 只是组变量的单热（二进制）编码（例如 as described here），因此 DataLongGenerator() 中的 Z 应该是 nxG（ 471x6) 稀疏矩阵没有？

如能说明如何使用实际数据指定 Z 参数，我们将不胜感激。

编辑

我的具体例子如下，我有一个响应变量（Y）。样本（用 id 标识）被随机分配到干预组（I、干预或无干预）。一组高维特征 (X)。在两个时间点（Time、基线和终点）测量特征和反应。我对使用 X 和 I 预测 Y 很感兴趣。我也有兴趣提取哪些特征对预测最重要Y（与 Capitaine 等人在他们的论文中对 HIV 所做的相同）。

我会调用REEMforest()如下

REEMforest(X=cbind(X,I), Y=Y, time=Time, id=id)

我应该为 Z 使用什么？

Answer 1

当函数DataLongGenerator()创建Z时，它是一个矩阵中的运行dom统一数据。实际编码为

Z <- as.matrix(cbind(rep(1, length(f)), 2 * runif(length(f))))

其中 f 表示表示每个元素的矩阵的长度。在您的示例中，您使用了 6 组，每组 50 名参与者，具有 6 个固定效果。这导致长度为 472.

据我所知，由于此函数旨在模拟纵向数据，因此这是对运行dom 对该数据的影响的模拟。如果您使用的是真实数据，我认为它会更容易理解。

虽然这个例子没有使用 RE-EM 森林，但我认为它很清楚，因为它使用有形元素作为例子。您可以在第 1.2.2 节“固定与随机效应”中阅读有关运行dom 效应的信息。 https://ademos.people.uic.edu/Chapter17.html#32_fixed_effects

查看第 3.2 节以查看运行dom 效果的示例，如果您使用的是真实数据，则可以有意对其进行建模。

另一个示例：您正在运行进行抗癌药物试验。您每周收集患者人口统计数据：体重、体温、CBC 面板和不同的给药组：每天 1 个单位、每天 2 个单位和每天 3 个单位。

在传统回归中，您会对这些变量建模以确定模型识别结果的准确程度。固定效应是 可解释方差 或 R²。所以如果你有 0.86 或 86%，那么 14% 是无法解释的。它可能是导致噪音的相互作用，完美与模型确定的结果之间无法解释的差异。

假设白细胞计数非常低和超重的患者对治疗的反应要好得多。或者也许红头发的患者反应更好；那不在你的数据中。就纵向数据而言，假设这种关系（交互关系）仅在经过一段时间后才出现。

您可以尝试对不同的关系建模，以评估数据中的运行dom 交互。不过，我认为与运行dom 尝试识别运行dom 效应相比，使用系统评估交互的众多方法之一会更好。

EDITED 我开始在@JustGettinStarted 的评论中写这个，但是太多了。

没有背景 - 最简单的方法是运行类似REEMtree::REEMtree()的方法，设置运行 random = ~1 | time / id) 的 dom 效果参数。在它运行s 之后，提取它计算的运行dom 效果。你可以这样做：

data2 <- data %>% mutate(oOrder = row_number()) %>%  # identify original order of the data
     arrange(time, id) %>%
     mutate(zOrder = row_number()) # because the random effects will be in order by time then id

extRE <- data.frame(time = attributes(fit$RandomEffects[2][["id"]])[["row.names"]]) %>% 
     separate(col = time,
              into = c("time", "id"),
              sep = "\/") %>%
     mutate(Z = fit$RandomEffects[[2]] %>% unlist(),
            id = as.integer(id),
            time = time))        # set data type to match dataset for time

data2 <- data2 %>% left_join(extRE) %>% arrange(oOrder) # return to original order

Z = cbind(rep(1, times = nrows(data2)), data2$Z)

或者，我建议您从运行dom 生成运行dom 效果开始。您开始使用的运行dom-effects 只是一个起点。运行最后的dom效果会不一样

无论我尝试使用多少种方法来处理 LongituRF::REEMforest() 真实数据，我运行都会出错。我每次都遇到不可逆的矩阵故障。

我注意到 DataLongGenerator() 生成的数据按 id，然后是时间顺序排列。我试图以这种方式对数据（和 Z）进行排序，但没有帮助。当我从包 LongituRF 中提取所有功能时，我使用 MERF（多重效果运行dom 森林）功能没有任何问题。即使在研究论文中，这种方法也是可靠的。只是觉得值得一提。

在 R 中使用 LongituRF 包实现纵向随机森林

Implementing Longitudinal Random Forest with LongituRF package in R

regression

r

random-forest

mixed-models

longitudinal