如何在 R 中模拟鞅过程问题？

Question

100个人在看戏。演出结束时，所有人都到衣帽间去拿他们的外套。衣帽间工作的人完全随机地把人们的外套还给他们。他们的参与者选对了外套就走。选错了的人把外套还给他，男的又随机还给了外套。当剧院的所有顾客都拿回了他们正确的外套时，这个过程就结束了。

我想在 R 中模拟这个鞅过程，以便找到这个过程结束的预期时间。 但我不知道如何。有什么帮助吗？ 比如：

# 100 customers
x = seq(1,100,by=1);x
# random sample from x 
y = sample(x,100,replace=FALSE)
x==y
# for the next iteration exclude those how are TRUE and run it again until everyone is TRUE

预期时间是需要多少次迭代。

或者类似的东西：

n = 100
X = seq(1,100,by=1)
martingale = rep(NA,n)

iterations = 0
accept     = 0
while (X != n) {
  iterations =  iterations + 1
  y = sample(1:100,100,replace=FALSE)
  if (X = y){ 
    accept = accept + 1
    X = X+1
    martingale [X] = y
  }
}
accept
iterations

Answer 1

一种方法如下（以10人为例，打印语句是不必要的，只是为了展示每次迭代做了什么）：

set.seed(0)
x <- 1:10
count <- 0
while(length(x) > 0){
  x <- x[x != sample(x)]
  print(x)
  count <- count + 1
}

# [1]  1  2  3  4  5  6  7  9 10
# [1] 3 4 5 6 7 9
# [1] 3 4 5 6 7
# [1] 3 4 5 6 7
# [1] 3 4 5 6 7
# [1] 3 4 5 6 7
# [1] 3 4 5 6 7
# [1] 3 4 5 6 7
# [1] 3 6
# 
count
# [1] 10

对于循环中的每一步，它都会删除客户被随机分配外套的 x 值，直到没有剩余。

要使用此代码获取 100 人的预期时间，您可以将其扩展为：

set.seed(0)
nits <- 1000 #simulate the problem 1000 times
count <- 0
for (i in 1:nits){
  x <- 1:100
  while(length(x) > 0){
    x <- x[x != sample(x)]
    count <- count + 1/nits
  } 
}
count
# [1] 99.901

我在没有证据的情况下假设 n 人的预期时间是 n 次迭代 - 当我尝试 50、100 或 200 人时，它似乎非常接近。

Answer 2

我没有关注你上面的讨论，我不完全确定这是否是你想要的，但我的理由如下：

• 您有 N 个人并排队。
• 在第一轮比赛中，第一人有 1/N 的机会获得正确的衣服。
• 此时您有两个选择。 Eitehr 人 1 的衣服穿对与否。
• 如果第 1 个人穿对了衣服，那么第 2 个人有 1/(N-1) 的机会穿对了衣服。如果第 1 个人没有穿上正确的衣服，第 1 个人会留在池中（最后），第 2 个人也有 1/N 的概率可以穿上正确的衣服。
• 您继续分配这些概率，直到所有 N 个人都见过店员一次。然后你挑选出那些穿对了衣服的人，并在第 1 步重复，直到每个人都穿对了。
• 出于模拟目的，您当然会重复整个过程 1000 或 10000 次。

如果我的理解正确，您对迭代次数感兴趣，即店员必须多久完成一次整个队列（或剩下的队列），直到每个人都穿好衣服。

图书馆（tidyverse）

people <- 100
results <- data.frame(people     = 1:people,
                      iterations = NA)

counter <- 0
finished <- 0

while (finished < people)
{
  loop_people <- results %>%
    filter(is.na(iterations)) %>%
    pull(people)

  loop_prob <- 1/length(loop_people)
  loop_correct <- 0

  for (i in 1:length(loop_people))
  {
    correct_clothes_i <- sample(c(0,1), size = 1, prob = c(1-loop_prob, loop_prob))
    if (correct_clothes_i == 1)
    {
      results[loop_people[i], 2] <- counter + 1
      loop_correct <- loop_correct + 1
      loop_prob <- 1/(length(loop_people) - loop_correct)
    }
  }
  counter <- counter + 1
  finished <- length(which(!is.na(results$iterations)))
}

max(results$iterations)

[1] 86

head(results)

  people iterations
1      1          7
2      2         42
3      3         86
4      4         67
5      5          2
6      6          9

results$iterations 列包含每个人正确穿衣服的迭代次数，因此max(results$iterations) 为您提供循环总数。

我没有证据，但根据经验和直觉，所需的迭代次数应该接近 N。