对R中向量的所有唯一枚举进行置换

|| 我正在尝试找到一个可以置换向量的所有唯一置换的函数，同时不计算同一元素类型的子集中的并置数。例如：

dat <- c(1,0,3,4,1,0,0,3,0,4)

已

factorial(10)
> 3628800

可能的排列，但只有10!/(2!*2!*4!*2!)

factorial(10)/(factorial(2)*factorial(2)*factorial(2)*factorial(4))
> 18900

忽略同一元素类型的子集中的并置时的唯一排列。我可以通过combinat包中的unique()和permn()函数获得此功能

unique( permn(dat) )

但这在计算上非常昂贵，因为它涉及枚举n!，这可能比我需要的排列高出一个数量级。有没有一种方法，而无需先计算n!？

已邀请:

5 个回复

翰冒绢县

编辑：这是一个更快的答案；再次基于Louisa Gray和Bryce Wagner的思想，但由于更好地使用了矩阵索引，因此R代码更快。它比我原来的要快很多：

> ddd <- c(1,0,3,4,1,0,0,3,0,4)
> system.time(up1 <- uniqueperm(d))
   user  system elapsed 
  0.183   0.000   0.186 
> system.time(up2 <- uniqueperm2(d))
   user  system elapsed 
  0.037   0.000   0.038

和代码：

