r-是否有一种方法可以将列表中的每个组件/元素除以一个对应于某个因子变量级别的值



metafor中运行rma.mv命令时,我当前遇到模型收敛问题。我已经设法弄清楚这是因为我的R脚本中早期的错误代码,这导致计算的方差非常大(而且是错误的(。因此,我已经到了知道问题的地步,但我不知道如何解决它。下面我将概述我正在努力做什么,我目前是如何做到这一点的,以及我如何隔离问题的原因:

我想做什么

我试图做的要点是为我的数据集中的每个预测模型取自举的SEE值,并将其缩放到相应研究的标准均值。(即实际测量值(。值得注意的是,一些研究有多个模型。缩放操作如下:

(SEE/标准平均值(*100。

数据

structure(list(Study = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L), .Label = c("1", "2", "3", "4", "5", "6", "7"), class = "factor"), 
Model = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 
5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 
6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 
7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 
8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L), .Label = c("1", "2", "3", 
"4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", 
"15", "16", "17", "18", "19", "20", "21"), class = "factor"), 
Residual = c(1.585838423, 10.53555743, 10.67781267, 1.22720582, 
2.598848981, 26.94888607, 1.440304678, 15.74086952, 24.27068334, 
11.74554852, 13.21437544, 12.81558569, 15.65125984, 0.25625278, 
7.941652407, 5.557930412, 5.359401505, 8.537299607, 2.671115348, 
17.34002169, 2.598006397, 1.011493647, 10.38934426, 23.80005223, 
5.655188679, 4.878735797, 11.28737632, 1.048913043, 6.670565717, 
6.100551281, 15.55386342, 0.452170844, 1.077426851, 4.152827648, 
13.16037736, 4.308822184, 3.383948815, 16.53474723, 5.323566515, 
2.386768718, 5.63280155, 3.577780725, 1.60254086, 15.74086952, 
0.490982433, 20.16709778, 4.622970061, 1.894674528, 14.60716285, 
1.353952437, 4.596126567, 6.579200694, 2.749886696, 7.002473325, 
6.999046711, 5.986306941, 15.80934315, 0.028470501, 8.299180328, 
2.372475627, 9.538286164, 4.878735797, 6.671197129, 7.172826087, 
0.19253775, 3.098225566, 14.50304585, 0.680903636, 15.64487166, 
2.742369838, 11.42322707, 2.365470852, 1.838235294, 9.9, 
1.470588235, 1.041666667, 0.454545455, 22.45614035, 6.196296296, 
9.181818182, 31.97674419, 6.607142857, 7.720588235, 15.88888889, 
14.09090909, 2.5, 1.875, 0, 7.75, 4.6875, 6.25, 14.8, 10.88235294, 
0, 2.121212121, 0, 9.824074074, 9.181818182, 9.186046512, 
13.03571429, 0, 5.333333333, 13, 1.25, 14.19642857, 1.346153846, 
5.583333333, 2.8125, -2.097428958, 1.745602165, 3.531799729, 
4.262516915, 1.826792963, 2.963464141, -2.963464141, -2.584573748, 
-4.50608931, -4.912043302, -5.209742896, -6.211096076, -3.775372124, 
-4.803788904, -8.944519621, -1.146245059, 2.938076416, 4.914361001, 
5.12516469, 4.65085639, 1.567852437, -2.437417655, -3.544137022, 
-4.677206851, -5.441370224, -6.442687747, -1.963109354, -3.675889328, 
-7.891963109, -4.466403162), Criterion.Mean = c(162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 
162, 162, 162, 162, 162, 162, 162, 162, 162, 162, 79.4, 79.4, 
79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 
79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 
79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4, 79.4)), row.names = c(NA, 
138L), class = "data.frame")

我所做的

首先,我创建了一个用户定义的函数,用于计算SEESEE<- function(x){sqrt((sum(x)/(length(x)-2))^2)}

然后,我为每个模型生成了1000个样本,并为每个样本计算SEE。我根据需要使用了tapply函数,以确保SEE是基于与特定模型相对应的残差值计算的。代码如下:

Bootstrapped_SEE<- tapply(new.meta$Residual, INDEX = new.meta$Model, FUN = function (x){
int<- lapply(1:1000, function(i) sample(x, replace = T))
Calc.SEE<- sapply(int, SEE)
})

这返回一个长度为21的列表(整个数据集中有21个模型(,其中每个元素的长度为1000(1000个SEE值(。为了将这些转换为SEE值,我使用了以下代码:

SEE.percentage.List<-lapply(SEE_Var_List, FUN = function(i) (i/new.meta$Crtierion.Mean)*100)

问题

正如我上面所说,我希望这样做的是,对于列表中的每个元素(即与特定模型相对应的SEE向量(,除以相应研究的标准均值(同样,每个研究可能有多个模型,因此也有列表元素,所有这些模型都具有相同的标准均值(。刚才发生的是,R只是沿着列向量Criterion.Mean循环,并除以第i个值。好消息?R按预期工作!坏消息?在第19个SEE值之后,它开始将对应于研究1的SEE值除以研究2的Criterion.Mean值。因此,我需要的是一种将SEE值除以相应研究的标准均值的方法,其中对于给定的研究可以有多个模型(即列表的多个元素具有相同的标准均值(。

解决方案

概念解决方案是以某种方式将列表中的每个元素与特定研究联系起来,然后使用代码将与研究相关的所有元素除以相应的Criterion.Mean值。我最初认为某种tapply函数可以在这里工作,但后来很快意识到列表与包含的数据帧无关,因此提供study作为索引可能不起作用。

有人知道怎么做吗?我很失落,因为这超出了我目前的编码能力。我理解这个问题以及为什么会发生,但我不知道如何从语法上解决这个问题。提前感谢您的帮助。

对于未来关注这一点的人来说,以下是我如何使用base R绕过它的。我知道可能有更优雅、更流畅的解决方案,但这可以做到:

步骤的口头描述

  1. 首先,我创建了一个新的"子集";仅包含模型及其相应标准的数据帧平均值一次(即,我删除了重复的观测值(

  2. 接下来,我按模型订购了新的数据帧,纯粹是为了便于解释和检查。因为我在一个";循环;按照打字的方式,我先手动完成了。这最终是一场胜利,因为它允许我将循环的结果与手动计算进行交叉引用,以确保它们是正确的。。。我只需要确保列表元素与正确的模型相对应(即,元素1=模型1(。

  3. 然后,我使用lapply将var函数应用于包含缩放SEE值的新列表的每个元素。

  4. 瓦拉,你现在有结果了。漂亮吗?可能不会。它能更优雅吗?当然,也许可以。它有效吗?对

代码

这是我使用的代码:

#Create new dataframe, and remove duplicates.
Temp_Dat<- full.data[!duplicated(full.data$Model), ]
#Order new dataframe by model 
Temp_Dat<- Temp_Dat[order(Temp_Dat$Model), ]
#Create list to save results of for loop.

Bootstrapped_SEE_Per <- list()

#Run for loop. [[i]] indicates sub-element of list and [i] indicates element of vector. Therefore this is applied to each sub-element of the list using the"ith" value in the criterion mean vector. This creates a new list, with the only change being that now all SEE values are scaled by criterion mean and expressed as a %
for( i in 1:length(Bootstrapped_SEE)){
Bootstrapped_SEE_Per[[i]] <- (Bootstrapped_SEE[[i]]/ Temp_Dat$Criterion.Mean[i])*100
}
#Use lapply function to calculate the variance of each element of the new list.

SEE_Per_Variance<- lapply(Bootstrapped_SEE_Per, var)

希望这能在未来的某个时候帮助其他处境相似的人。对于那些阅读本文并使用tidyverse的人,我相信您可能能够使用purr包中的map2函数。然而,作为一个整体,我对tidyverse非常陌生,所以我选择了基本的R解决方案。

感谢所有的贡献!

相关内容

最新更新