better visualization of uncertainty distributions

[markjrieke]

66/95% credible intervals are fine, but can get over-plotted for close states. May be better do do something like what ggblend does here:

If this requires all the draws for each day, then this is a non-starter, but I may be able to get around it by saving the results of ppoints() to a tibble per day.

[markjrieke]

Here's an example using ggblend that renders faster than ggdist's stat_lineribbon --- note that the rng-based example generates choppier distributions than the pure theoretical distribution Matthew Kay uses as an example in the ggblend docs.

library(tidyverse)
library(riekelib)
library(ggblend)

tibble(x = 1:180) %>%
  mutate(med = 0.2 * (x - 90)/180,
         sigma = 0.2 * (180 - x)/180 + 0.05,
         draws = pmap(list(med, sigma), ~rnorm(1e4, ..1, ..2))) %>%
  unnest(draws) %>%
  mutate(across(c(med, draws), expit)) %>%
  rename(d = draws) %>%
  mutate(r = 1 - d) %>%
  pivot_longer(c(d, r),
               names_to = "party",
               values_to = "p") %>%
  nest(data = p) %>%
  mutate(.width = map(data, ~seq(from = 0, to = 1, length.out = 101)),
         .lower = pmap(list(data, .width), ~quantile(..1$p, probs = (1 - ..2)/2)),
         .upper = pmap(list(data, .width), ~quantile(..1$p, probs = (1 - ..2)/2 + ..2))) %>%
  unnest(c(.width, .lower, .upper)) %>%
  ggplot(aes(x = x,
             y = med,
             ymin = .lower,
             ymax = .upper,
             fill = party,
             group = paste(party, .width))) + 
  geom_ribbon(alpha = 0.01) %>% partition(vars(party)) %>% blend("multiply") +
  scale_fill_manual(values = c("royalblue", "red2")) + 
  theme_rieke() +
  expand_limits(y = c(0, 1))

^{Created on 2024-10-25 with reprex v2.1.0}

[markjrieke]

For (potentially) multi-modal distributions, this is a good alternative:

library(tidyverse)
library(riekelib)
library(ggblend)

tibble(x = 1:180) %>%
  mutate(med1 = 0.2 * (x - 90)/180 - 0.2,
         med2 = 0.6 * (x - 90)/180 + 0.1,
         sigma = 0.2 * (180 - x)/180 + 0.1,
         draws1 = pmap(list(med1, sigma), ~rnorm(1e4, ..1, ..2)),
         draws2 = pmap(list(med2, sigma), ~rnorm(1e4, ..1, ..2)),
         draws = pmap(list(draws1, draws2), ~sample(c(..1, ..2), 1e4))) %>%
  select(-c(draws1, draws2)) %>%
  unnest(draws) %>%
  mutate(draws = expit(draws)) %>%
  rename(d = draws) %>%
  mutate(r = 1 - d) %>%
  pivot_longer(c(d, r),
               names_to = "party",
               values_to = "p") %>%
  mutate(e = round(p * 538)) %>%
  select(x, party, e) %>%
  count(x, party, e) %>%
  ggplot(aes(x = x,
             y = e,
             fill = party,
             alpha = n)) + 
  geom_raster(interpolate = TRUE) %>% partition(vars(party)) %>% blend("multiply") +
  scale_fill_manual(values = c("royalblue", "red2")) + 
  scale_alpha_continuous(range = c(0, 0.5)) + 
  expand_limits(x = c(0, 180),
                y = c(0, 538)) +
  theme_rieke()
#> Warning: Raster pixels are placed at uneven horizontal intervals and will be shifted
#> ℹ Consider using `geom_tile()` instead.

^{Created on 2024-10-25 with reprex v2.1.0}