生成固定值的非负(或正)随机整数总和

3
我想随机分配正整数G组,使它们加起来等于V
例如,如果G = 3V = 21,则有效的结果可能是(7, 7, 7)(10, 6, 5)等。
有没有直接的方法来做到这一点?

编辑注意事项(来自李哲源):

如果值不受整数限制,则问题很简单,已在选择固定和的n个数字中解决。

对于整数,有一个先前的Q&A:在R中生成总和为M的N个随机整数,但似乎更复杂且难以理解。那里的基于循环的解决方案也不尽如人意。

3个回答

5

非负整数

n 为样本大小:

x <- rmultinom(n, V, rep.int(1 / G, G))

这是一个G x n矩阵,其中每一列都是多项式样本,总和为V

通过将rep.int(1 / G, G)传递给参数prob,我假设每个组具有相等的“成功”概率。


正整数

正如Gregor所提到的,多项式样本可以包含0。如果不希望出现这种样本,则应该将其拒绝。因此,我们从截断的多项式分布中进行采样。

如何在拒绝准则下从分布中生成目标数量的样本中,我建议采用“过度采样”方法来实现截断采样的“向量化”。简单来说,知道接受概率后,我们可以估计预期的试验次数M来看到第一个“成功”(非零)。我们首先采样大约1.25 * M个样本,然后这些样本中至少会有一个“成功”。我们随机返回其中一个作为输出。

以下函数实现了这个想法,以生成没有0的截断多项式样本。

positive_rmultinom <- function (n, V, prob) {
  ## input validation
  G <- length(prob)
  if (G > V) stop("'G > V' causes 0 in a sample for sure!")
  if (any(prob < 0)) stop("'prob' can not contain negative values!")
  ## normalization
  sum_prob <- sum(prob)
  if (sum_prob != 1) prob <- prob / sum_prob
  ## minimal probability
  min_prob <- min(prob)
  ## expected number of trials to get a "success" on the group with min_prob
  M <- round(1.25 * 1 / min_prob)
  ## sampling
  N <- n * M
  x <- rmultinom(N, V, prob)
  keep <- which(colSums(x == 0) == 0)
  x[, sample(keep, n)]
  }

现在让我们尝试一下。
V <- 76
prob <- c(53, 13, 9, 1)

直接使用rmultinom绘制样本有时会导致其中一些样本为0:

## number of samples that contain 0 in 1000 trials
sum(colSums(rmultinom(1000, V, prob) == 0) > 0)
#[1] 355   ## or some other value greater than 0

但是使用positive_rmultinom不存在这样的问题:

## number of samples that contain 0 in 1000 trials
sum(colSums(positive_rmultinom(1000, V, prob) == 0) > 0)
#[1] 0

2
也许有更便宜的方法,但这种方法似乎可行。
G <- 3
V <- 21
m <- data.frame(matrix(rep(1:V,G),V,G))
tmp <- expand.grid(m) # all possibilities
out <- tmp[which(rowSums(tmp) == V),] # pluck those that sum to 'V'
out[sample(1:nrow(out),1),] # randomly select a column

不确定如何使用runif


@Gregor 哎呀!combn 的确是个好点子。至于递增顺序,随机抽取列数不就解决了吗? - Brian Davis
啊,我明白你的意思了。我已经编辑了我的回答。 - Brian Davis
1
看起来好多了!另一个建议,虽然有点复杂但可以显著提高效率。让你的 m 有数值 rep(1:(V - 1), G - 1)(只有 G - 1 列)。然后,在你展开网格之后,你可以添加最后一列以完成总和 tmp = cbind(tmp, V - rowSums(tmp))。你仍然需要取出添加的列小于1的行。但是,从 expand.grid 结果中有效地删除一个维度将有助于使它能够扩展到比原来更大的规模。 - Gregor Thomas

0

我找到了一个我认为更简单的解决方案。首先,从最小值到最大值范围内生成随机整数,计算它们的数量,然后创建一个包括零的计数向量。

请注意,即使最小值大于零,这个解决方案可能仍然包含零。

希望这对未来遇到这个问题的人有所帮助 :)

rand.vect.with.total <- function(min, max, total) {
  # generate random numbers
  x <- sample(min:max, total, replace=TRUE)
  # count numbers
  sum.x <- table(x)
  # convert count to index position
  out = vector()
  for (i in 1:length(min:max)) {
    out[i] <- sum.x[as.character(i)]
  }
  out[is.na(out)] <- 0
  return(out)
}

rand.vect.with.total(0, 3, 5)
# [1] 3 1 1 0

rand.vect.with.total(1, 5, 10)
#[1] 4 1 3 0 2

注意,我也在这里发布了在R中生成总和为M的N个随机整数,但这个答案对两个问题都适用。


网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接