logo标签列表

在R中查找向量中的元素

rmatrixvector
4
我有一个矩阵，恰好有2行n列，比如这样：
c(0,0,0,0,1,0,2,0,1,0,1,1,1,0,2)->a1
c(0,2,0,0,0,0,2,1,1,0,0,0,0,2,0)->a2
rbind(a1,a2)->matr

针对特定的列（在本例中为第9列，两行均为1），我需要找到左侧和右侧的第一个2/0或0/2实例 - 在本例中，左侧是2，右侧是14）

每一行的元素只能是0、1、2 - 没有其他。是否有一种方法可以快速地在大矩阵（有2行）上执行此操作？我需要执行600k次，因此速度可能是一个考虑因素。

你是否总是只对一个列感兴趣？还是你要寻找所有包含两个一的左右列？ - Heroka
@Heroka 不仅仅是特定的列。 - kutyw
这些数字总是在这样一个特定的低范围内吗？您是否只会搜索一个模式，而不考虑顺序，或者另一种搜索可能涉及“0/2和1/2”？ - alexis_laz
5个回答
2
library(compiler)
myfun <- cmpfun(function(m, cl) {
  li <- ri <- cl
  nc <- ncol(m)
  repeat {
    li <- li - 1
    if(li == 0 || ((m[1, li] != 1) && (m[1, li] + m[2, li] == 2))) {
      l <- li
      break
    }
  }
  repeat {
    ri <- ri + 1
    if(ri == nc || ((m[1, ri] != 1) && (m[1, ri] + m[2, ri] == 2))) {
      r <- ri
      break
    }
  }
  c(l, r)
})

在考虑了@Martin Morgan的意见之后，

set.seed(1)
N <- 1000000
test <- rbind(sample(0:2, N, replace = TRUE),
              sample(0:2, N, replace = TRUE))

library(microbenchmark)
microbenchmark(myfun(test, N / 2), fun(test, N / 2), foo(test, N / 2),
               AWebb(test, N / 2), RHertel(test, N / 2))
# Unit: microseconds
               expr         min          lq         mean      median          uq         max neval  cld
#    myfun(test, N/2)       4.658      20.033 2.237153e+01      22.536      26.022      85.567   100 a   
#      fun(test, N/2)   36685.750   47842.185 9.762663e+04   65571.546  120321.921  365958.316   100  b  
#      foo(test, N/2) 2622845.039 3009735.216 3.244457e+06 3185893.218 3369894.754 5170015.109   100    d
#    AWebb(test, N/2)  121504.084  142926.590 1.990204e+05  193864.670  209918.770  489765.471   100   c 
#  RHertel(test, N/2)   65998.733   76805.465 1.187384e+05   86089.980  144793.416  385880.056   100  b  

set.seed(123)
test <- rbind(sample(0:2, N, replace = TRUE, prob = c(5, 90, 5)),
              sample(0:2, N, replace = TRUE, prob = c(5, 90, 5)))
microbenchmark(myfun(test, N / 2), fun(test, N / 2), foo(test, N / 2),
               AWebb(test, N / 2), RHertel(test, N / 2))
# Unit: microseconds
#                expr         min          lq         mean      median         uq         max neval  cld
#    myfun(test, N/2)      81.805     103.732     121.9619     106.459     122.36     307.736   100 a   
#      fun(test, N/2)   26362.845   34553.968   83582.9801   42325.755  106303.84  403212.369   100  b  
#      foo(test, N/2) 2598806.742 2952221.561 3244907.3385 3188498.072 3505774.31 4382981.304   100    d
#    AWebb(test, N/2)  109446.866  125243.095  199204.1013  176207.024  242577.02  653299.857   100   c 
#  RHertel(test, N/2)   56045.309   67566.762  125066.9207   79042.886  143996.71  632227.710   100  b
2
对于第一个情况，将 all(...) 替换为 ((m[1, li] != 1) && (m[1, li] + m[2, li] == 2)) 可以提高大约2倍的速度，而 compiler::cmpfun() 又可以再次提高2倍的性能。对于第二种情况，加速比分别为3.5和17倍；编译使用 all() 的版本只有很小的影响。@A.Webb的答案似乎比这些答案中最快的版本慢了约500倍（我意识到它没有处理边缘情况），但当针对给定数据集进行100个以上的查询时，它开始领先。 - Martin Morgan
0

通过平方行并将它们相加来合并信息。正确的结果应该是4。然后，只需找到第一列小于9的列（rev(which())[1]）和第一列大于9的列（which()[1]）。

fun <- function(matr, col){
    valid <- which((matr[1,]^2 + matr[2,]^2) == 4)
    if (length(valid) == 0) return(c(NA,NA))

    left <- valid[rev(which(valid < col))[1]]
    right <- valid[which(valid > col)[1]]

    c(left,right)

    }

