与间隔匹配并提取两个矩阵 R 之间的值
match with an interval and extract values between two matrix R
我在列表中有 n 个矩阵和一个附加矩阵，其中包含我想在矩阵列表中找到的值。
要获取矩阵列表，我使用以下代码：
setwd("C:\~\Documents\R")
import.multiple.txt.files<-function(pattern=".txt",header=T)
{
list.1<-list.files(pattern=".txt")
list.2<-list()
for (i in 1:length(list.1))
{
list.2[[i]]<-read.delim(list.1[i])
}
names(list.2)<-list.1
list.2
}
txt.import.matrix<-cbind(txt.import)
我的列表是这样的：（我只展示了一个 n=2 的例子）。每个数组中的行数不同（这里我只取5行和6行来简化，但我的真实数据中有500多行）。
txt.import.matrix[1]
[[1]]
X. RT. Area. m.z.
1 1 1.01 2820.1 358.9777
2 2 1.03 9571.8 368.4238
3 3 2.03 6674.0 284.3294
4 4 2.03 5856.3 922.0094
5 5 3.03 27814.6 261.1299
txt.import.matrix[2]
[[2]]
X. RT. Area. m.z.
1 1 1.01 7820.1 358.9777
2 2 1.06 8271.8 368.4238
3 3 2.03 12674.0 284.3294
4 4 2.03 5856.6 922.0096
5 5 2.03 17814.6 261.1299
6 6 3.65 5546.5 528.6475
我想在矩阵列表中找到另一个值数组。该数组是通过将列表中的所有数组合并到一个数组中并删除重复项获得的。
reduced.list.pre.filtering
RT. m.z.
1 1.01 358.9777
2 1.07 368.4238
3 2.05 284.3295
4 2.03 922.0092
5 3.03 261.1299
6 3.56 869.4558
我想获得一个新矩阵，其中写入了列表中所有矩阵的匹配 RT. ± 0.02 和 m.z. ± 0.0002 的 Area. 结果。输出可能是这样的。
RT. m.z. Area.[1] Area.[2]
1 1.01 358.9777 2820.1 7820.1
2 1.07 368.4238 8271.8
3 2.05 284.3295 6674.0 12674.0
4 2.03 922.0092 5856.3
5 3.03 261.1299 27814.6
6 3.65 528.6475
我只知道如何在一个数组中只匹配一个精确值。这里的难点是在数组列表中找到值，需要找到值±一个区间。如果您有任何建议，我将不胜感激。
这是一个快速粗略的方法，如果我明白你想要做什么，可能会有所帮助。
取消列出两个矩阵的每个变量的值
areas <- unlist(lapply(txt.import.matrix, function(x) x$Area.))
rts <- unlist(lapply(txt.import.matrix, function(x) x$RT.))
mzs <- unlist(lapply(txt.import.matrix, function(x) x$m.z.))
找到 RT 和 m.z 的那些值的索引。最接近第三个值的 matrix/df:
rtmins <- lapply(reduced.list.pre.filtering$RT., function(x) which(abs(rts-x)==min(abs(rts-x))))
mzmins <- lapply(reduced.list.pre.filtering$m.z., function(x) which(abs(mzs-x)==min(abs(mzs-x))))
使用purrr快速计算出两者中有哪些索引（即每个索引的最小差异）：
inboth <- purrr::map2(rtmins,mzmins,intersect)
获取对应的面积值：
vals<-lapply(inboth, function(x) areas[x])
使用reshape2输入宽幅:
vals2 <- reshape2::melt(vals)
vals2$number <- ave(vals2$L1, vals2$L1, FUN = seq_along)
vals.wide <-reshape2::dcast(vals2, L1 ~ number, value.var="value")
cbind(reduced.list.pre.filtering, vals.wide)
# RT. m.z. L1 1 2
#1 1.01 358.9777 1 2820.1 7820.1
#2 1.07 368.4238 2 8271.8 NA
#3 2.05 284.3295 3 6674.0 12674.0
#4 2.03 922.0092 4 5856.3 NA
#5 3.03 261.1299 5 27814.6 NA
这可能会给您一些想法。可以轻松调整以检查共享最小值是否超过 +/- 值。
使用 data.table 当前开发版本 v1.9.7 中的 non-equi 连接（请参阅 installation instructions），它允许为连接提供非相等条件：
require(data.table) # v1.9.7
names(ll) = c("Area1", "Area2")
A = rbindlist(lapply(ll, as.data.table), idcol = "id") ## (1)
B = as.data.table(mat)
B[, c("RT.minus", "RT.plus") := .(RT.-0.02, RT.+0.02)]
B[, c("m.z.minus", "m.z.plus") := .(m.z.-0.0002, m.z.+0.0002)] ## (2)
ans = A[B, .(id, X., RT. = i.RT., m.z. = i.m.z., Area.),
on = .(RT. >= RT.minus, RT. <= RT.plus,
m.z. >= m.z.minus, m.z. <= m.z.plus)] ## (3)
dcast(ans, RT. + m.z. ~ id) ## (4)
# or dcast(ans, RT. + m.z. ~ id, fill = 0)
# RT. m.z. Area1 Area2
# 1: 1.01 358.9777 2820.1 7820.1
# 2: 1.07 368.4238 NA 8271.8
# 3: 2.03 922.0092 5856.3 NA
# 4: 2.05 284.3295 6674.0 12674.0
# 5: 3.03 261.1299 27814.6 NA
[1] 命名矩阵列表（此处称为 ll）并使用 lapply() 将每个矩阵转换为 data.table ，并使用 rbindlist 按行绑定它们，并将名称添加为额外的列 (idcol)。称之为 A.
[2] 将第二个矩阵（此处称为 mat）也转换为 data.table。添加与要搜索的 ranges/intervals 相对应的其他列（因为 on= 参数，正如我们将在下一步中看到的，无法处理表达式 yet).称之为 B.
[3] 执行条件 join/subset。对于 B 中的每一行，在 A 中找到与提供给 on= 参数的条件对应的匹配行，并提取这些匹配行索引的列 id, X., R.T. and m.z.。
[4]留在[3]比较好。但是，如果您希望它如您的答案所示，我们可以将其重塑为宽格式。 fill = 0 会将结果中的 NAs 替换为 0.
这是 Arun 使用 data.table 相当优雅的回答的另一种方法。我决定 post 它是因为它包含两个额外的方面，这两个方面在您的问题中是重要的考虑因素：
浮点数比较：比较浮点值是否在区间内需要考虑计算区间时的舍入误差。这是比较实数的浮点表示的一般问题。请参阅 R 上下文中的 this and this。以下在函数 in.interval.
多个匹配项： 如果间隔重叠，您的间隔匹配条件可能会导致多个匹配项。以下 假设 您只需要第一个匹配项（关于每个 txt.import.matrix 矩阵的增加行）。这是在函数 match.interval 中实现的，并在后面的注释中进行了解释。如果您想获得符合您标准的区域的平均值，则需要其他逻辑。
为了在矩阵 txt.import.matrix 中为矩阵 reduced.list.pre.filtering 中的每一行找到匹配行，以下代码将比较函数在 [=240= 上的应用向量化] reduced.list.pre.filtering 和来自 txt.import.matrix 的矩阵之间的所有枚举行对。从功能上讲，对于此应用程序，这与 Arun 使用 data.table 的 non-equi 连接的解决方案相同；然而，non-equi 连接功能更通用，而且 data.table 实现很可能针对内存使用和速度进行了更好的优化，即使对于此应用程序也是如此。
in.interval <- function(x, center, deviation, tol = .Machine$double.eps^0.5) {
return (abs(x-center) <= (deviation + tol))
}
match.interval <- function(r, t) {
r.rt <- rep(r[,1], each=nrow(t))
t.rt <- rep(t[,2], times=nrow(r))
r.mz <- rep(r[,2], each=nrow(t))
t.mz <- rep(t[,4], times=nrow(r)) ## 1.
ind <- which(in.interval(r.rt, t.rt, 0.02) &
in.interval(r.mz, t.mz, 0.0002))
r.ind <- floor((ind - 1)/nrow(t)) + 1 ## 2.
dup <- duplicated(r.ind)
r.ind <- r.ind[!dup]
t.ind <- ind[!dup] - (r.ind - 1)*nrow(t) ## 3.
return(cbind(r.ind,t.ind))
}
get.area.matched <- function(r, t) {
match.ind <- match.interval(r, t)
area <- rep(NA,nrow(r))
area[match.ind[,1]] <- t[match.ind[,2], 3] ## 4.
return(area)
}
res <- cbind(reduced.list.pre.filtering,
do.call(cbind,lapply(txt.import.matrix,
get.area.matched,
r=reduced.list.pre.filtering))) ## 5.
colnames(res) <- c(colnames(reduced.list.pre.filtering),
sapply(seq_len(length(txt.import.matrix)),
function(i) {return(paste0("Area.[",i,"]"))})) ## 6.
print(res)
## RT. m.z. Area.[1] Area.[2]
##[1,] 1.01 358.9777 2820.1 7820.1
##[2,] 1.07 368.4238 NA 8271.8
##[3,] 2.05 284.3295 6674.0 12674.0
##[4,] 2.03 922.0092 5856.3 NA
##[5,] 3.03 261.1299 27814.6 NA
##[6,] 3.56 869.4558 NA NA
备注：
这部分构造数据以实现比较函数在 reduced.list.pre.filtering 和矩阵 space 之间的所有枚举行对的应用的矢量化=19=]。要构造的数据是四个数组，它们是比较标准中使用的两列的复制（或扩展） reduced.list.pre.filtering 每个矩阵的行维度来自 txt.import.matrix 和复制reduced.list.pre.filtering 行维度中 txt.import.matrix 中每个矩阵的两列，用于比较标准。这里，术语数组指的是二维矩阵或一维向量。得到的四个数组是：
r.rt是reduced.list.pre.filtering的RT.列（即r[,1]）在t[=221=的行维度中的复制]
t.rt 是 rtxt.import.matrix 矩阵的 RT. 列的复制（即 t[,2]）rr.mz 是 reduced.list.pre.filtering 的 m.z. 列（即 r[,2]）在 tt.mz是r[=行维度中来自txt.import.matrix（即t[,4]）矩阵的m.z.列的复制221=]
重要的是，每个数组的索引都以相同的方式枚举了 r 和 t 中的所有行对。具体来说，将这些数组视为大小为 M by N 的二维矩阵，其中 M=nrow(t) 和 N=nrow(r)，行索引对应于 t 和列索引对应于 r 的行。因此，第 i 行和第 j 列（四个数组中的每一个）的数组值（在所有四个数组中）是 j - r 的第行和 i - t 的行。此复制过程的实现使用 R 函数 rep。例如，在计算 r.rt 时，使用 rep 和 each=M，其效果是将其数组输入 r[,1] 视为行向量并复制该行 M 次形成 M 行。结果是，对应于 r 中的一行的每一列都具有来自 r 对应行的 RT. 值，并且该值对于所有行（属于r.rt 的列），每个对应 t 中的一行。这意味着在将 r 中的该行与 t 中的任何行进行比较时，使用 r 中该行的 RT. 的值。相反，在计算 t.rt 时，使用 rep 和 times=N，其作用是将其数组输入视为列向量并将该列复制 N 次以形成一个N 列。结果是 t.rt 中的每一行对应于 t 中的一行，具有来自 t 对应行的 RT. 值，并且该值相同对于 t.rt 的所有（该行的）列，每一列对应于 r 中的一行。这意味着在将 t 中的该行与 r 中的任何行进行比较时，使用 t 中该行的 RT. 的值。同样，r.mz 和 t.mz 的计算分别使用 r 和 t 中的 m.z. 列。
这会执行向量化比较，生成 M by N 逻辑矩阵，其中第 i 行和第 j 行如果 r 的第 j 行与 t 的第 i 行匹配标准，则列为 TRUE，否则为 FALSE . which()的输出是这个逻辑比较结果矩阵的数组索引数组，其中它的元素是TRUE。我们想将这些数组索引转换为比较结果矩阵的行和列索引，以引用回 r 和 t 的行。下一行从数组索引中提取列索引。请注意，变量名称是 r.ind，表示这些对应于 r 的行。我们首先提取它，因为它对于检测 r.
这部分处理 r 中给定行的 t 中可能的多个匹配项。多个匹配项将在 r.ind 中显示为重复值。如上所述，这里的逻辑只保留 t 中增加行的第一个匹配项。函数duplicated returns 数组中所有重复值的索引。因此，删除这些元素将达到我们想要的效果。代码首先从r.ind中移除它们，然后从ind中移除它们，最后计算列索引到比较结果矩阵，它对应于t的行，使用修剪ind 和 r.ind。 match.interval 返回的是一个矩阵，其行是匹配的行索引对，其第一列是 r 的行索引，第二列是 t.[=164 的行索引=]
get.area.matched 函数只是使用 match.ind 的结果从 t 中提取所有匹配项的 Area。请注意，返回的结果是一个（列）向量，长度等于 r 中的行数并初始化为 NA。通过这种方式，r 中在 t 中没有匹配的行返回了 Area 的 NA。
这使用 lapply 在列表 txt.import.matrix 上应用函数 get.area.matched 并将返回的匹配 Area 结果附加到 reduced.list.pre.filtering 作为列向量。同样，相应的列名也被附加并设置在结果 res.
编辑： 使用 foreach 包的替代实现
事后看来，更好的实现方式是使用 foreach 包对比较进行矢量化。在此实现中，需要 foreach 和 magrittr 包
require("magrittr") ## for %>%
require("foreach")
然后match.interval中的代码用于向量化比较
r.rt <- rep(r[,1], each=nrow(t))
t.rt <- rep(t[,2], times=nrow(r))
r.mz <- rep(r[,2], each=nrow(t))
t.mz <- rep(t[,4], times=nrow(r)) # 1.
ind <- which(in.interval(r.rt, t.rt, 0.02) &
in.interval(r.mz, t.mz, 0.0002))
可以替换为
ind <- foreach(r.row = 1:nrow(r), .combine=cbind) %:%
foreach(t.row = 1:nrow(t)) %do%
match.criterion(r.row, t.row, r, t) %>%
as.logical(.) %>% which(.)
其中 match.criterion 定义为
match.criterion <- function(r.row, t.row, r, t) {
return(in.interval(r[r.row,1], t[t.row,2], 0.02) &
in.interval(r[r.row,2], t[t.row,4], 0.0002))
}
这更容易解析并反映正在执行的操作。请注意，嵌套的 foreach 与 cbind 组合返回的又是一个逻辑矩阵。最后，get.area.matched 函数在列表 txt.import.matrix 上的应用也可以使用 foreach:
执行
res <- foreach(i = 1:length(txt.import.matrix), .combine=cbind) %do%
get.area.matched(reduced.list.pre.filtering, txt.import.matrix[[i]]) %>%
cbind(reduced.list.pre.filtering,.)
使用foreach的完整代码如下：
require("magrittr")
require("foreach")
in.interval <- function(x, center, deviation, tol = .Machine$double.eps^0.5) {
return (abs(x-center) <= (deviation + tol))
}
match.criterion <- function(r.row, t.row, r, t) {
return(in.interval(r[r.row,1], t[t.row,2], 0.02) &
in.interval(r[r.row,2], t[t.row,4], 0.0002))
}
match.interval <- function(r, t) {
ind <- foreach(r.row = 1:nrow(r), .combine=cbind) %:%
foreach(t.row = 1:nrow(t)) %do%
match.criterion(r.row, t.row, r, t) %>%
as.logical(.) %>% which(.)
# which returns 1-D indices (row-major),
# convert these to 2-D indices in (row,col)
r.ind <- floor((ind - 1)/nrow(t)) + 1 ## 2.
# detect duplicates in r.ind and remove them from ind
dup <- duplicated(r.ind)
r.ind <- r.ind[!dup]
t.ind <- ind[!dup] - (r.ind - 1)*nrow(t) ## 3.
return(cbind(r.ind,t.ind))
}
get.area.matched <- function(r, t) {
match.ind <- match.interval(r, t)
area <- rep(NA,nrow(r))
area[match.ind[,1]] <- t[match.ind[,2], 3]
return(area)
}
res <- foreach(i = 1:length(txt.import.matrix), .combine=cbind) %do%
get.area.matched(reduced.list.pre.filtering, txt.import.matrix[[i]]) %>%
cbind(reduced.list.pre.filtering,.)
colnames(res) <- c(colnames(reduced.list.pre.filtering),
sapply(seq_len(length(txt.import.matrix)),
function(i) {return(paste0("Area.[",i,"]"))}))
希望对您有所帮助。
我在列表中有 n 个矩阵和一个附加矩阵，其中包含我想在矩阵列表中找到的值。
要获取矩阵列表，我使用以下代码：
setwd("C:\~\Documents\R")
import.multiple.txt.files<-function(pattern=".txt",header=T)
{
list.1<-list.files(pattern=".txt")
list.2<-list()
for (i in 1:length(list.1))
{
list.2[[i]]<-read.delim(list.1[i])
}
names(list.2)<-list.1
list.2
}
txt.import.matrix<-cbind(txt.import)
我的列表是这样的：（我只展示了一个 n=2 的例子）。每个数组中的行数不同（这里我只取5行和6行来简化，但我的真实数据中有500多行）。
txt.import.matrix[1]
[[1]]
X. RT. Area. m.z.
1 1 1.01 2820.1 358.9777
2 2 1.03 9571.8 368.4238
3 3 2.03 6674.0 284.3294
4 4 2.03 5856.3 922.0094
5 5 3.03 27814.6 261.1299
txt.import.matrix[2]
[[2]]
X. RT. Area. m.z.
1 1 1.01 7820.1 358.9777
2 2 1.06 8271.8 368.4238
3 3 2.03 12674.0 284.3294
4 4 2.03 5856.6 922.0096
5 5 2.03 17814.6 261.1299
6 6 3.65 5546.5 528.6475
我想在矩阵列表中找到另一个值数组。该数组是通过将列表中的所有数组合并到一个数组中并删除重复项获得的。
reduced.list.pre.filtering
RT. m.z.
1 1.01 358.9777
2 1.07 368.4238
3 2.05 284.3295
4 2.03 922.0092
5 3.03 261.1299
6 3.56 869.4558
我想获得一个新矩阵，其中写入了列表中所有矩阵的匹配 RT. ± 0.02 和 m.z. ± 0.0002 的 Area. 结果。输出可能是这样的。
RT. m.z. Area.[1] Area.[2]
1 1.01 358.9777 2820.1 7820.1
2 1.07 368.4238 8271.8
3 2.05 284.3295 6674.0 12674.0
4 2.03 922.0092 5856.3
5 3.03 261.1299 27814.6
6 3.65 528.6475
我只知道如何在一个数组中只匹配一个精确值。这里的难点是在数组列表中找到值，需要找到值±一个区间。如果您有任何建议，我将不胜感激。
这是一个快速粗略的方法，如果我明白你想要做什么，可能会有所帮助。
取消列出两个矩阵的每个变量的值
areas <- unlist(lapply(txt.import.matrix, function(x) x$Area.))
rts <- unlist(lapply(txt.import.matrix, function(x) x$RT.))
mzs <- unlist(lapply(txt.import.matrix, function(x) x$m.z.))
找到 RT 和 m.z 的那些值的索引。最接近第三个值的 matrix/df:
rtmins <- lapply(reduced.list.pre.filtering$RT., function(x) which(abs(rts-x)==min(abs(rts-x))))
mzmins <- lapply(reduced.list.pre.filtering$m.z., function(x) which(abs(mzs-x)==min(abs(mzs-x))))
使用purrr快速计算出两者中有哪些索引（即每个索引的最小差异）：
inboth <- purrr::map2(rtmins,mzmins,intersect)
获取对应的面积值：
vals<-lapply(inboth, function(x) areas[x])
使用reshape2输入宽幅:
vals2 <- reshape2::melt(vals)
vals2$number <- ave(vals2$L1, vals2$L1, FUN = seq_along)
vals.wide <-reshape2::dcast(vals2, L1 ~ number, value.var="value")
cbind(reduced.list.pre.filtering, vals.wide)
# RT. m.z. L1 1 2
#1 1.01 358.9777 1 2820.1 7820.1
#2 1.07 368.4238 2 8271.8 NA
#3 2.05 284.3295 3 6674.0 12674.0
#4 2.03 922.0092 4 5856.3 NA
#5 3.03 261.1299 5 27814.6 NA
这可能会给您一些想法。可以轻松调整以检查共享最小值是否超过 +/- 值。
使用 data.table 当前开发版本 v1.9.7 中的 non-equi 连接（请参阅 installation instructions），它允许为连接提供非相等条件：
require(data.table) # v1.9.7
names(ll) = c("Area1", "Area2")
A = rbindlist(lapply(ll, as.data.table), idcol = "id") ## (1)
B = as.data.table(mat)
B[, c("RT.minus", "RT.plus") := .(RT.-0.02, RT.+0.02)]
B[, c("m.z.minus", "m.z.plus") := .(m.z.-0.0002, m.z.+0.0002)] ## (2)
ans = A[B, .(id, X., RT. = i.RT., m.z. = i.m.z., Area.),
on = .(RT. >= RT.minus, RT. <= RT.plus,
m.z. >= m.z.minus, m.z. <= m.z.plus)] ## (3)
dcast(ans, RT. + m.z. ~ id) ## (4)
# or dcast(ans, RT. + m.z. ~ id, fill = 0)
# RT. m.z. Area1 Area2
# 1: 1.01 358.9777 2820.1 7820.1
# 2: 1.07 368.4238 NA 8271.8
# 3: 2.03 922.0092 5856.3 NA
# 4: 2.05 284.3295 6674.0 12674.0
# 5: 3.03 261.1299 27814.6 NA
[1] 命名矩阵列表（此处称为 ll）并使用 lapply() 将每个矩阵转换为 data.table ，并使用 rbindlist 按行绑定它们，并将名称添加为额外的列 (idcol)。称之为 A.
[2] 将第二个矩阵（此处称为 mat）也转换为 data.table。添加与要搜索的 ranges/intervals 相对应的其他列（因为 on= 参数，正如我们将在下一步中看到的，无法处理表达式 yet).称之为 B.
[3] 执行条件 join/subset。对于 B 中的每一行，在 A 中找到与提供给 on= 参数的条件对应的匹配行，并提取这些匹配行索引的列 id, X., R.T. and m.z.。
[4]留在[3]比较好。但是，如果您希望它如您的答案所示，我们可以将其重塑为宽格式。 fill = 0 会将结果中的 NAs 替换为 0.
这是 Arun 使用 data.table 相当优雅的回答的另一种方法。我决定 post 它是因为它包含两个额外的方面，这两个方面在您的问题中是重要的考虑因素：
浮点数比较：比较浮点值是否在区间内需要考虑计算区间时的舍入误差。这是比较实数的浮点表示的一般问题。请参阅 R 上下文中的 this and this。以下在函数
in.interval. 中实现此比较
多个匹配项： 如果间隔重叠，您的间隔匹配条件可能会导致多个匹配项。以下 假设 您只需要第一个匹配项（关于每个
txt.import.matrix矩阵的增加行）。这是在函数match.interval中实现的，并在后面的注释中进行了解释。如果您想获得符合您标准的区域的平均值，则需要其他逻辑。
为了在矩阵 txt.import.matrix 中为矩阵 reduced.list.pre.filtering 中的每一行找到匹配行，以下代码将比较函数在 [=240= 上的应用向量化] reduced.list.pre.filtering 和来自 txt.import.matrix 的矩阵之间的所有枚举行对。从功能上讲，对于此应用程序，这与 Arun 使用 data.table 的 non-equi 连接的解决方案相同；然而，non-equi 连接功能更通用，而且 data.table 实现很可能针对内存使用和速度进行了更好的优化，即使对于此应用程序也是如此。
in.interval <- function(x, center, deviation, tol = .Machine$double.eps^0.5) {
return (abs(x-center) <= (deviation + tol))
}
match.interval <- function(r, t) {
r.rt <- rep(r[,1], each=nrow(t))
t.rt <- rep(t[,2], times=nrow(r))
r.mz <- rep(r[,2], each=nrow(t))
t.mz <- rep(t[,4], times=nrow(r)) ## 1.
ind <- which(in.interval(r.rt, t.rt, 0.02) &
in.interval(r.mz, t.mz, 0.0002))
r.ind <- floor((ind - 1)/nrow(t)) + 1 ## 2.
dup <- duplicated(r.ind)
r.ind <- r.ind[!dup]
t.ind <- ind[!dup] - (r.ind - 1)*nrow(t) ## 3.
return(cbind(r.ind,t.ind))
}
get.area.matched <- function(r, t) {
match.ind <- match.interval(r, t)
area <- rep(NA,nrow(r))
area[match.ind[,1]] <- t[match.ind[,2], 3] ## 4.
return(area)
}
res <- cbind(reduced.list.pre.filtering,
do.call(cbind,lapply(txt.import.matrix,
get.area.matched,
r=reduced.list.pre.filtering))) ## 5.
colnames(res) <- c(colnames(reduced.list.pre.filtering),
sapply(seq_len(length(txt.import.matrix)),
function(i) {return(paste0("Area.[",i,"]"))})) ## 6.
print(res)
## RT. m.z. Area.[1] Area.[2]
##[1,] 1.01 358.9777 2820.1 7820.1
##[2,] 1.07 368.4238 NA 8271.8
##[3,] 2.05 284.3295 6674.0 12674.0
##[4,] 2.03 922.0092 5856.3 NA
##[5,] 3.03 261.1299 27814.6 NA
##[6,] 3.56 869.4558 NA NA
备注：
这部分构造数据以实现比较函数在
reduced.list.pre.filtering和矩阵 space 之间的所有枚举行对的应用的矢量化=19=]。要构造的数据是四个数组，它们是比较标准中使用的两列的复制（或扩展）reduced.list.pre.filtering每个矩阵的行维度来自txt.import.matrix和复制reduced.list.pre.filtering行维度中txt.import.matrix中每个矩阵的两列，用于比较标准。这里，术语数组指的是二维矩阵或一维向量。得到的四个数组是：r.rt是reduced.list.pre.filtering的RT.列（即r[,1]）在t[=221=的行维度中的复制]t.rt是rtxt.import.matrix矩阵的RT.列的复制（即t[,2]）rr.mz是reduced.list.pre.filtering的m.z.列（即r[,2]）在t 的行维度中的复制
t.mz是r[=行维度中来自txt.import.matrix（即t[,4]）矩阵的m.z.列的复制221=]
重要的是，每个数组的索引都以相同的方式枚举了
r和t中的所有行对。具体来说，将这些数组视为大小为MbyN的二维矩阵，其中M=nrow(t)和N=nrow(r)，行索引对应于t和列索引对应于r的行。因此，第i行和第j列（四个数组中的每一个）的数组值（在所有四个数组中）是j-r的第行和i-t的行。此复制过程的实现使用 R 函数rep。例如，在计算r.rt时，使用rep和each=M，其效果是将其数组输入r[,1]视为行向量并复制该行M次形成M行。结果是，对应于r中的一行的每一列都具有来自r对应行的RT.值，并且该值对于所有行（属于r.rt的列），每个对应t中的一行。这意味着在将r中的该行与t中的任何行进行比较时，使用r中该行的RT.的值。相反，在计算t.rt时，使用rep和times=N，其作用是将其数组输入视为列向量并将该列复制N次以形成一个N列。结果是t.rt中的每一行对应于t中的一行，具有来自t对应行的RT.值，并且该值相同对于t.rt的所有（该行的）列，每一列对应于r中的一行。这意味着在将t中的该行与r中的任何行进行比较时，使用t中该行的RT.的值。同样，r.mz和t.mz的计算分别使用r和t中的m.z.列。这会执行向量化比较，生成
MbyN逻辑矩阵，其中第i行和第j行如果r的第j行与t的第i行匹配标准，则列为TRUE，否则为FALSE.which()的输出是这个逻辑比较结果矩阵的数组索引数组，其中它的元素是TRUE。我们想将这些数组索引转换为比较结果矩阵的行和列索引，以引用回r和t的行。下一行从数组索引中提取列索引。请注意，变量名称是r.ind，表示这些对应于r的行。我们首先提取它，因为它对于检测r. 中一行的多个匹配很重要
这部分处理
r中给定行的t中可能的多个匹配项。多个匹配项将在r.ind中显示为重复值。如上所述，这里的逻辑只保留t中增加行的第一个匹配项。函数duplicatedreturns 数组中所有重复值的索引。因此，删除这些元素将达到我们想要的效果。代码首先从r.ind中移除它们，然后从ind中移除它们，最后计算列索引到比较结果矩阵，它对应于t的行，使用修剪ind和r.ind。match.interval返回的是一个矩阵，其行是匹配的行索引对，其第一列是r的行索引，第二列是t.[=164 的行索引=]get.area.matched函数只是使用match.ind的结果从t中提取所有匹配项的Area。请注意，返回的结果是一个（列）向量，长度等于r中的行数并初始化为NA。通过这种方式，r中在t中没有匹配的行返回了Area的NA。这使用
lapply在列表txt.import.matrix上应用函数get.area.matched并将返回的匹配Area结果附加到reduced.list.pre.filtering作为列向量。同样，相应的列名也被附加并设置在结果res.
编辑： 使用 foreach 包的替代实现
事后看来，更好的实现方式是使用 foreach 包对比较进行矢量化。在此实现中，需要 foreach 和 magrittr 包
require("magrittr") ## for %>%
require("foreach")
然后match.interval中的代码用于向量化比较
r.rt <- rep(r[,1], each=nrow(t))
t.rt <- rep(t[,2], times=nrow(r))
r.mz <- rep(r[,2], each=nrow(t))
t.mz <- rep(t[,4], times=nrow(r)) # 1.
ind <- which(in.interval(r.rt, t.rt, 0.02) &
in.interval(r.mz, t.mz, 0.0002))
可以替换为
ind <- foreach(r.row = 1:nrow(r), .combine=cbind) %:%
foreach(t.row = 1:nrow(t)) %do%
match.criterion(r.row, t.row, r, t) %>%
as.logical(.) %>% which(.)
其中 match.criterion 定义为
match.criterion <- function(r.row, t.row, r, t) {
return(in.interval(r[r.row,1], t[t.row,2], 0.02) &
in.interval(r[r.row,2], t[t.row,4], 0.0002))
}
这更容易解析并反映正在执行的操作。请注意，嵌套的 foreach 与 cbind 组合返回的又是一个逻辑矩阵。最后，get.area.matched 函数在列表 txt.import.matrix 上的应用也可以使用 foreach:
res <- foreach(i = 1:length(txt.import.matrix), .combine=cbind) %do%
get.area.matched(reduced.list.pre.filtering, txt.import.matrix[[i]]) %>%
cbind(reduced.list.pre.filtering,.)
使用foreach的完整代码如下：
require("magrittr")
require("foreach")
in.interval <- function(x, center, deviation, tol = .Machine$double.eps^0.5) {
return (abs(x-center) <= (deviation + tol))
}
match.criterion <- function(r.row, t.row, r, t) {
return(in.interval(r[r.row,1], t[t.row,2], 0.02) &
in.interval(r[r.row,2], t[t.row,4], 0.0002))
}
match.interval <- function(r, t) {
ind <- foreach(r.row = 1:nrow(r), .combine=cbind) %:%
foreach(t.row = 1:nrow(t)) %do%
match.criterion(r.row, t.row, r, t) %>%
as.logical(.) %>% which(.)
# which returns 1-D indices (row-major),
# convert these to 2-D indices in (row,col)
r.ind <- floor((ind - 1)/nrow(t)) + 1 ## 2.
# detect duplicates in r.ind and remove them from ind
dup <- duplicated(r.ind)
r.ind <- r.ind[!dup]
t.ind <- ind[!dup] - (r.ind - 1)*nrow(t) ## 3.
return(cbind(r.ind,t.ind))
}
get.area.matched <- function(r, t) {
match.ind <- match.interval(r, t)
area <- rep(NA,nrow(r))
area[match.ind[,1]] <- t[match.ind[,2], 3]
return(area)
}
res <- foreach(i = 1:length(txt.import.matrix), .combine=cbind) %do%
get.area.matched(reduced.list.pre.filtering, txt.import.matrix[[i]]) %>%
cbind(reduced.list.pre.filtering,.)
colnames(res) <- c(colnames(reduced.list.pre.filtering),
sapply(seq_len(length(txt.import.matrix)),
function(i) {return(paste0("Area.[",i,"]"))}))
希望对您有所帮助。