---
title: "Visualizing Proportions"
output:
  html_document:
    toc: yes
    code_folding: show
    code_download: true
---

<link rel="stylesheet" href="stat4580.css" type="text/css" />
<style type="text/css"> .remark-code { font-size: 85%; } </style>
<!-- title based on Wilke's chapter -->

```{r setup, include = FALSE, message = FALSE}
source(here::here("setup.R"))
knitr::opts_chunk$set(collapse = TRUE, message = FALSE,
                      fig.height = 5, fig.width = 6, fig.align = "center")

set.seed(12345)
library(dplyr)
library(ggplot2)
library(lattice)
library(gridExtra)
library(patchwork)
source(here::here("datasets.R"))
```


## Categorical Data

Categorical data can be

* nominal, qualitative

* ordinal

For visualization, the main difference is that ordinal data suggests a
particular display order.

Purely categorical data can come in a range of formats.

The most common are

* raw data: individual observations;

* aggregated data: counts for each unique combination of levels;

* cross-tabulated data.


### Raw Data

```{r, echo = FALSE}
ah <- as.data.frame(HairEyeColor)
raw <- ah[rep(seq_len(nrow(ah)), times = ah$Freq), ][-4]
raw <- raw[sample(seq_len(nrow(raw))), ]
row.names(raw) <- NULL
```

Raw data for a survey of individuals that records hair color, eye
color, and gender of `r nrow(raw)` individuals might look like this:

```{r}
head(raw)
```


### Aggregated Data

One way to aggregate raw categorical data is to use `count` from `dplyr`:

```{r}
library(dplyr)
agg <- count(raw, Hair, Eye, Sex) |>
    as_tibble()
agg
```

<!--
The `count_` function from `dplyr` allows the variables to use to be
read from the data:

```{r, eval = FALSE}
agg <- count_(raw, names(raw))
head(agg)
```

Apparently the "modern" way to do this is

```{r}
count(raw, !!! syms(names(raw)))
```
-->


### Cross-Tabulated Data

Cross-tabulated data can be produced from aggregate data using `xtabs`:

```{r}
xtabs(n ~ Hair + Eye + Sex, data = agg)
```

Cross-tabulated data can be produced from raw data using `table`:

```{r}
xtb <- table(raw)
xtb
```

Both raw and aggregate data in this example are in _tidy_ form; the
cross-tabulated data is not.

Cross-tabulated data on $p$ variables is arranged in a $p$-way array.

The cross-tabulated data can be converted to the tidy aggregate form
using `as.data.frame`:

```{r}
class(xtb)
head(as.data.frame(xtb))
```

The variable `xtb` corresponds to the data set `HairEyeColor` in the
`datasets` package,


### Working With Categorical Variables

Categorical variables are usually represented as:

* character vectors

* factors.

Some advantages of factors:

* more control over ordering of levels

* levels are preserved when forming subsets

* levels can reflect possible values not present in the data

Most plotting and modeling functions will convert character vectors to
factors with levels ordered alphabetically.

Some standard R functions for working with factors include

* `factor` creates a factor from another type of variable
* `levels` returns the levels of a factor
* `reorder` changes level order to match another variable
* `relevel` moves a particular level to the first position as a base line
* `droplevels` removes levels not in the variable.

The `tidyverse` package `forcats` adds some more tools, including

* `fct_inorder` creates a factor with levels ordered by first appearance
* `fct_infreq` orders levels by decreasing frequency
* `fct_rev` reverses the levels
* `fct_recode` changes factor levels
* `fct_relevel` moves one or more levels
* `fct_c` merges two or more factors
* `fct_collapse` merge some factor levels


## Bar Charts For Frequencies


### Basics

A bar chart is often used to show the frequencies of a categorical
variable.

By default, `geom_bar` uses `stat = "count"` and maps its result to
the `y` aesthetic.

This is suitable for raw data:

```{r, class.source = "fold-hide"}
thm <- theme_minimal() +
    theme(text = element_text(size = 16))
ggplot(raw) +
    geom_bar(aes(x = Hair),
             fill = "deepskyblue3") +
    thm
```

For a nominal variable it is often better to order the bars by
decreasing frequency:

```{r, class.source = "fold-hide"}
library(forcats)
ggplot(mutate(raw,
              Hair = fct_infreq(Hair))) +
    geom_bar(aes(x = Hair),
             fill = "deepskyblue3") +
    thm
```

If the data have already been aggregated, then you need to either

* specify `stat = "identity"` as well as the variable containing the
  counts as the `y` aesthetic, or

* use `geom_col`:

```{r, class.source = "fold-hide"}
ggplot(agg) +
    geom_col(aes(x = Hair,
                 y = n),
             fill = "deepskyblue3") +
    thm
```

For aggregated data, reordering can be based on the computed counts
using

```{r}
agg_ord <-
    mutate(agg,
           Hair = reorder(Hair, -n, sum))
```

* `-n` is used to order largest to smallest;

* the default summary used by `reorder` is `mean`; `sum` is better here.

```{r, class.source = "fold-hide"}
ggplot(agg_ord) +
    geom_col(aes(x = Hair,
                 y = n),
             fill = "deepskyblue3") +
    thm
```


### Adding a Grouping Variable

Mapping the `Eye` variable to `fill` in `ggplot` produces a _stacked
bar chart_.

An alternative, specified with `position = "dodge"`, is a _side by
side_ bar chart, or a _clustered_ bar chart.

For the side by side chart in particular it may be useful to also
reorder the `Eye` color levels.

```{r, fig.width = 8, class.source = "fold-hide"}
ecols <- c(Brown = "brown4",
           Blue = "blue2",
           Hazel = "darkgoldenrod3",
           Green = "green4")
agg_ord <-
    mutate(agg,
           Hair = reorder(Hair, -n, sum),
           Eye = reorder(Eye, -n, sum))
p1 <- ggplot(agg_ord) +
    geom_col(aes(x = Hair,
                 y = n,
                 fill = Eye)) +
    scale_fill_manual(values = ecols) +
    thm
p2 <- ggplot(agg_ord) +
    geom_col(aes(x = Hair,
                 y = n,
                 fill = Eye),
             position = "dodge") +
    scale_fill_manual(values = ecols) +
    thm
(p1 + guides(fill = "none")) | p2
```

Faceting can be used to bring in additional variables:

```{r, fig.width = 8, class.source = "fold-hide"}
p1 + facet_wrap(~ Sex)
```

The counts shown here may not be the most relevant features for
understanding the joint distributions of these variables.


## Pie Charts and Doughnut Charts