uniqueperm2 <- function(d) {
  dat <- factor(d)
  N <- length(dat)
  n <- tabulate(dat)
  ng <- length(n)
  if(ng==1) return(d)
  a <- N-c(0,cumsum(n))[-(ng+1)]
  foo <- lapply(1:ng, function(i) matrix(combn(a[i],n[i]),nrow=n[i]))
  out <- matrix(NA, nrow=N, ncol=prod(sapply(foo, ncol)))
  xxx <- c(0,cumsum(sapply(foo, nrow)))
  xxx <- cbind(xxx[-length(xxx)]+1, xxx[-1])
  miss <- matrix(1:N,ncol=1)
  for(i in seq_len(length(foo)-1)) {
    l1 <- foo[[i]]
    nn <- ncol(miss)
    miss <- matrix(rep(miss, ncol(l1)), nrow=nrow(miss))
    k <- (rep(0:(ncol(miss)-1), each=nrow(l1)))*nrow(miss) + 
               l1[,rep(1:ncol(l1), each=nn)]
    out[xxx[i,1]:xxx[i,2],] <- matrix(miss[k], ncol=ncol(miss))
    miss <- matrix(miss[-k], ncol=ncol(miss))
  }
  k <- length(foo)
  out[xxx[k,1]:xxx[k,2],] <- miss
  out <- out[rank(as.numeric(dat), ties=\"first\"),]
  foo <- cbind(as.vector(out), as.vector(col(out)))
  out[foo] <- d
  t(out)
}

它不会返回相同的顺序，但是排序之后，结果是相同的。

up1a <- up1[do.call(order, as.data.frame(up1)),]
up2a <- up2[do.call(order, as.data.frame(up2)),]
identical(up1a, up2a)

对于我的第一次尝试，请参阅编辑历史记录。

悲帽慑彤电

下面的函数（就像您在问题中手动完成的那样，它实现了重复排列的经典公式）对我来说似乎非常快：

upermn <- function(x) {
    n <- length(x)
    duplicates <- as.numeric(table(x))
    factorial(n) / prod(factorial(duplicates))
}

它确实计算n!，但不像permn函数那样先生成所有排列。实际观看：

> dat <- c(1,0,3,4,1,0,0,3,0,4)
> upermn(dat)
[1] 18900
> system.time(uperm(dat))
   user  system elapsed 
  0.000   0.000   0.001

更新：我刚刚意识到问题是关于生成所有唯一排列的问题，而不仅仅是指定它们的数量-抱歉！您可以通过为较少的一个元素指定唯一的排列，然后在其前面添加uniqe元素，来改进ѭ17部分。好吧，我的解释可能会失败，所以让消息人士说：

uperm <- function(x) {
u <- unique(x)                    # unique values of the vector
result <- x                       # let\'s start the result matrix with the vector
for (i in 1:length(u)) {
    v <- x[-which(x==u[i])[1]]    # leave the first occurance of duplicated values
    result <- rbind(result, cbind(u[i], do.call(rbind, unique(permn(v)))))
}
return(result)
}

这样您可以提高速度。我懒于在您提供的向量上运行代码（花了很多时间），这是在较小的向量上进行的比较：

> dat <- c(1,0,3,4,1,0,0)
> system.time(unique(permn(dat)))
   user  system elapsed 
  0.264   0.000   0.268 
> system.time(uperm(dat))
   user  system elapsed 
  0.147   0.000   0.150

我认为您可以通过将此函数重写为递归来获得更多收益！更新（再次）：我试图用我有限的知识来构造一个递归函数：

uperm <- function(x) {
    u <- sort(unique(x))
    l <- length(u)
    if (l == length(x)) {
        return(do.call(rbind,permn(x)))
    }
    if (l == 1) return(x)
    result <- matrix(NA, upermn(x), length(x))
    index <- 1
    for (i in 1:l) {
        v <- x[-which(x==u[i])[1]]
        newindex <- upermn(v)
        if (table(x)[i] == 1) {
            result[index:(index+newindex-1),] <- cbind(u[i], do.call(rbind, unique(permn(v))))
            } else {
                result[index:(index+newindex-1),] <- cbind(u[i], uperm(v))
            }
        index <- index+newindex
    }
    return(result)
}

这有很大的收获：

> system.time(unique(permn(c(1,0,3,4,1,0,0,3,0))))
   user  system elapsed 
 22.808   0.103  23.241 

> system.time(uperm(c(1,0,3,4,1,0,0,3,0)))
   user  system elapsed 
  4.613   0.003   4.645

如果适合您，请举报！

徐百晴墓斜

此处未提及的一个选项是multicool包中的allPerm函数。可以很容易地使用它来获取所有唯一的排列：

library(multicool)
perms <- allPerm(initMC(dat))
dim(perms)
# [1] 18900    10
head(perms)
#      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
# [1,]    4    4    3    3    1    1    0    0    0     0
# [2,]    0    4    4    3    3    1    1    0    0     0
# [3,]    4    0    4    3    3    1    1    0    0     0
# [4,]    4    4    0    3    3    1    1    0    0     0
# [5,]    3    4    4    0    3    1    1    0    0     0
# [6,]    4    3    4    0    3    1    1    0    0     0

在基准测试中，我发现它在ѭ25上比OP和daroczig的解决方案要快，但比Aaron的解决方案要慢。

倾坞髓

我实际上并不知道R，但是这是解决问题的方法：查找每种元素类型有多少，即

按频率排序（上面已经是）。从最频繁的值开始，该值占10个点中的4个。确定10个可用点内4个值的唯一组合。（0,1,2,3），（0,1,2,4），（0,1,2,5），（0,1,2,6） ...（0,1,2,9），（0,1,3,4），（0,1,3,5） ...（6,7,8,9）转到第二个最频繁的值，它占用6个可用点中的2个，并确定它是6个中的2个的唯一组合。（0,1），（0,2），（0,3），（0,4），（0,5），（1,2），（1,3）...（4,6），（5,6）然后是2之4：（0,1），（0,2），（0,3），（1,2），（1,3），（2,3）剩下的值2之2：（0,1）然后，您需要将它们组合成每种可能的组合。这是一些伪代码（我确信这是一种更有效的算法，但这应该不会太糟）：

lookup = (0,1,3,4)
For each of the above sets of combinations, example: input = ((0,2,4,6),(0,2),(2,3),(0,1))
newPermutation = (-1,-1,-1,-1,-1,-1,-1,-1,-1,-1)
for i = 0 to 3
  index = 0
  for j = 0 to 9
    if newPermutation(j) = -1
      if index = input(i)(j)
        newPermutation(j) = lookup(i)
        break
      else
        index = index + 1

超可林

另一个选择是iterpc封装，我相信它是现有方法中最快的。更重要的是，结果按字典顺序排列（可能更可取）。

dat <- c(1, 0, 3, 4, 1, 0, 0, 3, 0, 4)
library(iterpc)
getall(iterpc(table(dat), order=TRUE))

基准测试表明，iterpc比此处描述的所有其他方法快得多

library(multicool)
library(microbenchmark)
microbenchmark(uniqueperm2(dat), 
               allPerm(initMC(dat)), 
               getall(iterpc(table(dat), order=TRUE))
              )

Unit: milliseconds
                                     expr         min         lq        mean      median
                         uniqueperm2(dat)   23.011864   25.33241   40.141907   27.143952
                     allPerm(initMC(dat)) 1713.549069 1771.83972 1814.434743 1810.331342
 getall(iterpc(table(dat), order = TRUE))    4.332674    5.18348    7.656063    5.989448
          uq        max neval
   64.147399   74.66312   100
 1855.869670 1937.48088   100
    6.705741   49.98038   100

要回复问题请先登录或注册

对R中向量的所有唯一枚举进行置换

5 个回复

发起人

combinatorics

algorithm

r

permutation

问题状态

对R中向量的所有唯一枚举进行置换

与内容相关的链接

5 个回复

发起人

combinatorics

algorithm

r

permutation

问题状态