7.2 homework

library(purrrfect)

## 
## Attaching package: 'purrrfect'

## The following objects are masked from 'package:base':
## 
##     replicate, tabulate

library(tidyverse)

## ── 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

order_stat_sim <- parameters(
  ~n,
  c(4, 8, 12, 16)
) %>%
  add_trials(10000)

order_stat_sim <- order_stat_sim %>%
  mutate(
    sample = map(n, ~ runif(.x)),
    y1 = map_dbl(sample, min),
    y3 = map_dbl(sample, ~ sort(.x)[3]),
    yn = map_dbl(sample, max)
  )

order_stat_long <- order_stat_sim %>%
  select(n, y1, y3, yn) %>%
  pivot_longer(
    cols = c(y1, y3, yn),
    names_to = "stat",
    values_to = "value"
  )

beta_df <- expand_grid(
  n = c(4, 8, 12, 16),
  stat = c("y1", "y3", "yn"),
  x = seq(0, 1, length.out = 400)
) %>%
  mutate(
    k = case_when(
      stat == "y1" ~ 1,
      stat == "y3" ~ 3,
      stat == "yn" ~ n
    ),
    density = dbeta(x, shape1 = k, shape2 = n - k + 1)
  )

ggplot(order_stat_long, aes(x = value)) +
  geom_histogram(aes(y = after_stat(density)),
                 bins = 40,
                 fill = "lavender",
                 alpha = 0.6) +
  geom_line(
    data = beta_df,
    aes(x = x, y = density),
    color = "purple",
    linewidth = 1
  ) +
  facet_grid(stat ~ n, scales = "free_y") +
  theme_classic(base_size = 14) +
  labs(
    title = "Simulated vs analytic densities for order statistics",
    x = "value",
    y = "density"
  )