如何根据拒绝标准从分布中生成目标样本数
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" 部分以获得详细说明。
我正在尝试 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
背后的思想是几何分布理论。