── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr 1.1.4 ✔ readr 2.1.6
✔ forcats 1.0.1 ✔ stringr 1.6.0
✔ ggplot2 4.0.1 ✔ tibble 3.3.1
✔ lubridate 1.9.4 ✔ tidyr 1.3.2
✔ purrr 1.2.1
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
sim_sn_pois <- function(n, lambda){
sum(rpois(n, lambda = lambda))
}
You can add options to executable code like this
n_grid <- c(2, 5, 10)
lambda_grid <- c(0.5, 1.0, 1.5)
grid <- expand.grid(n = n_grid, lambda = lambda_grid) %>%
as_tibble()
pois_sims <- grid %>%
mutate(mu = n * lambda) %>%
slice(rep(1:n(), each = 10000)) %>% # repeat each (n,lambda) row 10,000 times
mutate(trial = rep(1:10000, times = nrow(grid))) %>%
mutate(sn = map2_dbl(n, lambda, sim_sn_pois))
pois_sims
# A tibble: 90,000 × 5
n lambda mu trial sn
<dbl> <dbl> <dbl> <int> <dbl>
1 2 0.5 1 1 0
2 2 0.5 1 2 0
3 2 0.5 1 3 1
4 2 0.5 1 4 2
5 2 0.5 1 5 1
6 2 0.5 1 6 1
7 2 0.5 1 7 1
8 2 0.5 1 8 0
9 2 0.5 1 9 1
10 2 0.5 1 10 0
# ℹ 89,990 more rows
pois_cdf <- pois_sims %>%
arrange(n, lambda, sn) %>%
group_by(n, lambda, mu) %>%
mutate(
Fhat = row_number() / n(),
F = ppois(sn, lambda = mu)
) %>%
ungroup()
head(pois_cdf)
# A tibble: 6 × 7
n lambda mu trial sn Fhat F
<dbl> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
1 2 0.5 1 1 0 0.0001 0.368
2 2 0.5 1 2 0 0.0002 0.368
3 2 0.5 1 8 0 0.0003 0.368
4 2 0.5 1 10 0 0.0004 0.368
5 2 0.5 1 11 0 0.0005 0.368
6 2 0.5 1 12 0 0.0006 0.368
ggplot(pois_cdf, aes(x = sn)) +
geom_step(aes(y = F, color = "analytic cdf")) +
geom_step(aes(y = Fhat, color = "empirical cdf")) +
facet_grid(n ~ lambda, labeller = label_both) +
labs(
title = "simulated vs analytic cdfs for sn = sum of poisson variables",
x = "sn",
y = "cdf",
color = ""
)
if you notice you’ll see the empirical and analytic cdf are the exact same supporting the result that from the question before
