R expand.grid用于组中向量的重复组合?

cwxwcias  于 2023-10-13  发布在  其他
关注(0)|答案(2)|浏览(78)

我已经研究了两天了,没有任何进展。假设我有20个号码,我想(不替换)

  • 得到一组唯一的10个数字
  • 得到两组唯一的3个数
  • 得到两组唯一的2个数

共有(20选10)(10选3)(7选3)(4选2)(2选2)= 4,655,851,200个组。
我想要一个嵌套列表或data.frame,我可以评估这些组。在R中可以快速地实现吗?
下面我有一个非常笨拙的解决方案,当我只是试图从10个数字中得到一个唯一的5组,一个唯一的3组和一个唯一的2组时。
这段代码实际上是有效的,但很快就不适用于更大的群体。我想知道是否有一个快速的方法来创建这个网格?看起来应该很容易,但我不知道怎么做。提前感谢!!

library(data.table)
library(gtools)

n <- 10
group_sizes <- c(5,3,2)

# initialize first group
group_3_values <-
  combinations(
    n = n,
    r = group_sizes[1], # size of the first group!
    repeats.allowed = T,
    v = 1:n
  )
group_3_values <- as.data.table(group_3_values)
group_3_values[, row := 1:nrow(group_3_values)]
output <- as.data.table(group_3_values)

# run for loop for the remaining groups
start_time <- Sys.time()
for (group_size in group_sizes[-1]) { # size of the second, etc. groups
  output_list <- list()
  for (i in 1:nrow(output)) {
    remaining_options <-
      setdiff(1:n, output[row == i, .SD, .SDcols = !c('row')])
    second_group <- combinations(
      n = length(remaining_options),
      r = group_size,
      repeats.allowed = F,
      v = remaining_options
    )
    second_group <- as.data.table(second_group)
    second_group[, row := i]
    out <-
      merge(
        output[row == i,],
        second_group,
        by = "row",
        all.x = T,
        all.y = T,
        allow.cartesian = T
      )
    out$row <- NULL
    output_list[[i]] <- as.data.table(out)
    out <- NULL
  }
  output <- as.data.table(rbindlist(output_list))
  output[, row := 1:nrow(output)]
}
end_time <- Sys.time()
duration <- end_time - start_time
print(duration)
ctehm74n

ctehm74n1#

我认为实体化一个表,甚至一个包含4,655,851,200个元素的向量都不是一个好的选择。也许您能做的最好的事情就是使用combn来生成每个组合,并使用回调来运行您想要的代码。
这个例子生成并打印每个组合的元素,并在全局范围内增加一个计数器n。我使用callCC在10次迭代中提前退出。

x <- 1:20

n<-0
callCC(function(exit)
  combn(x, 10, function(i)
    combn(setdiff(x, c(i)), 3, function(j)
      combn(setdiff(x, c(i,j)), 3, function(k)
        combn(setdiff(x, c(i,j,k)), 2, function(l)
          combn(setdiff(x, c(i,j,k,l)), 2, function(m){
            n<<-n+1
            # YOUR CODE HERE
            print(paste0(sapply(list(i,j,k,l,m), paste0, collapse=" "), collapse=" | "))
            if(n>=10) exit("early exit")
          }, simplify=F),
        simplify=F),
      simplify=F),
    simplify=F),
  simplify=F))

#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 16 | 17 18 | 19 20"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 16 | 17 19 | 18 20"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 16 | 17 20 | 18 19"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 16 | 18 19 | 17 20"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 16 | 18 20 | 17 19"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 16 | 19 20 | 17 18"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 17 | 16 18 | 19 20"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 17 | 16 19 | 18 20"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 17 | 16 20 | 18 19"
#> [1] "1 2 3 4 5 6 7 8 9 10 | 11 12 13 | 14 15 17 | 18 19 | 16 20"
#> [1] "early exit"
dauxcl2d

dauxcl2d2#

这个问题与Create Combinations in R by Groups非常相似。唯一的区别是,链接问题要求每个组的大小相同。虽然这是一个很小的变化,但从算法的Angular 来看,它具有巨大的影响。
正如@Ric所指出的,这在硬件和计算时间方面将是困难的。
在我们开始之前,重要的是要注意实际上有1,163,962,800总结果,因为OP给出的公式没有考虑组的顺序。作为一个小例子,考虑生成2组2。
为了得到正确的数字,我们必须考虑到相同大小的组。在我们的例子中,我们有2组3和2组2,所以我们必须将原始结果除以2! * 2! = 4

(choose(20, 10) * choose(10, 3) *
    choose(7, 3) * choose(4, 2) *choose(2, 2)) /
(factorial(2) * factorial(2))
#> [1] 1163962800

我们减少了输出的大小,这很好,但我们仍然在谈论超过10亿个结果。

输入RcppAlgos

我是RcppAlgos软件包的作者,从2.8.0版本开始,我们可以处理任何大小的组的分区,从2.8.2版本开始,我们可以通过comboGroupsIter来处理它们。
如果你的机器可以处理它,我们可以很快地用函数comboGroups生成所有的结果,并使用参数nThreads并行生成结果。不用说,你需要一个像样的内存块:

