As with length, it’s probably not really interesting to look at the distribution of penis circumferences regardless of state, but well, we’re a thorough bunch.

ggplot(penis_long, aes(x = circumf)) +
  geom_histogram(binwidth = .5, alpha = .7, position = "dodge") +
  scale_x_continuous(breaks = seq(0, 100, 1), 
                     minor_breaks = seq(0, 100, .5),
                     labels = label_cm,
                     sec.axis = sec_axis(trans = ~./2.54, 
                                         labels = label_in)) +
  labs(title = "World Penis Data", subtitle = "Circumference",
       x = "Circumference", y = "Frequency", 
       fill = "State", color = "State",
       caption = plot_caption) +
  theme(legend.position = "top")

By State

ggplot(penis_long, aes(x = circumf, fill = state, color = state)) +
  geom_histogram(binwidth = .5, alpha = .7, position = "dodge") +
  scale_x_continuous(breaks = seq(0, 100, 1), 
                     minor_breaks = seq(0, 100, .5),
                     labels = label_cm,
                     sec.axis = sec_axis(trans = ~./2.54, 
                                         labels = label_in)) +
  scale_fill_brewer(palette = "Paired") +
  scale_color_brewer(palette = "Paired") +
  labs(title = "World Penis Data", subtitle = "Circumference by State",
       x = "Circumference", y = "Frequency", 
       fill = "State", color = "State",
       caption = plot_caption) +
  theme(legend.position = "top")

By State and Region

penis_long %>%
  group_by(Region) %>%
  mutate(order_circumf = mean(circumf)) %>%
  {
    ggplot(., aes(x = reorder(Region, order_circumf), y = circumf, 
                  fill = state, color = state)) +
    geom_boxplot(alpha = .7) +
    coord_flip() +
    scale_y_continuous(breaks = seq(0, 100, 2), 
                       minor_breaks = seq(0, 100, .5),
                       labels = label_cm,
                       sec.axis = sec_axis(trans = ~./2.54, 
                                           labels = label_in)) +
    scale_fill_brewer(palette = "Paired") +
    scale_color_brewer(palette = "Paired") +
    labs(title = "World Penis Data", subtitle = "Circumference by State",
         x = "", y = "Circumference", fill = "State", color = "State",
         caption = plot_caption) +
    theme(axis.text.y = element_text(size = rel(1.2)),
          legend.position = "top")
  }

# CIs
penis_long %>%
  group_by(Region) %>%
  mutate(order_circumf = mean(circumf)) %>%
  {
    ggplot(data = ., aes(x = reorder(Region, order_circumf), y = circumf, color = state)) +
    stat_summary(fun.data = mean_cl_normal, geom = "errorbar") +
    stat_summary(fun.y = mean, geom = "point", size = 2) +
    coord_flip() +
    scale_y_continuous(breaks = seq(0, 100, 1), 
                       minor_breaks = seq(0, 100, .5),
                       labels = label_cm,
                       sec.axis = sec_axis(trans = ~./2.54, 
                                           labels = label_in)) +
    scale_color_brewer(palette = "Paired") +
    labs(title = "World Penis Data", subtitle = "Circumference by State",
         x = "", y = "Circumference", fill = "State", color = "State",
         caption = plot_caption) +
    theme(axis.text.y = element_text(size = rel(1.2)),
          legend.position = "top")
  }

By Length

Is length a good predictor of circumference?

ggplot(data = penis_long, aes(x = length, y = circumf, color = state)) +
  geom_smooth(method = lm, se = F, color = "gray") +
  geom_smooth(method = lm, se = F) +
  geom_point(size = 3, color = "black") +
  geom_point(size = 2) +
  scale_x_continuous(labels = label_cm,
                     sec.axis = sec_axis(trans = ~./2.54, 
                                         labels = label_in)) +
  scale_y_continuous(labels = label_cm,
                     sec.axis = sec_axis(trans = ~./2.54, 
                                         labels = label_in)) +
  scale_color_brewer(palette = "Paired") +
  labs(title = "World Penis Data",
       subtitle = "Circumference by Length",
       x = "Length", y = "Circumference", color = "State",
       caption = plot_caption) +
  theme(legend.position = "top")

Choropleth-Penismap

A worldmapwise look at penile circumference distribution.

highchart() %>% 
  hc_add_series_map(worldgeojson, map, value = "circumf_erect", joinBy = "iso3",
                    name = "Penis Circumference (erect)") %>% 
  hc_title(text = "Erect Penis Circumference by Country") %>% 
  hc_tooltip(valueDecimals = 2, valuePrefix = "<b>", valueSuffix = " cm</b>", 
             useHTML = TRUE) %>% 
  hc_colorAxis(stops = idk, min = 10, max = 14) %>% 
  hc_legend(enabled = TRUE) %>% 
  hc_mapNavigation(enabled = TRUE)