fun(matr,9)
# [1]  2 14

fun(matr,1)
# [1] NA  2

fun(matrix(0,nrow=2,ncol=100),9)
# [1] NA NA

基准测试

set.seed(1)
test <- rbind(sample(0:2,1000000,replace=T),
              sample(0:2,1000000,replace=T))

microbenchmark::microbenchmark(fun(test,9))
# Unit: milliseconds
#         expr     min       lq     mean   median       uq      max neval
# fun(test, 9) 22.7297 27.21038 30.91314 27.55106 28.08437 51.92393   100

编辑：感谢 @MatthewLundberg 指出了很多错误。

我没有给你点踩，但是楼主特别要求使用矩阵来解决问题。 - Heroka
@Heroka，说得对，我只是用a1和a2作为输入数据，但现在我已经改写成了rbind(vec1，vec2)。 - slamballais
@Laterow 有没有可能添加一个功能，如果左侧或右侧没有类似于2/0或0/2这样的东西，只显示NA或0以及dim（matr）[2]（当在右侧时）？ - kutyw
@Laterow，您能否在您的解决方案中添加一些内容，使其更加完整和有用，以便其他人也能受益 - 否则这段代码将被埋没在这里 - 再次感谢您的回复！ - kutyw
fun(matr,10) 应该返回什么？ - Matthew Lundberg
@MatthewLundberg 我误读了，只有最接近的列应该被打印出来... 再次感谢。我要去睡觉了，因为这太草率了。 - slamballais
0

我比 @Laterow 慢，但无论如何，这是类似的方法

foo  <- function(mtr, targetcol) {
  matr1  <-  colSums(mtr)
  matr2  <- apply(mtr, 2, function(x) x[1]*x[2])
  cols  <-  which(matr1 == 2 & matr2 == 0) - targetcol
  left  <-   cols[cols < 0]
  right  <-  cols[cols > 0]
  c(ifelse(length(left) == 0, NA, targetcol + max(left)),
    ifelse(length(right) == 0, NA, targetcol + min(right)))
}

foo(matr,9) #2 14
有没有可能添加一个功能，如果左侧或右侧没有2/0或0/2这样的内容，只显示NA或0和dim（matr）[2]（当在右侧时）？ - kutyw
如果没有2/0或0/2，它将在左侧实际显示“-Inf”，在右侧显示“Inf”。当然，这可以轻松地重新编码为您需要的任何内容，但我不确定您希望在这种情况下显示什么。左侧只有NA，右侧是dim（matr）[2]吗？如果最后一列确实有2/0呢？ - sparrow
不要返回NA或dim(matr)[2]右侧的0左侧的翻译，只返回翻译后的文本。 - kutyw
好的，我进行了编辑，如果两侧没有2/0或0/2，则返回NA。 - sparrow
0

这是一个有趣的问题。以下是我的解答。

首先定义一个向量，其中包含每列的乘积：

a3 <- matr[1,]*matr[2,]

然后，我们可以很容易地找到具有（0/2）或（2/0）对的列，因为我们知道矩阵只能包含值0、1和2：

the02s <- which(colSums(matr)==2 & a3==0)

接下来我们想要找到距离特定列数最近的(0/2)或(2/0)的对，分别在该列的左侧和右侧。例如，该列的编号可能是9：

thecol <- 9

现在我们基本上已经拥有了所有需要的东西，以找到最接近列thecol的(0/2)或(2/0)组合的索引（矩阵中的列数）。我们只需要使用findInterval()的输出即可：
pos <- findInterval(thecol,the02s)
pos <- c(pos, pos+1)
pos[pos==0] <- NA # output NA if no column was found on the left

结果如下：

the02s[pos]
#  2 14

因此，在满足所需条件的thecol两侧最接近的列的索引在此情况下为2和14，我们可以确认这些列号都包含其中一个相关组合：

matr[,14]
#a1 a2 
# 0  2
matr[,2]
#a1 a2 
# 0  2

编辑：我更改了答案，以便在矩阵左侧和/或右侧不存在满足所需条件的列时，返回NA。

但是这个想法是，给定一个输入列（比如9），找到第一个同时包含（0,2）或（2,0）的左侧和右侧的列。 - slamballais
@Laterow 再次感谢您指出我误解了 OP。我编辑了答案，现在我认为它提供了所需的输出。 - RHertel
0
如果您需要多次执行此操作，请预先计算所有位置。
loc <- which((a1==2 & a2==0) | (a1==0 & a2==2))

然后您可以使用findInterval函数找到左右两侧的第一个值。

i<-findInterval(9,loc);loc[c(i,i+1)]
# [1]  2 14

请注意，如果您想指定多个目标列，则findInterval是矢量化的。
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接