library(RcppAlgos)

desired_grps <- c(10, 3, 3, 2, 2)
v <- seq_len(sum(desired_grps))
comboGroupsCount(v, grpSizes = desired_grps)
#> [1] 1163962800

## If you really wanted to generate them all,
## you will need at least 87 GBs free. The formula used is:
## ((# rows) * (# columns) * (bytes per cell)) / (# bytes in Gibabyte)
(comboGroupsCount(v, grpSizes = desired_grps) * 20 * 4) / 2^30
#> [1] 86.72199

或者,我们可以使用参数lowerupper,就像我们在Create Combinations in R by Groups的答案中所做的那样,以保持内存较低,但这可能会变得很麻烦。
一个更好的选择是使用迭代器来降低内存。这些迭代器非常灵活,提供了各种遍历结果的方法,甚至通过[[提供随机访问:

## nThreads is optional, however we choose 4 threads for increased efficiency
it <- comboGroupsIter(v, grpSizes = desired_grps, nThreads = 4)

## To see the next result, call the method nextIter
it@nextIter()
#> Grp1 Grp1 Grp2 Grp2 Grp3 Grp3 Grp3 Grp4 Grp4 Grp4 Grp5 Grp5 Grp5 Grp5 Grp5 Grp5 
#>    1    2    3    4    5    6    7    8    9   10   11   12   13   14   15   16 
#> Grp5 Grp5 Grp5 Grp5 
#>   17   18   19   20

## To see the next n results, call nextNIter
it@nextNIter(5)
#>      Grp1 Grp1 Grp2 Grp2 Grp3 Grp3 Grp3 Grp4 Grp4 Grp4 Grp5 Grp5 Grp5 Grp5 Grp5
#> [1,]    1    2    3    4    5    6    7    8    9   11   10   12   13   14   15
#> [2,]    1    2    3    4    5    6    7    8    9   12   10   11   13   14   15
#> [3,]    1    2    3    4    5    6    7    8    9   13   10   11   12   14   15
#> [4,]    1    2    3    4    5    6    7    8    9   14   10   11   12   13   15
#> [5,]    1    2    3    4    5    6    7    8    9   15   10   11   12   13   14
#>      Grp5 Grp5 Grp5 Grp5 Grp5
#> [1,]   16   17   18   19   20
#> [2,]   16   17   18   19   20
#> [3,]   16   17   18   19   20
#> [4,]   16   17   18   19   20
#> [5,]   16   17   18   19   20

## Note, the state is updated
it@summary()
#> $description
#> [1] "Partition of v of length 20 into 5 groups of sizes: 2, 2, 3, 3, 10"
#> 
#> $currentIndex
#> [1] 6
#> 
#> $totalResults
#> [1] 1163962800
#> 
#> $totalRemaining
#> [1] 1163962794

## Quickly jump to a specific point, m, without generating the m - 1 results before it
it[[1e8]]
#> Grp1 Grp1 Grp2 Grp2 Grp3 Grp3 Grp3 Grp4 Grp4 Grp4 Grp5 Grp5 Grp5 Grp5 Grp5 Grp5 
#>    1   10    3   12    5   14   19    6   18   20    2    4    7    8    9   11 
#> Grp5 Grp5 Grp5 Grp5 
#>   13   15   16   17

## Again, the state is updated
it@summary()
#> $description
#> [1] "Partition of v of length 20 into 5 groups of sizes: 2, 2, 3, 3, 10"
#> 
#> $currentIndex
#> [1] 100000000
#> 
#> $totalResults
#> [1] 1163962800
#> 
#> $totalRemaining
#> [1] 1063962800

## To see the last result
it@back()
#> Grp1 Grp1 Grp2 Grp2 Grp3 Grp3 Grp3 Grp4 Grp4 Grp4 Grp5 Grp5 Grp5 Grp5 Grp5 Grp5 
#>   17   20   18   19   11   15   16   12   13   14    1    2    3    4    5    6 
#> Grp5 Grp5 Grp5 Grp5 
#>    7    8    9   10

这些算法也非常高效,可以在合理的时间内生成所有1,163,962,800,特别是当我们使用nThreads

system.time({
    ## Reset the iterator
    it@startOver()

    ## Generate one million at a time
    while (!is.null(it@nextNIter(1e6))) {}
    print(it@summary())
})
#> $description
#> [1] "Partition of v of length 20 into 5 groups of sizes: 2, 2, 3, 3, 10"
#> 
#> $currentIndex
#> [1] 1163962801
#> 
#> $totalResults
#> [1] 1163962800
#> 
#> $totalRemaining
#> [1] -1
#>    user  system elapsed 
#>  25.271   3.668   8.771

我们不仅能够在10秒内遍历所有1,163,962,800结果,而且每次迭代只需要round((1e6 * 20 * 4) / 2^20) ~= 76兆字节,这也让我们的RAM很满意。

相关问题