0

我有一个关于 R 的问题。

我有一些按顺序编号的矩阵(所有相同的维度),我想搜索它们并生成一个最终矩阵,其中包含(对于每个矩阵元素)超过定义阈值的次数。

例如,我可以选择 0.7 的阈值,我可以有以下三个矩阵。

matrix1
    [,1] [,2] [,3]
[1,] 0.38 0.72 0.15
[2,] 0.58 0.37 0.09
[3,] 0.27 0.55 0.22

matrix2
    [,1] [,2] [,3]
[1,] 0.19 0.78 0.72
[2,] 0.98 0.65 0.46
[3,] 0.72 0.57 0.76

matrix3
     [,1] [,2] [,3]
[1,] 0.39 0.68 0.31
[2,] 0.40 0.05 0.92
[3,] 1.00 0.43 0.21

我想要的输出将是

      [,1] [,2] [,3]
[1,]    0    2    1
[2,]    1    0    1
[3,]    2    0    1

如果我这样做:

test <- matrix1 >= 0.7
test[test==TRUE] = 1

然后我得到一个矩阵,其中 1 表示超过阈值,0 表示未超过阈值。所以这是我想做的关键一步:

test=
      [,1] [,2] [,3]
[1,]    0    1    0
[2,]    0    0    0
[3,]    0    0    0

我的想法是做一个循环,所以我对每个矩阵执行这个计算并添加“测试”的每个结果,这样我就得到了我想要的最终矩阵。但我不确定两件事:如何在变量名“矩阵”中使用计数器,其次是否有比使用循环更有效的方法。

所以我在想这样的事情:

output = matrix(0,3,3)

for i in 1:3 {

test <- matrixi >= 0.7        
test[test==TRUE] = 1
output = output + test }

当然,这不起作用,因为 matrixi 不会转换为 matrix1、matrix2 等。

我真的很感谢你的帮助!!!

4

2 回答 2

2

如果您将矩阵存储在列表中,您会发现操作更容易:

lst <- list(matrix(c(0.38, 0.58, 0.27, 0.72, 0.37, 0.55, 0.15, 0.09, 0.22), nrow=3),
            matrix(c(0.19, 0.98, 0.72, 0.78, 0.65, 0.57, 0.72, 0.46, 0.76), nrow=3),
            matrix(c(0.39, 0.40, 1.00, 0.68, 0.05, 0.43, 0.31, 0.92, 0.21), nrow=3))
Reduce("+", lapply(lst, ">=", 0.7))
#      [,1] [,2] [,3]
# [1,]    0    2    1
# [2,]    1    0    1
# [3,]    2    0    1

在这里,lapply(lst, ">=", 0.7)返回一个列表,调用存储在 中x >= 0.7的每个矩阵。然后调用 with将它们全部加起来。xlstReduce+

如果你只有三个矩阵,你可以做类似的事情lst <- list(matrix1, matrix2, matrix3)。但是,如果您有更多(比如说 100,编号为 1 到 100),则可能更容易执行lst <- lapply(1:100, function(x) get(paste0("matrix", x)))or lst <- mget(paste0("matrix", 1:100))

对于 100 个矩阵,每个矩阵的大小为 100 x 100(根据您的评论,这大致是您的用例的大小),使用Reduce列表的方法似乎比rowSums使用数组的方法快一点,尽管两者都很快:

# Setup test data
set.seed(144)
for (i in seq(100)) {
    assign(paste0("matrix", i), matrix(rnorm(10000), nrow=100))
}

all.equal(sum.josilber(), sum.gavin())
# [1] TRUE
library(microbenchmark)
microbenchmark(sum.josilber(), sum.gavin())
# Unit: milliseconds
#            expr       min       lq   median       uq      max neval
#  sum.josilber()  6.534432 11.11292 12.47216 17.13995 160.1497   100
#     sum.gavin() 11.421577 16.54199 18.62949 23.09079 165.6413   100
于 2014-06-12T23:17:07.873 回答
0

如果将矩阵放入数组中,则无需循环即可轻松完成。这是一个例子:

## dummy data
set.seed(1)
m1 <- matrix(runif(9), ncol = 3)
m2 <- matrix(runif(9), ncol = 3)
m3 <- matrix(runif(9), ncol = 3)

将这些粘贴到数组中

arr <- array(c(m1, m2, m3), dim = c(3,3,3))

现在每个矩阵就像一个盘子,阵列是这些盘子的堆叠。

照你做的把数组转换成指标数组(你不需要保存这一步,它可以在下一次调用中内联完成)

ind <- arr > 0.7

这给出了:

> ind
, , 1

      [,1]  [,2]  [,3]
[1,] FALSE  TRUE  TRUE
[2,] FALSE FALSE FALSE
[3,] FALSE  TRUE FALSE

, , 2

      [,1]  [,2]  [,3]
[1,] FALSE FALSE FALSE
[2,] FALSE FALSE  TRUE
[3,] FALSE  TRUE  TRUE

, , 3

      [,1]  [,2]  [,3]
[1,] FALSE FALSE FALSE
[2,]  TRUE FALSE FALSE
[3,]  TRUE FALSE FALSE

现在使用该rowSums()函数来计算您想要的值

> rowSums(ind, dims = 2)
     [,1] [,2] [,3]
[1,]    0    1    1
[2,]    1    0    1
[3,]    1    2    1

请注意,汇总的rowSums()内容是(有点令人困惑!)维度dims + 1。在这种情况下,我们将每个 3*3 单元格的板堆栈(数组)中的值相加,因此输出中有 9 个值。

如果您需要将对象放入数组形式,您可以通过

arr2 <- do.call("cbind", mget(c("m1","m2","m3")))
dim(arr2) <- c(3,3,3) # c(nrow(m1), ncol(m1), nmat)

> all.equal(arr, arr2)
[1] TRUE

对于更大的问题(更多矩阵)使用类似的东西

nmat <- 200 ## number matrices
matrices <- paste0("m", seq_len(nmat))
arr <- do.call("cbind", mget(matrices))
dim(arr) <- c(dim(m1), nmat)
于 2014-06-12T23:17:14.513 回答