- Even though the problem give us our lambda
- Lets build function to find lambda
- Then build function to find probability
library(kableExtra)
library(knitr)
## Build functions
find_lambda <- function(probability,N){
N*(probability)
}
poisson_prob <- function(lambda,x_occurences){
poisson_prob <- (lambda**(x_occurences)*exp(1)**(-lambda))/(factorial(x_occurences))
return(poisson_prob)
}
## calculate value at 5
probability <- 1/500
N=600
x_occurences <- 5
my_lambda <- find_lambda(probability,N)
my_poisson_probability <- poisson_prob(lambda=my_lambda,x_occurences)
paste("My poisson lambda is 1.2"," and the percent chance ",x_occurences," Occurences is ",my_poisson_probability*100,"%")
## [1] "My poisson lambda is 1.2 and the percent chance 5 Occurences is 0.624556317821142 %"
But wait a minute
- Our question wants to know the probability of 5 or more
- this is equivalent to the probability of 1 - x={0,1,2,3,4}
probability <- 1/500
N=600
x_occurences <- 5
my_lambda <- find_lambda(probability,N)
cumulative_range <- 0:4
total_probability <- 0
for (number in cumulative_range){
probability <- poisson_prob(my_lambda,number)
total_probability <- probability+total_probability
}
via_poisson_formula <- 1-total_probability
paste("The percent chance of getting 5 or more via poisson formula is",via_poisson_formula*100,"%")
## [1] "The percent chance of getting 5 or more via poisson formula is 0.774578827644146 %"
What if we simulate?
#library(yuima)
samples <- rpois(10000000, 1.2)
#hist_btimes_c <- replicate(100000,(runif(1)*runif(1)))
hist(samples,probability = TRUE)

P = ecdf(samples)
plot(P)

via_simulation <- 1-P(4)
paste("The percent chance of getting 5 or more via simulation is",via_simulation*100,"%")
## [1] "The percent chance of getting 5 or more via simulation is 0.77623 %"
compare results
round(via_poisson_formula,3)==round(via_simulation,3)
## [1] TRUE
Lucky cookie eaters
round(500*via_poisson_formula,0)
## [1] 4
- about 4 lucky(unlucky?) cookie eaters will get 5 or more raisins,