I keep getting error messages saying that I am creating hundreds of NA's in my for loop. Where do those NA come from? Any help would be greatly appreciated!
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()
You get integer overflow. Try
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
As soon as lambda is too large, 32-bit representation of integer is insufficient and NA is returned. (Recall that Poisson random variables are integers).
Your loop has a dynamic growth on rate (lambda), which can eventually become too big. Running your function with a smaller minutes, say 10, is fine.
By contrast, ppois and dpois which produce double-precision floating point numbers are fine with large 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%, wherexis a random value from the uniform distribution on the interval[0, 5].
The rate increases by x% not x. So you should use
rate <- rate * (1 + runif(1, 0, 5) / 100)