为什么我在模拟泊松随机变量的 "for" 循环中得到很多 NA

Question

我不断收到错误消息，说我在 for 循环中创建了数百个 NA。那些 NA 来自哪里？任何帮助将不胜感激！

drip <- function(rate = 1, minutes = 120) {
  count <- 0
  for(i in 1:(minutes)) {
    count <- count + rpois(1, rate)
    rate <- rate * runif(1, 0, 5)
  }
  count
}
drip()

Answer 1

你得到了整数溢出。尝试

set.seed(0)
rpois(1, 1e+8)
#[1] 100012629
rpois(1, 1e+9)
#[1] 999989683
rpois(1, 1e+10)
#[1] NA
#Warning message:
#In rpois(1, 1e+10) : NAs produced

一旦lambda过大，整数的32位表示不足，返回NA。（回想一下，泊松随机变量是整数）。

您的循环在 rate (lambda) 上动态增长，最终会变得太大。运行你的函数有一个较小的 minutes，比如 10，没问题。

相比之下，产生双精度浮点数的ppois和dpois可以很好地处理大lambda。

dpois(1e+8, 1e+8)
#[1] 3.989423e-05
dpois(1e+9, 1e+9)
#[1] 1.261566e-05
dpois(1e+10, 1e+10)
#[1] 3.989423e-06
dpois(1e+11, 1e+11)
#[1] 1.261566e-06

ppois(1e+8, 1e+8)
#[1] 0.5000266
ppois(1e+9, 1e+9)
#[1] 0.5000084
ppois(1e+10, 1e+10)
#[1] 0.5000027
ppois(1e+11, 1e+11)
#[1] 0.5000008

With each passing minute, the rate parameter increases by x%, where x is a random value from the uniform distribution on the interval [0, 5].

rate 增加 x% 而不是 x。所以你应该使用

rate <- rate * (1 + runif(1, 0, 5) / 100)

为什么我在模拟泊松随机变量的 "for" 循环中得到很多 NA

Why do I get many NA's in a "for" loop that simulates Poisson random variables

loops

for-loop

r

poisson

na