r - 生成具有非重复元素的向量的多个排列

Question

我有一个向量:

seq1<-c('a','b','c','b','a','b','c','b','a','b','c')

我希望置换此向量的元素以创建多个(理想情况下最多 5000 个)向量，条件是置换的向量不能在向量中的连续元素中包含重复元素。例如“abbca ....” 是不允许的，因为 'b-b' 是重复的。

我意识到对于这个小例子，可能没有 5000 个解决方案。我通常处理更大的向量。我也愿意考虑更换 sample ，尽管目前我正在研究无需更换的解决方案。

我正在寻找比我目前的想法更好的解决方案。

选项 1. - 蛮力。

在这里，我只是反复采样并检查是否有任何连续的元素是重复的。

set.seed(18)
seq1b <-  sample(seq1a)
seq1b
#[1] "b" "b" "a" "a" "c" "b" "b" "c" "a" "c" "b"
sum(seq1b[-length(seq1b)]==seq1b[-1])  #3

这不是解决方案，因为有 3 个重复的连续元素。我也意识到 lag可能是检查重复元素的更好方法，但由于某种原因它很挑剔(我认为它被我加载的另一个包所掩盖)。

set.seed(1000)
res<-NULL
for (i in 1:10000){res[[i]]<-sample(seq1a)}
res1 <- lapply(res, function(x) sum(x[-length(x)]==x[-1]))
sum(unlist(res1)==0) #228

这会在 10000 次迭代中产生 228 个选项。但让我们看看有多少独特的:

res2 <- res[which(unlist(res1)==0)]
unique(unlist(lapply(res2, paste0, collapse="")))  #134

在 10000 次尝试中，我们只能从这个简短的示例向量中获得 134 次独特的尝试。

以下是生成的 134 个示例序列中的 3 个:

# "bcbabcbabca" "cbabababcbc" "bcbcababacb"

事实上，如果我尝试超过 500,000 个样本，我只能得到 212 个符合我的非重复标准的独特序列。这可能接近可能的上限。

选项 2. - 迭代

我的第二个想法是对方法进行更多的迭代。

seq1a
table(seq1a)
#a b c 
#3 5 3

我们可以采样这些字母之一作为我们的起点。然后从剩余的中抽取另一个，检查它是否与先前选择的相同，如果不相同，则将其添加到末尾。等等等等……

set.seed(10)
newseq <- sample(seq1a,1)  #b
newseq #[1] "b"

remaining <-seq1a[!seq1a %in% newseq | duplicated(seq1a)]
table(remaining)
#a b c 
#3 4 3 

set.seed(10)
newone <- sample(remaining,1) #c

#check if newone is same as previous one.
newone==newseq[length(newseq)] #FALSE
newseq <- c(newseq, newone) #update newseq
newseq #[1] "b" "c"

remaining <-seq1a[!seq1a %in% newseq | duplicated(seq1a)] #update remaining
remaining
table(remaining)

#a b c 
#3 4 2 

这可能有效，但我也可以看到它遇到了很多问题 - 例如我们可以去:

# "a" "c" "a" "c" "a" "b"  ...

然后再留下 3 个 'b's，因为它们是重复的，所以最后不能去。

当然，如果我允许带替换采样，这会容易得多，但现在我正在尝试不替换采样。

Answer 1

您可以使用 iterpc包以使用组合和迭代。在尝试回答这个问题之前我没有听说过它，所以可能还有更有效的方法来使用相同的包。

这里我用过 iterpc设置迭代器，以及 getall根据该迭代器查找向量的所有组合。这似乎只是报告独特的组合，这比使用 expand.grid 查找所有组合要好一些。 .

#install.packages("iterpc")
require("iterpc")

seq1 <- c('a','b','c','b','a','b','c','b','a','b','c')

I <- iterpc(n = table(seq1), ordered=TRUE)

all_seqs <- getall(I)

