如何在R中生成等效于neg_binomial_2_rng的(mu,phi)参数化负二项分布随机数？

当然可以实现。Stan的neg_binomial_2_rng(mu, phi)基于Gamma-Poisson混合的参数化：先从Gamma分布采样速率参数λ，再从泊松分布生成最终的负二项随机数，该分布满足均值为mu，方差为mu + mu²/phi。在R里有两种简单的实现方式：

方法1：直接使用基础包的rnbinom函数 #

R的stats::rnbinom()函数支持以size和mu为参数的负二项分布采样，其中size参数恰好对应Stan中的phi（两者满足相同的方差公式：Var(Y) = mu + mu²/size）。直接调用即可：

# 生成1000个参数为mu=5，phi=2的负二项随机数
n <- 1000
mu <- 5
phi <- 2

y_rnbinom <- rnbinom(n = n, size = phi, mu = mu)

# 验证均值和方差（理论值：均值=5，方差=5 + 5²/2=17.5）
mean(y_rnbinom)
var(y_rnbinom)

方法2：手动实现Gamma-Poisson混合（完全对齐Stan底层逻辑） #

Stan的neg_binomial_2_rng本质是先采样λ ~ Gamma(shape=phi, rate=phi/mu)，再采样y ~ Poisson(λ)。你可以手动复现这个过程：

# 第一步：从Gamma分布采样λ
lambda <- rgamma(n = n, shape = phi, rate = phi / mu)
# 第二步：从泊松分布生成最终随机数
y_mixed <- rpois(n = n, lambda = lambda)

# 两种方法生成的随机数分布一致（随机性导致数值不完全相同）
hist(y_rnbinom, col = rgb(1,0,0,0.5), main = "两种方法生成的随机数分布")
hist(y_mixed, col = rgb(0,0,1,0.5), add = TRUE)
legend("topright", legend = c("rnbinom", "Gamma-Poisson混合"), fill = c(rgb(1,0,0,0.5), rgb(0,0,1,0.5)))

关键参数约束 #

• mu：负二项分布的均值，必须大于0
• phi：离散参数，值越大，分布越接近泊松分布（离散程度越低），必须大于0

内容的提问来源于stack exchange，提问作者Gled Guri