代码之家 › 专栏 › 技术社区 › compbiostats

生成与固定值相加的非负(或正)随机整数

random r

compbiostats · 技术社区 · 7 年前

到 G V .

例如,如果 G = 3 V = 21 ,有效结果可能是 (7, 7, 7) , (10, 6, 5)

有没有直接的方法?

编辑通知(从 æå²æº ):

如果值不局限于整数,那么问题就很简单,在 Choosing n numbers with fixed sum .

对于整数,前面有一个Q&答: Generate N random integers that sum to M in R 但它看起来更复杂,也很难理解。基于循环的解决方案也不令人满意。

2 回复 | 直到 7 年前

Zheyuan Li 7 年前

非负整数

让 n

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

是一个 G x n 矩阵,其中每列都是 multinomial 总计为 V

路过 rep.int(1 / G, G) prob 我假设每个小组“成功”的概率相等。

作为 Gregor 提到,多项式样本可以包含0。如果这些样品是不需要的,则应拒收。因此,我们从截断多项式分布中取样。

How to generate target number of samples from a distribution under a rejection criterion 我建议使用“过采样”方法来实现截断采样的“矢量化”。简单地说,知道接受概率我们就可以估计出预期的试验次数 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

## 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

Brian Davis 7 年前

可能是一种比较便宜的方法,但这似乎是可行的。

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

Jakey 5 年前

我找到了一个更简单的解决办法。首先生成从最小值到最大值的随机整数,对它们进行计数,然后生成一个计数向量(包括零)。

希望这能帮助未来的研究人员解决这个问题:)

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

注意,我也在这里贴了这个 Generate N random integers that sum to M in R