如何根据拒绝标准从分布中生成目标样本数

How to generate target number of samples from a distribution under a rejection criterion

我正在尝试 rnbinom 如下

x<- rnbinom(500, mu = 4, size = .1)
xtrunc <- x[x>0]

然后我只得到 125 个观察值。

但是,我想在相同条件 (mu = 4, size =.1) 下进行 500 次不包括 0(零)的观察。

这样就可以了:

N <- 500    ## target number of samples

## set seed for reproducibility
set.seed(0)
## first try
x <- rnbinom(N, mu = 4, size = .1)
p_accept <- mean(success <- x > 0)  ## estimated probability of accepting samples
xtrunc <- x[success]
## estimated further sampling times
n_further <- (N - length(xtrunc)) / p_accept
## further sampling
alpha <- 1.5   ## inflation factor
x_further <- rnbinom(round(alpha * n_further), mu = 4, size = .1)
## filter and combine
xtrunc <- c(xtrunc, (x_further[x_further > 0])[seq_len(N - length(xtrunc))])

## checking
length(xtrunc)
# [1] 500

summary(xtrunc)
#   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
#   1.00    2.00    5.00   12.99   16.00  131.00 

在上面,采样分为两个阶段。初始阶段的结果用于估计接受率的概率,以指导第二阶段的抽样。

但是,由于基础分布是明确已知的,因此接受率的理论概率是已知的。因此,在这种情况下无需执行两阶段方法。尝试:

p <- 1 - pnbinom(0, mu = 4, size = .1)  ## theoretical probability
alpha <- 1.5
n_try <- round(alpha * N / p)
set.seed(0)
x <- rnbinom(n_try, mu = 4, size = .1)
xtrunc <- (x[x > 0])[1:N]

## checking
length(xtrunc)
# [1] 500

summary(xtrunc)
#   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
#   1.00    2.00    5.00   12.99   16.00  131.00

背后的思想是几何分布理论。 密切相关。阅读 "More efficient vectorized method" 部分以获得详细说明。