_Pie charts_ go by many different names (from a [Twitter
thread](https://twitter.com/ElephantEating/status/1361039771414319106)):

```{r, echo = FALSE, out.width = "80%"}
knitr::include_graphics(IMG("pienames.jpeg"))
```

Pie charts can be viewed as stacked bar charts in polar coordinates:

```{r, fig.width = 8, class.source = "fold-hide"}
hcols <- c(Black = "black", Brown = "brown4",
           Red = "brown1", Blond = "lightgoldenrod1")
p1 <- ggplot(agg_ord) +
    geom_col(aes(x = 1, y = n, fill = Hair), position = "fill") +
    coord_polar(theta = "y") +
    scale_fill_manual(values = hcols) +
    thm
p2 <- ggplot(agg_ord) +
    geom_col(aes(x = Hair, y = n, fill = Hair)) +
    scale_fill_manual(values = hcols) +
    thm
(p1 + guides(fill = "none")) | p2
```

The axes and grid lines are not helpful for the pie chart and can be
removed with some _theme_ settings.

Using faceting we can also separately show the distributions for men
and women:

```{r, fig.width = 8, class.source = "fold-hide"}
pie_thm <- thm +
    theme(axis.title = element_blank(),
          axis.text = element_blank(),
          axis.ticks = element_blank(),
          panel.grid.major = element_blank(),
          panel.grid.minor = element_blank(),
          panel.border = element_blank())

p3 <- p1 + facet_wrap(~ Sex) + pie_thm
p3
```

_Doughnut charts_ are a variant that has recently become popular in
the media:

```{r, fig.width = 8, class.source = "fold-hide"}
p4 <- p3 + xlim(0, 1.5)
p4
```

The center is often used for annotation:

```{r, fig.width = 8, class.source = "fold-hide"}
p4 + geom_text(aes(x = 0, y = 0, label = Sex), size = 5) +
    theme(strip.background = element_blank(),
          strip.text = element_blank())
```

An alternative to the polar coordinates approach uses
`geom_arc_bar` and `stat_pie` from package `ggforce`:

```{r, fig.width = 8, class.source = "fold-hide"}
#| warning: false
library(ggforce)
arrange(agg_ord, desc(Hair)) |>
    ggplot(aes(x0 = 0, y0 = 0, r0 = 0, r = 1, amount = n, fill = Hair)) +
    geom_arc_bar(stat = "pie", color = NA) +
    coord_fixed() +
    scale_fill_manual(values = hcols) +
    pie_thm +
    facet_wrap(~ Sex)
```

For doughnut charts:

```{r, fig.width = 8, class.source = "fold-hide"}
arrange(agg_ord, desc(Hair)) |>
    ggplot(aes(x0 = 0, y0 = 0, r0 = 0.4, r = 1, amount = n, fill = Hair)) +
    geom_arc_bar(stat = "pie", color = NA) +
    geom_text(aes(x = 0, y = 0, label = Sex), size = 5) +
    coord_fixed() +
    scale_fill_manual(values = hcols) +
    pie_thm +
    theme(strip.background = element_blank(),
          strip.text = element_blank()) +
    facet_wrap(~ Sex)
```


## Some Notes

Pie charts are effective for judging part/whole relationships.

Pie charts can be effective for comparing proportions to

* one half
* one quarter

Pie charts are not very effective for comparing proportions to each other.

3D pie charts are popular and a very bad idea. An example
([Fig. 6.61](https://www.dropbox.com/s/tlehzi3kb6ikbsz/6.61.3DIllustration.png?dl=0))
from Andy Kirk's book (2016),
[_Data Visualization: A Handbook for Data Driven Design_](http://book.visualisingdata.com/home):

```{r, echo = FALSE, out.width = "50%"}
knitr::include_graphics(IMG("badpie.png"))
```

Pie charts are widely used for political data.

* With the right ordering, pie charts are very good at showing
  which coalitions of parties can form a majority.

* When no one candidate earns a majority of the votes, pie charts
  do not show which candidate has earned a plurality very well.

* Good orientation and factor ordering can help.

```{r, fig.width = 8, class.source = "fold-hide"}
elect <- geofacet::election |>
    group_by(candidate) |>
    summarize(votes = sum(votes))

p1 <- ggplot(elect) +
    geom_col(aes(x = 1, y = votes, fill = candidate), position = "fill") +
    coord_polar(theta = "y", start = -1) +
    xlim(c(-0.5, 1.5)) +
    scale_fill_manual(values = c(Trump = scales::muted("red", 50, 80),
                                 Clinton = scales::muted("blue", 50, 70),
                                 Other = "grey")) +
    pie_thm

p2 <- mutate(elect,
             candidate = factor(candidate,
                                c("Clinton", "Other", "Trump"))) |>
    ggplot() +
    geom_col(aes(x = 1, y = votes, fill = candidate), position = "fill") +
    coord_polar(theta = "y") +
    xlim(c(-0.5, 1.5)) +
    scale_fill_manual(values = c(Trump = scales::muted("red", 50, 80),
                                 Clinton = scales::muted("blue", 50, 70),
                                 Other = "grey")) +
    pie_thm

p3 <- ggplot(elect) +
    geom_col(aes(x = candidate,
                 y = 100 * (votes / sum(votes)),
                 fill = candidate)) +
    scale_fill_manual(values = c(Trump = scales::muted("red", 50, 80),
                                 Clinton = scales::muted("blue", 50, 70),
                                 Other = "grey")) +
    labs(y = "percent") +
    thm +
    theme(axis.text.x = element_blank(),
          axis.title.x = element_blank())

(p1 + guides(fill = "none")) + (p2 + guides(fill = "none")) + p3
```


## Some Alternatives


### Stacked Bar Charts

Stacked bar charts with equal heights, or filled bar charts, are an
alternative for representing part-whole relationships.

* Top and bottom proportions are easy to compare.

* Comparing proportions to one half and one quarter is harder.

```{r, class.source = "fold-hide"}
ggplot(agg) +
    geom_col(aes(x = Sex, y = n, fill = Hair), position = "fill") +
    scale_fill_manual(values = hcols) +
    thm
```


### Waffle Charts

Another alternative is a _waffle chart_, sometimes also called a
_square pie chart_.

```{r, echo = FALSE, out.width = "50%"}
knitr::include_graphics(IMG("waffle.png"))
```

The [`waffle`](https://github.com/hrbrmstr/waffle) package is one R
implementation of this idea.

Currently the development version on GitHub is needed for the
following examples.

Showing the counts:

```{r, fig.width = 7, class.source = "fold-hide"}
library(waffle)
stopifnot(packageVersion("waffle") >= "1.0.1")
ggplot(arrange(agg, Hair), aes(values = n, fill = Hair)) +
    geom_waffle(n_rows = 18, flip = TRUE, color = "white", size = 0.33,
                na.rm = FALSE) +
    coord_equal() +
    facet_wrap(~ Sex) +
    scale_fill_manual(values = hcols) +
    theme_minimal() +
    theme_enhance_waffle()
```

```{r, eval = TRUE, echo = FALSE}
ggplot(arrange(agg, Hair), aes(values = n, fill = Hair)) +
    geom_waffle(n_rows = 10, flip = TRUE, color = "white", size = 0.33,
                na.rm = FALSE) +
    coord_equal() +
    facet_grid(Sex ~ Eye) +
    scale_fill_manual(values = hcols) +
    theme_minimal() +
    theme_enhance_waffle()
```

Showing the proportions:

```{r, fig.width = 7, class.source = "fold-hide"}
round_pct <- function(n) {
    pct <- 100 * (n / sum(n))
    nn <- floor(pct)
    if (sum(nn) < 100) {
        rem <- pct - nn
        idx <- sort(order(rem), decreasing = TRUE)[seq_len(100 - sum(nn))]
        nn[idx] <- nn[idx] + 1
    }
    nn
}

group_by(agg, Sex) |>
    mutate(pct = round_pct(n)) |>
    ungroup() |>
    ggplot(aes(values = pct, fill = Hair)) +
    geom_waffle(n_rows = 10, flip = TRUE, color = "white", size = 0.33,
                na.rm = FALSE) +
    coord_equal() +
    facet_wrap(~ Sex) +
    scale_fill_manual(values = hcols) +
    theme_minimal() +
    theme_enhance_waffle()
```

```{r, eval = FALSE, echo = FALSE}
group_by(agg, Sex, Eye) |>
    mutate(pct = round_pct(n)) |>
    ungroup() |>
    arrange(pct, Hair) |>
    ggplot(aes(values = pct, fill = Hair)) +
    geom_waffle(n_rows = 10, flip = TRUE, color = "white", size = 0.33,
                na.rm = FALSE) +
    coord_equal() +
    facet_grid(Sex ~ Eye) +
    scale_fill_manual(values = hcols) +
    theme_minimal() +
    theme_enhance_waffle()
```


## Population Pyramids

Bar charts for two groups can be shown back to back.

```{r, class.source = "fold-hide"}
mutate(agg, Hair = reorder(Hair, n, sum)) |>
    ggplot(aes(x = ifelse(Sex == "Male", n, -n),
               y = Hair,
               fill = Sex)) +
    geom_col() +
    xlab("Count") +
    thm
```

This is often used for showing age distributions by sex for
populations; the result is called a [_population
pyramid_](https://www.visualcapitalist.com/us-population-pyramid-1980-2050/).

Age distribution data for many countries and years is available from a
[Census Bureau
website](https://www.census.gov/data-tools/demo/idb/).

Data files for 2026 for
[Germany](https://stat.uiowa.edu/~luke/data/germany-2026.csv) and
[Nigeria](https://stat.uiowa.edu/~luke/data/nigeria-2026.csv) are
available locally.

```{r, warning = FALSE, class.source = "fold-hide"}
library(dplyr)
library(forcats)
library(readr)
if (! file.exists("germany-2026.csv"))
    download.file("https://stat.uiowa.edu/~luke/data/germany-2026.csv",
                  "germany-2026.csv")
if (! file.exists("nigeria-2026.csv"))
    download.file("https://stat.uiowa.edu/~luke/data/nigeria-2026.csv",
                  "nigeria-2026.csv")

gm_pop <- read_csv("germany-2026.csv") |>
    select(Year, Age = GROUP, Population,
           Male.Population = "Male Population",
           Female.Population = "Female Population") |>
    filter(! is.na(Age)) |>
    mutate(Age = fct_inorder(Age))

ni_pop <- read_csv("nigeria-2026.csv") |>
    select(Year, Age = GROUP, Population,
           Male.Population = "Male Population",
           Female.Population = "Female Population") |>
    filter(! is.na(Age)) |>
    mutate(Age = fct_inorder(Age))
```

Combining the data sets allows a side by side comparison of the counts:

```{r, class.source = "fold-hide"}
library(tidyr)
pop2 <-
    bind_rows(mutate(gm_pop, Country = "Germany"),
              mutate(ni_pop, Country = "Nigeria")) |>
    select(Age,
           Male = Male.Population,
           Female = Female.Population, Country) |>
    pivot_longer(Male : Female,
                 names_to = "Sex",
                 values_to = "n")

ggplot(pop2) +
    geom_col(aes(x = ifelse(Sex == "Male", n, -n),
                 y = Age,
                 fill = Sex)) +
    facet_wrap(~ Country) +
    scale_x_continuous(
        labels = function(n)
        scales::label_number(scale_cut =
                                 scales::cut_short_scale())(abs(n))) +
    xlab("Count") +
    thm +
    theme(legend.position = "top")
```

The different shapes are evident, but are harder to see than
they could be because of the difference in total population:

```{r, class.source = "fold-hide"}
group_by(pop2, Country) |>
    summarize(Population = sum(n)) |>
    ungroup() |>
    mutate(Population = scales::comma(Population)) |>
    knitr::kable(format = "html", align = "lr") |>
    kableExtra::kable_styling(full_width = FALSE)
```

Using a group mutate we can compute sex/age group percentages within
each country:

```{r, class.source = "fold-hide"}
group_by(pop2, Country) |>
    mutate(pct = 100 * n / sum(n)) |>
    ungroup() |>
    ggplot() +
    geom_col(aes(x = ifelse(Sex == "Male", pct, -pct),
                 y = Age,
                 fill = Sex)) +
    facet_wrap(~ Country) +
    scale_x_continuous(
        labels = function(x) scales::percent(abs(x / 100))) +
    xlab("Percent") +
    thm +
    theme(legend.position = "top")
```

```{r, eval = FALSE, echo = FALSE}
p <- ggplot(gm_pop) +
    geom_col(aes(x = Male.Population, y = Age, fill = "Male")) +
    geom_col(aes(x = -Female.Population, y = Age, fill = "Female")) +
    thm

p | (p + ni_pop)
```


## Multiple Categorical Variables

Visualizing the distribution of multiple categorical variables
involves visualizing counts and proportions.

Distributions can be viewed as

* joint distributions;

* conditional distributions.

When one variable (or several) can be viewed as a response and others
as predictors then it is common to focus on the conditional
distribution of the response given the predictors.

The most common approaches use variants of bar and area charts.

The resulting plots are often called [_mosaic
plots_](https://en.wikipedia.org/wiki/Mosaic_plot).


## Two Data Sets


### Hair and Eye Color

```{r}
HairEyeColorDF <-
    as.data.frame(HairEyeColor)
head(HairEyeColorDF)
```

Marginal distributions of the variables:

```{r, fig.height = 3, fig.width = 9, class.source = "fold-hide"}
p1 <- ggplot(HairEyeColorDF) +
    geom_col(aes(Sex, Freq), fill = "deepskyblue3") +
    thm
p2 <- mutate(HairEyeColorDF, Hair = reorder(Hair, -Freq, sum)) |>
    ggplot() +
    geom_col(aes(Hair, Freq), fill = "deepskyblue3") +
    thm
p3 <- ggplot(HairEyeColorDF) +
    geom_col(aes(Eye, Freq), fill = "deepskyblue3") +
    thm
p1 | p2 | p3
```


### Arthritis Data

The `vcd` package includes the data frame `Arthritis` with several
variables for 84 patients in a clinical trial for a treatment for
rheumatoid arthritis.

```{r, message = FALSE}
data(Arthritis, package = "vcd")
head(Arthritis)
```

* The `Improved` variable is the response.

* The predictors are `Treatment`, `Sex`, and `Age`.

Counts for the categorical predictors:
```{r}
xtabs(~ Sex, Arthritis)
```

```{r}
xtabs(~ Treatment, Arthritis)
```

```{r}
xtabs(~ Treatment + Sex, data = Arthritis)
```

Joint distribution of the predictors:

```{r, class.source = "fold-hide"}
ggplot(Arthritis) +
    geom_histogram(aes(x = Age),
                   binwidth = 10,
                   fill = "deepskyblue3",
                   color = "black") +
    facet_grid(Treatment ~ Sex) +
    thm
```

Conditional distribuiton of age, given sex and treatment:

```{r, class.source = "fold-hide"}
ggplot(Arthritis) +
    geom_histogram(aes(x = Age,
                       y = after_stat(density)),
                   binwidth = 10,
                   fill = "deepskyblue3",
                   color = "black") +
    facet_grid(Treatment ~ Sex) +
    thm
```


## Bar Charts


### Hair and Eye Color

Default bar charts show the individual count or joint proportions.

For the hair-eye color aggregated data counts:

```{r, class.source = "fold-hide"}
ggplot(HairEyeColorDF) +
    geom_col(aes(x = Eye, y = Freq, fill = Sex)) +
    facet_wrap(~ Hair) +
    thm
```

Joint proportions:

```{r, class.source = "fold-hide"}
ggplot(mutate(HairEyeColorDF, Prop = Freq / sum(Freq))) +
    geom_col(aes(x = Eye, y = Prop, fill = Sex)) +
    facet_wrap(~ Hair) +
    thm
```

* Differing frequencies of the hair colors are visible.

* Conditional distributions of eye color within hair color are
  harder to compare.

Showing conditional distributions requires computing proportions
within groups.

For the joint conditional distribution of sex and eye color given hair color:

```{r, class.source = "fold-hide"}
group_by(HairEyeColorDF, Hair) |>
    mutate(Prop = Freq / sum(Freq)) |>
    ungroup() |>
    ggplot() +
    geom_col(aes(x = Eye, y = Prop, fill = Sex)) +
    facet_wrap(~ Hair) +
    thm
```

* It is easier to compare the skewness of the eye color
  distributions for black, brown, and red hair.

* Assessing the proportion of females or males withing the different
  groups is possible but challenging since it requires relative
  length comparisons.

To more clearly see the that the proportion of females among subjects
with blond hair and blue eyes is higher than for other hair/eye color
combinations we can look at the conditional distribution of sex given
hair and eye color.

```{r, class.source = "fold-hide"}
group_by(HairEyeColorDF, Hair, Eye) |>
    mutate(Prop = Freq / sum(Freq)) |>
    ungroup() |>
    ggplot() +
    geom_col(aes(x = Eye,
                 y = Prop,
                 fill = Sex)) +
    facet_wrap(~ Hair, nrow = 1) +
    thm +
    theme(axis.text.x =
              element_text(angle = 45,
                           hjust = 1))
```

This plot can also be obtained using `position = "fill"`.

```{r, eval = FALSE, class.source = "fold-hide"}
ggplot(HairEyeColorDF) +
    geom_col(aes(x = Eye,
                 y = Freq,
                 fill = Sex),
             position = "fill") +
    facet_wrap(~ Hair, nrow = 1) +
    thm +
    theme(axis.text.x =
              element_text(angle = 45,
                           hjust = 1))
```

One drawback: This visualization no longer shows that some of the
hair/eye color combinations are more common than others.


### Arthritis Data

For the raw arthritis data, `geom_bar` computes the aggregate counts
and produces a stacked bar chart by default:

```{r}
p <- ggplot(Arthritis, aes(x = Sex,
                           fill = Improved)) +
    facet_wrap(~ Treatment)
p + geom_bar() +
    scale_fill_brewer(palette = "Blues") +
    thm

```

Specifying `position = "dodge"` produces a side-by-side plot:

```{r, class.source = "fold-hide"}
p + geom_bar(position = "dodge") +
    scale_fill_brewer(palette = "Blues") +
    thm
```

There are no cases of male patients on placebo reporting `Some`
improvement, resulting in wider bars for the other options.

One way to produce a zero height bar:

* aggregate with `count`, and

* use `complete` from `tidyr`

```{r, class.source = "fold-hide"}
library(tidyr)
comp_counts <-
    count(Arthritis,
          Treatment, Sex, Improved) |>
    complete(Treatment, Sex, Improved,
             fill = list(n = 0))
ggplot(comp_counts,
       aes(x = Sex, y = n, fill = Improved)) +
    geom_col(position = "dodge") +
    facet_wrap(~ Treatment) +
    scale_fill_brewer(palette = "Blues") +
    thm
```

Another option is to use the `preserve = "single"` option with
`position_dodge`.

```{r, class.source = "fold-hide"}
p + geom_bar(position =
                 position_dodge(
                     preserve = "single")) +
    scale_fill_brewer(palette = "Blues") +
    thm
```

Showing conditional distributions of `Improved` given different levels
of `Treatment` and `Sex`:

```{r, class.source = "fold-hide"}
group_by(comp_counts, Treatment, Sex) |>
    mutate(prop = n / sum(n)) |>
    ungroup() |>
    ggplot() +
    geom_col(aes(x = Sex,
                 y = prop,
                 fill = Improved),
             position = "dodge") +
    facet_wrap(~ Treatment) +
    scale_fill_brewer(palette = "Blues") +
    thm
```

Stacked bar charts with height one are another option to make
these conditional distributions easier to compare:

```{r, class.source = "fold-hide"}
p + geom_bar(position = "fill") +
    scale_fill_brewer(palette = "Blues") +
    thm
```

Ordering of variables affects which comparisons are easier.

* A researcher might want to emphasize the differential response among
  males and females.

* A patient might prefer to be able to focus on whether the treatment
  is effective for them:

```{r, class.source = "fold-hide"}
ggplot(Arthritis, aes(x = Treatment, fill = Improved)) +
    geom_bar(position = "fill") +
    scale_fill_brewer(palette = "Blues") +
    thm +
    facet_wrap(~ Sex)
```

Some notes:

* The stacked bar chart is effective for two categories, and a few
  more if they are ordered.

* Providing a visual indication of uncertainty in the estimates is a
  challenge. The standard errors in this case are around 0.1.

* The proportions of each treatment group that are male or female
  could be encoded in the bar widths.

* The resulting plot is called a _spine plot_.

* Basic `ggplot2` does not seem to make this easy.


## Spine Plots

_Spine plots_ are a special case of [_mosaic
plots_](https://en.wikipedia.org/wiki/Mosaic_plot), and can be seen as
a generalization of stacked bar plots.

For a spine plot the proportions for the categories of a predictor
variable are encoded in the bar widths.

The `ggmosaic` package provides support for mosaic plots in the
`ggplot` framework. (It can be a little rough around the edges.)

Spine plots are provided by the base graphics function `spineplot` and
the `vcd` function `spine`.

`vcd` plots are built on the `grid` graphics system, like `lattice`
and `ggplot2` graphics.

<!-- spine plot in the wild:
https://t.co/91xqkWXIfD;
in img/energy-spine.png -->

A spine plot for the distribution of `Improved` given `Sex` in the
`Treated` group:

```{r, warning = FALSE, class.source = "fold-hide"}
library(ggmosaic)
filter(Arthritis, Treatment == "Treated") |>
    mutate(Improved = fct_rev(Improved)) |>
    ggplot() +
    geom_mosaic(aes(x = Sex,
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    facet_wrap(~ Treatment) +
    thm + labs(x = "", y = "Improved")
```

Spine plots for `Treatment` groups using faceting:

```{r, class.source = "fold-hide"}
library(ggmosaic)
mutate(Arthritis,
       Improved = fct_rev(Improved)) |>
    ggplot() +
    geom_mosaic(aes(x = Sex,
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    facet_wrap(~ Treatment) +
    thm + labs(x = "", y = "Improved")
```

Spine plots for the arthritis data, faceted on `Sex`:

```{r, class.source = "fold-hide"}
library(ggmosaic)
mutate(Arthritis,
       Improved = fct_rev(Improved)) |>
    ggplot() +
    geom_mosaic(aes(x = Treatment,
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    facet_wrap(~ Sex) +
    thm + labs(x = "", y = "Improved")
```
This no longer shows the Female/Male imbalance.

For aggregate counts use the weight aesthetic:

```{r, class.source = "fold-hide"}
mutate(HairEyeColorDF, Sex = fct_rev(Sex)) |>
    ggplot() +
    geom_mosaic(aes(weight = Freq,
                    x = Hair,
                    fill = Sex)) +
    thm + labs(x = "Hair", y = "")
```

Spine plots of `Sex` within `Eye` color, faceted on `Hair` color:

```{r, class.source = "fold-hide"}
mutate(HairEyeColorDF, Sex = fct_rev(Sex)) |>
    ggplot() +
    geom_mosaic(aes(weight = Freq,
                    x = Eye,
                    fill = Sex)) +
    thm + labs(x = "Eye", y = "") +
    facet_wrap(~ Hair,
               nrow = 1,
               scales = "free_x") +
    theme(legend.position = "top",
          axis.text.y = element_blank(),
          axis.text.x =
              element_text(angle = 45,
                           hjust = 1)) +
    scale_y_continuous(expand = c(0, 0))
```

The relative sizes of the groups on the `x` (eye color) axis are shown
within the facets.

The sizes of the faceted variable (hair color) groups are not reflected.

_Double decker plots_ try to address this.


## Doubledecker Plots

_Doubledecker plots_ can be viewed as a generalization of spine plots
to multiple predictors.

Package `vcd` provides the `doubledecker` function.

This function can use a formula interface.

```{r, class.source = "fold-hide"}
arth_pal <-
    RColorBrewer::brewer.pal(3, "Blues")
arth_gp <- grid::gpar(fill = arth_pal)
vcd::doubledecker(Improved ~ Treatment + Sex,
                  data = Arthritis,
                  gp = arth_gp,
                  margins = c(2, 5, 4, 2))
```

```{r, class.source = "fold-hide"}
vcd::doubledecker(Improved ~ Sex + Treatment,
                  data = Arthritis,
                  gp = arth_gp,
                  margins = c(2, 5, 4, 2))
```

Using `ggmosaic`:

```{r, class.source = "fold-hide"}
mutate(Arthritis,
       Improved = fct_rev(Improved)) |>
    ggplot() +
    geom_mosaic(
        aes(x = product(Sex, Treatment),
            fill = Improved),
        divider = ddecker()) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    thm +
    theme(axis.text.x =
              element_text(angle = 15,
                           hjust = 1)) +
    labs(x = "", y = "")
```

```{r, class.source = "fold-hide"}
mutate(Arthritis,
       Improved = fct_rev(Improved)) |>
    ggplot() +
    geom_mosaic(
        aes(x = product(Treatment, Sex),
            fill = Improved),
        divider = ddecker()) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    thm +
    theme(axis.text.x =
              element_text(angle = 15,
                           hjust = 1)) +
    labs(x = "", y = "")
```

Another approach using the
[`ggh4x`](https://teunbrand.github.io/ggh4x/) package:

```{r, class.source = "fold-hide"}
library(dplyr)
library(ggplot2)
library(ggh4x)
data(Arthritis, package = "vcd")
ggplot(Arthritis, aes(x = 0, fill = Improved)) +
    geom_bar(position = "fill") +
    scale_fill_brewer(palette = "Blues") +
    facet_nested(~ Treatment + Sex) +
    force_panelsizes(cols = count(Arthritis, Treatment, Sex)$n) +
    labs(x = NULL, y = NULL) +
    scale_x_continuous(expand = c(0, 0), breaks = NULL, labels = NULL) +
    scale_y_continuous(expand = c(0, 0), breaks = NULL, labels = NULL)
```

```{r, class.source = "fold-hide"}
library(dplyr)
library(ggplot2)
library(ggh4x)
data(Arthritis, package = "vcd")
ggplot(Arthritis, aes(x = 0, fill = Improved)) +
    geom_bar(position = "fill") +
    scale_fill_brewer(palette = "Blues") +
    facet_nested(~ Sex + Treatment) +
    force_panelsizes(cols = count(Arthritis, Sex, Treatment)$n) +
    labs(x = NULL, y = NULL) +
    scale_x_continuous(expand = c(0, 0), breaks = NULL, labels = NULL) +
    scale_y_continuous(expand = c(0, 0), breaks = NULL, labels = NULL)
```

## Mosaic Plots

_Mosaic plots_ recursively partition the axes to represent counts of
categorical variables as rectangles.

* Base graphics provides `mosaicplot`;

* `vcd` provides `mosaic`.

Both support a formula interface.

A Mosaic plot for the predictors `Sex` and `Treatment`:

```{r, class.source = "fold-hide"}
vcd::mosaic(~ Sex + Treatment,
            data = Arthritis)
```

Adding `Improved` to the joint distribution:

```{r, class.source = "fold-hide"}
vcd::mosaic(~ Sex + Treatment + Improved,
            data = Arthritis)
```

Identifying `Improved` as the response:

<!-- ##vcd::mosaic(Improved ~ Sex + Treatment, data = Arthritis, gp = arth_gp)-->

```{r}
vcd::mosaic(Improved ~ Sex + Treatment,
            data = Arthritis)
```

Matching the doubledecker plots:

```{r, class.source = "fold-hide"}
vcd::mosaic(
         Improved ~ Treatment + Sex,
         data = Arthritis,
         split_vertical = c(TRUE, TRUE, FALSE))
```

```{r, class.source = "fold-hide"}
vcd::mosaic(
         Improved ~ Sex + Treatment,
         data = Arthritis,
         split_vertical = c(TRUE, TRUE, FALSE))
```

Some variants using `ggmosaic`:

```{r, class.source = "fold-hide"}
ggplot(mutate(Arthritis, Sex = fct_rev(Sex))) +
    geom_mosaic(
        aes(x = product(Treatment, Sex))) +
    coord_flip() +
    labs(x = "", y = "")
```

```{r, class.source = "fold-hide"}
ggplot(mutate(Arthritis,
              Sex = fct_rev(Sex),
              Improved = fct_rev(Improved))) +
    geom_mosaic(aes(x = product(Improved,
                                Treatment,
                                Sex))) +
    coord_flip()
```

A mosaic plot for all bivariate marginals:

```{r, class.source = "fold-hide"}
pairs(xtabs(~ Sex + Treatment + Improved, data = Arthritis))
```


## Spinograms and CD Plots

<!-- building a spinogram from scratch, more or less:

```r
Arth <- mutate(Arthritis,
               AgeBin = cut(Age, seq(20, by = 10, len = 7)),
               Improved = fct_rev(Improved))
d <- filter(Arth, Treatment == "Treated") |>
    count(AgeBin) |>
    mutate(brk = (cumsum(lag(n, default = 0)) + 0.5 * n) / sum(n))
p <- ggplot(d)

p + geom_col(aes(x = AgeBin, y = n))

p2 <- p +
    geom_col(aes(x = 0.5, y = n, color = fct_rev(AgeBin)),
             position = "fill", fill = NA) +
    guides(color = "none") +
    scale_y_continuous(breaks = d$brk,
                       labels = d$AgeBin)
p2

p2 + coord_flip()
```

```r
props <- filter(Arth, Treatment == "Treated") |>
    count(AgeBin, Improved) |>
    mutate(AgeBin = fct_drop(AgeBin)) |>
    complete(AgeBin, Improved, fill = list(n = 0)) |>
    group_by(AgeBin) |>
    mutate(prop = n / sum(n)) |>
    ungroup() |>
    select(-n)
props
```
-->

_Spinograms_ and _CD plots_ show the conditional distribution of a
categorical variable given the value of a numeric variable.

* Spinograms use the same binning as a histogram and then create a
  spine plot.

* CD plots use a smoothing or density estimation approach.

A spinogram for `Improved` against `Age`:

```{r, class.source = "fold-hide"}
ArthT <- filter(Arthritis,
                Treatment == "Treated") |>
    mutate(Improved = fct_rev(Improved))
arthT_gp <-
    grid::gpar(fill = rev(arth_gp$fill))
vcd::spine(Improved ~ Age,
           data = ArthT,
           gp = arthT_gp,
           breaks = 5)
```

An analogous plot created with `ggmosaic` by binning the `Age` variable:

```{r, class.source = "fold-hide"}
Arth <-
    mutate(Arthritis,
           AgeBin = cut(Arthritis$Age,
                        seq(20, by = 10,
                            len = 7)),
           Improved = fct_rev(Improved))
filter(Arth, Treatment == "Treated") |>
    count(Improved, AgeBin) |>
    ggplot() +
    geom_mosaic(aes(weight = n,
                    x = AgeBin,
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    theme_minimal() +
    theme(axis.title = element_blank())
```

A facet grid can be used to create spinograms for each of the
`Sex`/`Treatment` combinations:

```{r, fig.width = 8, class.source = "fold-hide"}
ggplot(count(Arth, Improved, Sex, Treatment, AgeBin)) +
    geom_mosaic(aes(weight = n,
                    x = AgeBin,
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues",
                      direction = -1,
                      guide = guide_legend(reverse = TRUE)) +
    theme_minimal() +
    facet_grid(Treatment ~ Sex) +
    theme(axis.title = element_blank()) +
    theme(axis.text.x = element_text(angle = 35,
                                     hjust = 1),
          axis.text.y = element_blank())
```

A [spinogram](https://flowingdata.com/2021/08/02/decline-of-u-s-vaccination-rate-compared-against-europes/) in the media (NYT, August 2021):

```{r, echo = FALSE, out.width = "70%"}
knitr::include_graphics(IMG("Vaccination-rates-with-Marimekko-1536x925.png"))
```

Some plots in a [Twitter thread](https://twitter.com/jburnmurdoch/status/1503420660869214213?s=20&t=RZxOXiHZ0eXkV6WXMIuVVA):

```{r, echo = FALSE, out.width = "90%"}
knitr::include_graphics(IMG("covid-spinogram.png"))
```

```{r, echo = FALSE, out.width = "90%"}
knitr::include_graphics(IMG("covid-mort.png"))
```

CD plots estimate the conditional density of the `x` variable given the
levels of `y`, weighted by the marginal proportions of `y` and use
these to estimate cumulative probabilities.

* The slice at a particular `x` level visualizes the conditional
  distribution of `y` given `x` at that level.

* `geom_density` with `position = fill` is one way to create a CD
  plot.

* The `cd_plot` function from the `vcd` package produces a CD plot
  using `grid` graphics.

* The `cdplot` function from the base `graphics` package provides the
  same plots using base graphics.

<!-- Building a CD plot in steps:
```r
d <- filter(Arthritis, Treatment == "Treated", Sex == "Female")

## use y = after_stat(count) so area is number of rows
p0 <- ggplot(d, aes(x = Age, y = after_stat(count))) +
    geom_density(bw = 5, fill = "grey")
p0

## separate count-weighted densities for each group with alpha blending
p1 <- ggplot(d, aes(x = Age, y = after_stat(count), fill = Improved)) +
    geom_density(bw = 5, alpha = 0.5) +
    scale_fill_brewer(palette = "Blues") +
    guides(fill = "none")
p1

## easier to see the separate densities with faceting
p1 + facet_wrap(~ Improved)

## stack the densities
ggplot(d, aes(x = Age, y = after_stat(count), fill = Improved)) +
    geom_density(position = "stack", bw = 5) +
    scale_fill_brewer(palette = "Blues") +
    guides(fill = "none")

## rescale to height 1 with position = "fill"
p2 <- ggplot(d, aes(x = Age, y = after_stat(count), fill = Improved)) +
    geom_density(position = "fill", bw = 5) +
    scale_fill_brewer(palette = "Blues") +
    guides(fill = "none")
p2

p2 + geom_vline(xintercept = 40, lty = 2)

p2 + geom_vline(xintercept = 60, lty = 2)
```
-->

CD plots for the `Treated` group:

```{r, fig.height = 6, class.source = "fold-hide"}
filter(Arthritis, Treatment == "Treated") |>
    ggplot(aes(x = Age, fill = Improved)) +
    geom_density(position = "fill", bw = 5) +
    scale_fill_brewer(palette = "Blues") +
    facet_wrap(~ Sex, ncol = 1) +
    thm
```

CD plots for all combinations end up with one group of size one and
one of size zero, which produces a non-useful plot for one
combination:

```{r, fig.width = 8, fig.height = 6, class.source = "fold-hide"}
count(Arthritis, Treatment, Sex, Improved) |>
    complete(Treatment, Sex, Improved,
             fill = list(n = 0)) |>
    filter(n < 2)

ggplot(Arthritis,
       aes(x = Age, fill = Improved)) +
    geom_density(position = "fill", bw = 5) +
    scale_fill_brewer(palette = "Blues") +
    facet_grid(Treatment ~ Sex) +
    thm
```


## Uncertainty Representation

Categorical data are often analyzed by fitting models representing
conditional independence structures.

* Plotting residuals from these models can help assess how well they
  fit.

* `vcd::mosaic` supports using color to represent magnitude of residuals
  for comparing to a simple independence model.

For the `Arthritis` data, observed counts and expected counts under an
independence model assuming `Treatment` and `Improved` are independent
can be visualized as mosaic plots:

```{r, fig.width = 8, fig.height = 4, class.source = "fold-hide"}
## there are easier ways do do this ...
v <- count(Arthritis, Treatment, Improved)
pT <- group_by(v, Treatment) |>
    summarize(n = sum(n)) |>
    mutate(pT = n / sum(n)) |>
    select(-n)
pI <- group_by(v, Improved) |>
    summarize(n = sum(n)) |>
    mutate(pI = n / sum(n)) |>
    select(-n)
v <- left_join(v, pT, "Treatment") |>
    left_join(pI, "Improved") |>
    mutate(p = pT * pI,
           Treatment = fct_rev(Treatment))

po <- ggplot(v) +
    geom_mosaic(aes(weight = n, x = product(Improved, Treatment),
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues") +
    guides(fill = "none") +
    labs(title = "Observed Proportions") +
    thm +
    coord_flip() +
    theme(axis.text.y = element_text(angle = 90, hjust = 0))

pe <- ggplot(v) +
    geom_mosaic(aes(weight = p, x = product(Improved, Treatment),
                    fill = Improved)) +
    scale_fill_brewer(palette = "Blues") +
    guides(fill = "none") +
    labs(title = "Expected Proportions") +
    thm +
    coord_flip() +
    theme(axis.text.y = element_text(angle = 90, hjust = 0))

po + pe
```

A plot for assessing the fit of the residuals between the observed and
expected data under a model assuming independence of `Treatment` and
`Improved` produces:

```{r, class.source = "fold-hide"}
vcd::mosaic(~ Treatment + Improved,
            data = Arthritis,
            gp = vcd::shading_max)
```

Another visualization of the residuals is the _association plot_
produced by `assoc`:

```{r, class.source = "fold-hide"}
vcd::assoc(~ Treatment + Improved,
           data = Arthritis,
           gp = vcd::shading_max)
```


## References

> The vignette [_Residual-Based Shadings in
> vcd_](https://cran.r-project.org/package=vcd/vignettes/residual-shadings.pdf)
> in the `vcd` package.

> Zeileis, Achim, David Meyer, and Kurt Hornik. "Residual-based
> shadings for visualizing (conditional) independence." Journal of
> Computational and Graphical Statistics 16, no. 3 (2007): 507-525.

> The vignette [_Working with categorical data with R and the vcd and
  vcdExtra
  packages_](https://www.datavis.ca/courses/VCD/vcd-tutorial.pdf)
  in the `vcdExtra` package.

Several other experimental mosaic plot implementations are available
for `ggplot`.


## Some Other Visualizations


### Tree Maps

Tree maps show hierarchically structured (or tree-tructured) data.

* Each branch is represented by a rectangle.

* Leaf node tiles have areas proportional to the value of a variable.

* Tiles are often colored to reflect the value of another variable.

The package `treemapify` provides a `ggplot`-based implementation.

The data set `G20` includes some variables on the G-20 member countries:

```{r, class.source = "fold-hide"}
library(treemapify)
select(G20, region,
       country, gdp_mil_usd, hdi) |>
    knitr::kable(format = "html") |>
    kableExtra::kable_styling(
                    full_width = FALSE)
```

A simple tree with only one level, the individual countries:

```{r, include = FALSE}
library(nomnoml)
```
<!--

# nolint start
-->
<center>
```{nomnoml, echo = FALSE}
#padding: 25
#fontsize: 18
#linewidth: 2
#direction: down

[Root] -> [Germany]
[Root] -> [France]
[Root] -> [Japan]
[Root] -> [China]
[Root] -> [...]
```
</center>
<!--

# nolint end
-->

A corresponding tree map based on `gdp_mil_usd`:

```{r, class.source = "fold-hide"}
ggplot(G20, aes(area = gdp_mil_usd)) +
    geom_treemap() +
    geom_treemap_text(aes(label = country),
                      color = "white")
```

A tree grouping by region:

<center>
```{nomnoml, echo = FALSE}
#padding: 25
#fontsize: 18
#linewidth: 2
#direction: down

[Root] -> [Europe]
[Root] -> [Asia]
[Root] -> [...]
[Europe] -> [Germany]
[Europe] -> [France]
[Europe] -> [....]
[Asia] -> [Japan]
[Asia] -> [China]
[Asia] -> [.....]

```
</center>

A corresponding tree map:

```{r, class.source = "fold-hide"}
ggplot(G20, aes(area = gdp_mil_usd,
                subgroup = region)) +
    geom_treemap() +
    geom_treemap_text(aes(label = country),
                      color = "white") +
    geom_treemap_subgroup_border(
        color = "red") +
    geom_treemap_subgroup_text(color = "red")
```

A tree map showing GDP values for the G-20 members, grouped by region,
with fill mapped to the country's Human Development Index:

```{r, class.source = "fold-hide"}
ggplot(G20, aes(area = gdp_mil_usd,
                fill = hdi,
                subgroup = region)) +
    geom_treemap() +
    geom_treemap_text(aes(label = country),
                      color = "white") +
    geom_treemap_subgroup_border() +
    geom_treemap_subgroup_text(
        color = "lightgrey")
```

A treemap representing the distribution of eye color within hair color:

```{r, class.source = "fold-hide"}
group_by(agg, Eye, Hair) |>
    summarize(n = sum(n)) |>
    ungroup() |>
    ggplot(aes(area = n,
               subgroup = Hair)) +
    geom_treemap(aes(fill = Eye),
                 color = "white") +
    geom_treemap_subgroup_text() +
    geom_treemap_subgroup_border(
        color = "black", size = 6) +
    geom_treemap_text(aes(label = Eye),
                      color = "grey90") +
    scale_fill_manual(values = ecols) +
    guides(fill = "none")
```

A treemap representing proportions for `Improved` within `Treatment`
within `Sex` for the Arthritis data:

```{r, class.source = "fold-hide"}
count(Arth, Treatment, Improved, Sex) |>
    ggplot(aes(area = n,
               subgroup = Sex, fill = Improved,
               subgroup2 = Treatment)) +
    geom_treemap() +
    geom_treemap_subgroup_text() +
    scale_fill_brewer(palette = "Blues",
                      direction = -1) +
    geom_treemap_subgroup_border() +
    geom_treemap_subgroup2_text(place = "top",
                                size = 20)
```


### Alluvial plots

These are also known as

* _parallel sets_, or

* _Sankey diagrams_.

They can be viewed as a parallel coordinates plot for categorical data.

Several implementations are available, including:

<!-- * `alluvial` using base graphics; -->

* `geom_parallel_sets` from `ggforce`;

* `geom_sankey` from [`ggsankey`](https://github.com/davidsjoberg/ggsankey);

* `geom_alluvium` from `ggalluvial`.

<!-- Hair/Eye color using the `alluvial` package: -->

```{r, echo = FALSE, eval = FALSE}
HDF <- mutate(HairEyeColorDF, Sex = fct_rev(Sex))
library(alluvial)
pal <- RColorBrewer::brewer.pal(3, "Set1")
with(HDF,
     alluvial(Hair, Eye, Sex, freq = Freq, col = pal[as.numeric(Sex)]))
```

Hair/Eye color using the `ggforce` package:

```{r, class.source = "fold-hide"}
pal <- RColorBrewer::brewer.pal(3, "Set1")
HDF <- mutate(HairEyeColorDF,
              Sex = fct_rev(Sex))
library(ggforce)
sHDF <- gather_set_data(HDF, 3 : 1)
sHDF <- mutate(sHDF, x = fct_inorder(as.factor(x))) #**** simplify this?
ggplot(sHDF, aes(x, id = id,
                 split = y,
                 value = Freq)) +
    geom_parallel_sets(aes(fill = Sex),
                       alpha = 0.5,
                       axis.width = 0.1) +
    geom_parallel_sets_axes(
        axis.width = 0.1) +
    geom_parallel_sets_labels(
        colour = 'white') +
    scale_fill_manual(
        values = c(Male = pal[2],
                   Female = pal[1])) +
    theme_void() + guides(fill = "none")
```

<!-- Arthritis data with `alluvial`: -->

```{r, echo = FALSE, eval = FALSE}
with(count(Arth, Improved, Treatment, Sex),
     alluvial(Improved, Treatment, Sex, freq = n, col = pal[as.numeric(Sex)]))
```

Arthritis data with `ggforce`:

```{r, class.source = "fold-hide"}
sArth <- mutate(Arth,
                Improved = factor(Improved,
                                  ordered = FALSE)) |>
    count(Improved, Treatment, Sex) |>
    gather_set_data(3 : 1)
sArth <- mutate(sArth,
                x = fct_inorder(factor(x)),
                Sex = fct_rev(Sex))
ggplot(sArth, aes(x,
                  id = id,
                  split = y,
                  value = n)) +
    geom_parallel_sets(aes(fill = Sex),
                       alpha = 0.5,
                       axis.width = 0.1) +
    geom_parallel_sets_axes(axis.width = 0.1) +
    geom_parallel_sets_labels(
        colour = 'white') +
    scale_fill_manual(
        values = c(Male = pal[2],
                   Female = pal[1])) +
    theme_void() + guides(fill = "none")
```

```{r, eval = FALSE, echo = FALSE}
## refugee data (needs adjusting text sizes or labels)
## could also try chord diagram
R <- count(ref.dtl, Destination, Nationality, wt = Cases)
nats <- (count(R, Nationality, wt = n) |> slice_max(n, n = 10))$Nationality
dsts <- (count(R, Destination, wt = n) |> slice_max(n, n = 10))$Destination
mutate(R, Nationality = fct_other(Nationality, keep = nats),
       Destination = fct_other(Destination, keep = dsts)) |>
    gather_set_data(1 : 2) |>
    ggplot(aes(x, id = id, split = y, value = n)) +
    geom_parallel_sets(aes(fill = Nationality), alpha = 0.5, axis.width = 0.1) +
    geom_parallel_sets_axes(axis.width = 0.1) +
    geom_parallel_sets_labels(colour = 'white', size = 3) +
    theme_void() + guides(fill = "none")
```


### Stream Graphs

[Stream
  graphs](https://www.visualisingdata.com/2010/08/making-sense-of-streamgraphs/)
  are a generalization of stacked bar charts plotted against a numeric
  variable.

In some cases the origins of the bars are shifted to improve some
aspect of the overall visualization.

An early example is the [Baby Name
Voyager](https://www.bewitched.com/namevoyager.html).
(A more recent variant is also
[available](https://namerology.com/baby-name-grapher/).)

A [NY Times
visualization](https://archive.nytimes.com/www.nytimes.com/interactive/2008/02/23/movies/20080223_REVENUE_GRAPHIC.html?_r=0)
of movie box office results is another example. ([Blog post with a
static
version](https://flowingdata.com/2008/02/25/ebb-and-flow-of-box-office-receipts-over-past-20-years/)).

Some R implementations on GitHub:

* [`ggTimeSeries`](https://github.com/AtherEnergy/ggTimeSeries)
* [`streamgraph`](https://hrbrmstr.github.io/streamgraph/) (uses D3)
* [`ggstream`](https://github.com/davidsjoberg/ggstream).

A stream graph for movie genres (these are not mutually exclusive):

```{r, warning = FALSE, class.source = "fold-hide"}
## install with: remotes::install_github("hrbrmstr/streamgraph")
library(streamgraph)
library(tidyverse)
genres <- c("Action", "Animation", "Comedy",
            "Drama", "Documentary", "Romance")
mymovies <- select(ggplot2movies::movies,
                   year, one_of(genres))
mymovies_long <- pivot_longer(
    mymovies, -year,
    names_to = "genre",
    values_to = "value")
movie_counts <- count(mymovies_long,
                      year, genre)
streamgraph(movie_counts, "genre", "n", "year")
```

<!--
nice example with Xmen characters from TidyTuesday
https://github.com/Z3tt/TidyTuesday/blob/master/R/2020_27_ClaremontRunXMen.Rmd
-->


## Reading

Chapters [_Visualizing
  proportions_](https://clauswilke.com/dataviz/visualizing-proportions.html)
  and [_Visualizing nested
  proportions_](https://clauswilke.com/dataviz/nested-proportions.html)
  in [_Fundamentals of Data
  Visualization_](https://clauswilke.com/dataviz/).


## Interactive Tutorial

An interactive [`learnr`](https://rstudio.github.io/learnr/) tutorial
for these notes is [available](`r WLNK("tutorials/proportions.Rmd")`).

You can run the tutorial with

```{r, eval = FALSE}
STAT4580::runTutorial("proportions")
```

You can install the current version of the `STAT4580` package with

```{r, eval = FALSE}
remotes::install_gitlab("luke-tierney/STAT4580")
```

You may need to install the `remotes` package from CRAN first.


## Exercises

1. Figure A shows a bar char of the flights leaving NYC airports in
   2013 for each day of the week. Figure B shows the market share of
   five major internet browsers in 2015.

    ```{r, message = FALSE, echo = FALSE, fig.height = 4, fig.width = 8}
    library(lubridate)
    library(dplyr)
    library(nycflights13)
    library(ggplot2)
    library(patchwork)
    thm <- theme_minimal() + theme(text = element_text(size = 15))
    p1 <- mutate(flights,
                 date = make_date(year, month, day),
                 wday = wday(date, label = TRUE, abbr = TRUE)) |>
        ggplot(aes(x = wday)) +
        geom_bar(fill = "deepskyblue3") +
        scale_y_continuous(expand = c(0, 0)) +
        labs(title = "Flights from NYC in 2013",
             subtitle = "By Day of the Week",
             caption = "Figure A",
             x = NULL,
             y = "Number of Flights") +
        thm

    browsers2015 <-
        data.frame(Browser = c("Opera", "Safari", "Firefox", "Chrome", "IE"),
                   share = c(2, 22, 21, 27, 29))
    p2 <- ggplot(browsers2015, aes(x = Browser, y = share)) +
        geom_col(fill = "deepskyblue3") +
        scale_y_continuous(expand = c(0, 0)) +
        labs(title = "Browser Market Share",
             subtitle = "2015",
             caption = "Figure B",
             x = NULL,
             y = "Percent") +
        thm

    p1 | p2
    ```

    For which of these bar charts would it be better to reorder the
    categories so the bars are ordered from largest to smallest?
    
    a. Yes for Figure A. No for Figure B.
    b. No for Figure A. Yes for Figure B.
    c. Yes for both.
    d. No for both.

2.  Consider the stacked bar chart `p1` and the spine plot `p2` for
    the hair and eye color data produced by the following code:

    ```{r, eval = FALSE}
    library(dplyr)
    library(ggplot2)
    library(ggmosaic)
    ecols <- c(Brown = "brown4", Blue = "blue2",
               Hazel = "darkgoldenrod3", Green = "green4")
    HairEyeColorDF <- as.data.frame(HairEyeColor)
    p0 <- ggplot(HairEyeColorDF) +
        scale_fill_manual(values = ecols) +
        theme_minimal()

    p1 <- p0 + geom_col(aes(x = Hair, y = Freq / sum(Freq), fill = Eye))

    p2 <- p0 + geom_mosaic(aes(x = Hair, fill = Eye, weight = Freq))
    ```

    Use the two plots to answer: Which hair color has the highest
    proportion of individuals with green eyes?

    a. Black
    b. Brown
    c. Red
    d. Blond

    Which plot makes it easiest to answer this question?

3.  Use the plots of the previous question to answer: The proportion
    of individuals with red hair is closest to:

    a. 5%
    b. 8%
    c. 12%
    d. 20%

    Which plot makes it easiest to answer this question?

<!--
pareto chart

  - pie charts
  - bar charts
  - stacked bar charts
  - grouped bar charts
  - population pyramids
  - waffle charts, square pie charts

  - joint and conditional distributions
  - spine plots
  - spinograms?
  - double decker plots
  - cd_plots
  - mosaic plots
  - tree maps
  - sankey diagrams, alluvial charts, parallel sets
  - chord diagrams
  - stream graphs
-->