# result is a matrix with permutations as rows:
head(all_seqs)
#     [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
#[1,] "a"  "a"  "a"  "b"  "b"  "b"  "b"  "b"  "c"  "c"   "c"  
#[2,] "a"  "a"  "a"  "b"  "b"  "b"  "b"  "c"  "b"  "c"   "c"  
#[3,] "a"  "a"  "a"  "b"  "b"  "b"  "b"  "c"  "c"  "b"   "c"  
#[4,] "a"  "a"  "a"  "b"  "b"  "b"  "b"  "c"  "c"  "c"   "b"  
#[5,] "a"  "a"  "a"  "b"  "b"  "b"  "c"  "b"  "b"  "c"   "c"  
#[6,] "a"  "a"  "a"  "b"  "b"  "b"  "c"  "b"  "c"  "b"   "c" 

rle函数告诉我们向量中彼此相等的连续值。 lengths输出的分量告诉我们 values的每个元素的次数。重复:

rle(c("a", "a", "b", "b", "b", "c", "b"))

# Run Length Encoding
#   lengths: int [1:3] 2 3 1 1
#   values : chr [1:3] "a" "b" "c" "b"

values的长度或 lengths仅对于没有连续重复的组合才等于原始向量的长度。

因此，您可以申请 rle到每一行，计算 values的长度或 lengths并保留来自 all_seqs 的行其中计算值与 seqs1的长度相同.

#apply the rle function 
all_seqs_rle <- apply(getall(I), 1, function(x) length(rle(x)$values))

# keep rows which have an rle with a length equal to length(seq1)
all_seqs_good <- all_seqs[which(all_seqs_rle == length(seq1)), ]

all_seqs_good有一个 nrow 212，表明您确实找到了示例向量的所有可能组合。

nrow(all_seqs_good)
# 212 

从技术上讲，这仍然是蛮力(除了它不会计算所有可能的组合 - 只有唯一的组合)，但对于您的示例来说相当快。我不确定它会如何处理更大的向量......

编辑:对于较大的向量，这似乎确实失败了。一种解决方案是将较大的向量分解为较小的块，然后按上述方式处理这些块并组合它们 - 仅保留符合您标准的组合。

例如，将长度为 24 的向量分解为长度为 12 的两个向量，然后组合结果可以为您提供 200,000 多个符合您的标准并且非常快的组合(对我来说大约 1 分钟):

# function based on the above solution
seq_check <- function(mySeq){
I = iterpc(n = table(mySeq), ordered=TRUE)
all_seqs <- getall(I)
all_seqs_rle <- apply(getall(I), 1, function(x) length(rle(x)$values))
all_seqs_good <- all_seqs[which(all_seqs_rle == length(mySeq)), ]
return(all_seqs_good)
}

set.seed(1)
seq1<-sample(c(rep("a", 8), rep("b", 8), rep("c", 8)),24)

seq1a <- seq1[1:12]
seq1b <- seq1[13:24]

#get all permutations with no consecutive repeats
seq1a = apply(seq_check(seq1a), 1, paste0, collapse="")
seq1b = apply(seq_check(seq1b), 1, paste0, collapse="")

#combine seq1a and seq1b: 
combined_seqs <- expand.grid(seq1a, seq1b)
combined_seqs <- apply(combined_seqs, 1, paste0, collapse="") 

#function to calculate rle lengths
rle_calc <- function(x) length(rle(unlist(strsplit(x, "")))$values)

#keep combined sequences which have rle lengths of 24
combined_seqs_rle <- sapply(combined_seqs, rle_calc)
passed_combinations <- combined_seqs[which(combined_seqs_rle == 24)]

#find number of solutions
length(passed_combinations)
#[1] 245832
length(unique(passed_combinations))
#[1] 245832

您可能需要对起始向量重新排序以获得最佳结果。例如，如果 seq1在上面的例子中，连续八次以“a”开头，没有通过的解决方案。例如，尝试使用 seq1 <- c(rep("a", 8), rep("b", 8), rep("c", 8)) 的拆分解决方案并且您没有得到任何解决方案，即使对于随机序列确实有相同数量的解决方案。

看起来您不需要找到所有可能的传递组合，但是如果您这样做了，那么对于更大的向量，您可能需要遍历 I使用 getnext函数来自 iterpc ，并在循环中检查每一个，这会非常慢。