12-Graphs.Rmd

# Graphical summaries {#Graphs}




<!-- Introductions; easier to separate by format -->

```{r, child = if (knitr::is_html_output()) {'./introductions/12-Graphs-HTML.Rmd'} else {'./introductions/12-Graphs-LaTeX.Rmd'}}
```


## Introduction {#GraphsIntro}

To answer a RQ, a study is designed to collect data, because *data are the means by which the research question is answered*.
But before analysis begins, understanding, describing and summarising the collected data are important.
This chapter discusses the use of *graphs* to summarise data.
Graphs are produced to aid understanding of the data, so the message in the graphs should always be explained, especially regarding what it reveals about answering the RQ.


::: {.importantBox .important data-latex="{iconmonstr-warning-8-240.png}"}
The purpose of a graph is to display the information in the clearest, simplest possible way, to help the reader understand the message(s) in the data.
:::

In this chapter, graphs are discussed for

:::::: {.cols data-latex=""}

::: {.col data-latex="{0.48\textwidth}"}

* One quantitative variable (Sect.\ \@ref(GraphsOneQuant));
* One qualitative variable (Sect.\ \@ref(GraphsOneQual));
* One qualitative and quantitative variable (Sect.\ \@ref(GraphOneQualOneQuant));

:::

::: {.col data-latex="{0.01\textwidth}"}
\ 
<!-- an empty Div (with a white space), serving as
a column separator -->
:::

::: {.col data-latex="{0.51\textwidth}"}

* Two quantitative variables (Sect.\ \@ref(Scatterplots));
* Two qualitative variables (Sect.\ \@ref(TwoQualVars));
* Other types of graphs (Sect.\ \@ref(GraphsOther)).

:::
::::::


## One quantitative variable {#GraphsOneQuant}

For quantitative data, a graph shows the *distribution* of the data: the values taken by the variable, and how often they appear.


::: {.definition #Distribution name="Distribution"}
The *distribution* of a quantitative variable describes what values are present in the data, and how often those values appear.
:::


The graphs discussed in this section are appropriate for *continuous* quantitative data, but may also be useful for discrete quantitative data if many possible values are present.
Sometimes, *discrete* quantitative data with very few values are best graphed using graphs for qualitative data (Sect.\ \@ref(GraphsOneQual)).

Three different types of graphs can be used to show the distribution of one quantitative variable:

* *Stemplots* (or *stemplots*) (Sect.\ \@ref(StemAndLeafPlots)):
  Best for small amounts of data; only sometimes useful.
* *Dot chart* (Sect.\ \@ref(DotChartsOneVar)):
  Best for small amounts of data; good for moderate amounts of data.
* *Histogram* (Sect.\ \@ref(Histograms)):
  Best for moderate to large amounts of data.



### Stemplots {#StemAndLeafPlots}

*Stemplots* (or *stem-and-leaf plots*) are best described and explained using an example.
Consider the data in
`r if( knitr::is_html_output() ) {
  'Fig.\\ \\@ref(fig:BabyBoomDataHTML),'
} else
{
  'Table\\ \\@ref(tab:BabyBoomDataLATEX) (which shows the first 10 of the total 44 observations).'
}`
The data give the weights (in kg) of babies born in a Brisbane hospital on one day [@mypapers:Dunn:dataset:1999; @data:Steele:BabyBoom].
The dataset also includes the gender of each baby, and the number of minutes after midnight of the birth.
The data are given in the order in which the births occurred.



<div style="float:right; width: 222x; border: 1px; padding:10px">
<img src="Illustrations/pexels-nappy-3569597.jpg" width="200px"/>
</div>





```{r}
data(BabyBoom)
BB <- BabyBoom

BB2 <- data.frame(Gender = factor(BB$Gender),
                  "Weight (in kg)" = round( BB$Weight/1000, 1),
                   Mins.Since.Midnight = BB$Mins.Since.Midnight)
```

```{r BabyBoomDataLATEX}
if( knitr::is_latex_output() ) {
  T1 <- kable( BB2[1:5, ],
               format = "latex",
               longtable = FALSE,
               row.names = FALSE,
               digits = c(0, 1, 0),
               escape = FALSE,
               col.names = c("Gender", 
                             "Weight (kg)", 
                             "Time of birth"),
               booktabs = TRUE) %>%
    row_spec(0, bold = TRUE)
#    column_spec(1, width = "10mm") %>%
#    column_spec(2, width = "10mm") %>%
#    column_spec(3, width = "10mm")
  T2 <- kable( BB2[6:10, ],
               format = "latex",
               longtable = FALSE,
               row.names = FALSE,
               digits = c(0, 1, 0),
               col.names = c("Gender", 
                             "Weight (kg)", 
                             "Time of birth"),
               booktabs = TRUE) %>%
    row_spec(0, bold = TRUE)
  
  out <- knitr::kables(list(T1, T2),
                       format = "latex",
                       label = "BabyBoomDataLATEX",
                       caption = "The first ten observations (of 44) of the baby-births data. The 'time of birth is the number of minutes since midnight.") %>% 
    kable_styling(font_size = 10)
  out2 <- prepareSideBySideTable(out, 
                                 numberOfTables = 2,
                                 gap = "\\quad") 
  out2
  
}
```
```{r BabyBoomDataHTML, fig.cap="The baby-births data"}
if( knitr::is_html_output() ) {
      DT::datatable( BB2,
                     fillContainer = FALSE, # Make more room, so we don't just have ten values
                     options = list(searching = FALSE), # Remove searching: See: https://stackoverflow.com/questions/35624413/remove-search-option-but-leave-search-columns-option
                     colnames = c("Gender", 
	                                "Weight (in kg)", 
		                              "Minutes since midnight"),
                     caption = "The baby-births data"
        )
}
```


In a stemplot, part of each number is placed to the left of a vertical line (the *stem*), and the rest of each number to the right (the *leaf*).
The weights (quantitative) 
`r if( knitr::is_html_output() ) {
  'in Fig.\\ \\@ref(fig:BabyBoomDataHTML)'
} else
{
  'in Table\\ \\@ref(tab:BabyBoomDataLATEX)'
}`
are given to one decimal place of a kilogram, so the whole number of kilograms is placed to the left (as a *stem*), and the first decimal place is placed on the right (as a *leaf*).
`r if (knitr::is_latex_output()) {
   'Figure\\ \\@ref(fig:BBStem1) shows the stemplot starting to be built, and Fig.\\ \\@ref(fig:BBStem44) shows the final stemplot. (The online version has an animation.)'
} else {
   'The animation below shows how the stemplot is constructed.'
}`
The first weight, of 1.7kg, is entered with the 1 to the right of the stem, and the 7 to the left: `1 | 7`. 
Similarly, 2.1kg is entered as `2 | 1` and 2.2kg is entered as `2 | 2`, sharing the same stem as for 2.1kg.
The plots shows that most birthweights are 3-point-something kilograms.


```{r animation.hook="gifski", interval=0.85, dev=if (is_latex_output()){"pdf"}else{"png"}}
# BABYBOOM
if (knitr::is_html_output()){
  Wts <- sort(BB$Weight)
    
  DP2 <- function(x){format(round(x, 1), nsmall = 1) }
  
  par(mfrow = c(2, 1), 
      mar = c(2, 3, 4, 1) + 0.1 )
  for (i in (1:length(Wts))){
    # Plot the data
    datagrid.x <- 11
    datagrid.y <- 4
    
    plot( expand.grid(1:datagrid.x, 
                      1:datagrid.y),
          type = "n",
          ylim = c(0.7, 4.3),
          main = "The birthweights (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    cex.data <- 0.9
    text(  1:datagrid.x, 
           rep(4, 11), 
	   DP2(Wts[1:11]),  
	   cex = cex.data )
    text(  1:datagrid.x, 
           rep(3, 11), 
	   DP2(Wts[12:22]), 
	   cex = cex.data )
    text(  1:datagrid.x, 
           rep(2, 11), 
	   DP2(Wts[23:33]), 
	   cex = cex.data )
    text(  1:datagrid.x, 
           rep(1, 11), 
	   DP2(Wts[34:44]), 
	   cex = cex.data )
    
    # Now highlight the observation being used
    if (i <= 11) {
      points(  i, 4, 
               pch = 1, 
	       col = plotDark, 
	       lwd = 3, 
	       cex = 4 )
    }
    if ( (i > 11) & (i <= 22)){
      points(  i - 11, 3, 
               pch = 1, 
	       col = plotDark, 
	       lwd = 3, 
	       cex = 4 )
    }
    if ( (i > 22) & (i <= 233)){
      points(  i - 22, 2, 
               pch = 1, 
	       col = plotDark, 
	       lwd = 3, 
	       cex = 4 )
    }
    if ( (i > 33) ){
      points(  i - 33, 1, 
               pch = 1, 
	       col = plotDark, 
	       lwd = 3, 
	       cex = 4 )
    }
    
    # Produce the stemplot
    stemgrid.x <- 34
    stemgrid.y <- 4
    stemcounts <- array(1, dim = stemgrid.y)
    
    plot( expand.grid(1:stemgrid.x, 1:stemgrid.y), 
          type = "n",
          xlab = "",
          ylab = "",
          main = "Stemplot of birthweights\n(2 | 6 means 2.6 kg)",
          xlim = c(0, 34),
          ylim = c(0.8, 4.2),
          axes = FALSE)
    
    if (Wts[i] < 2 ) {
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 4, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("4", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1)
    }
    if ( (Wts[i] >= 2) & (Wts[i] < 3) ){
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 4)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("4", 
             side = 2, 
	     line = 1, 
	     las = 1, 
	     at = 1)
    }
    if ( (Wts[i] >= 3) & (Wts[i] < 4) ){
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 4)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
      mtext("4", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1)
    }
    if ( Wts[i] >= 4 ){
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 4)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("4", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
    }
    
    # Divide the stems and leaves
    abline(v = 0, 
           lwd = 2)
    
    # Plot leaves
    for (j in (1:i) ){
      if ( Wts[j] < 2) {
        if (j != i) points( stemcounts[1], 4, 
	                    pch = substr( sprintf("%1.1f", Wts[j]), 3, 3))
        if (j == i ) text( stemcounts[1], 4, 
	                   substr( sprintf("%1.1f", Wts[j]), 3, 3), 
			   col = plotDark, 
			   cex = 1.5, 
			   font = 2)
        stemcounts[1] <- stemcounts[1] + 1
      }
      
      if ( (Wts[j] >= 2) & (Wts[j] < 3) ) {
        if (j != i ) points( stemcounts[2], 3, 
	                     pch = substr( sprintf("%1.1f", Wts[j]), 3, 3) )
        if (j == i ) text( stemcounts[2], 3, 
	                   substr( sprintf("%1.1f", Wts[j]), 3, 3), 
			   col = plotDark, 
			   cex = 1.5, 
			   font = 2)
        stemcounts[2] <- stemcounts[2] + 1
      }
      
      if ( (Wts[j] >= 3) & (Wts[j] < 4) ) {
        if (j != i ) points( stemcounts[3], 2, 
	                     pch = substr( sprintf("%1.1f", Wts[j]), 3, 3))
        if (j == i ) text( stemcounts[3], 2, 
	                   substr( sprintf("%1.1f", Wts[j]), 3, 3), 
			   col = plotDark, 
			   cex = 1.5, 
			   font = 2)
        stemcounts[3] <- stemcounts[3] + 1
      }
      if ( Wts[j] >= 4 ) {
        if (j != i ) points( stemcounts[4], 1, 
	                     pch = substr( sprintf("%1.1f", Wts[j]), 3, 3))
        if (j == i ) text( stemcounts[4], 1, 
	                   substr( sprintf("%1.1f", Wts[j]), 3, 3), 
			   col = plotDark, 
			   cex = 1.5,  
			   font = 2)
        stemcounts[4] <- stemcounts[4] + 1
      }
    }
  }
}
```


```{r BBStem1, fig.align="center", fig.width=5, fig.cap="Starting to make the stemplot for the baby-weight data: the first 4 observations added. The data are on the left; the stemplot during construction on the right.", fig.width=10, fig.height=3, out.width='95%'}
if (knitr::is_latex_output()){
  Wts <- sort(BB$Weight)
  
  DP2 <- function(x){ format(round(x, 1), 
                             nsmall = 1) }
  
  #  par(width=5, height=5, mar=c(0, 0, 2, 0)+0.1, mfrow=c(2, 1))
  par(mfrow = c(1, 2), 
      mar = c(3.85, 2, 3.5, 1) + 0.1)
  i <- 4
    # Plot the data
    datagrid.x <- 11
    datagrid.y <- 4
    
    plot( expand.grid(1:datagrid.x, 1:datagrid.y),
          type = "n",
          ylim = c(0.7, 4.3),
          xlim = c(1, 13),
          main = "The birthweights (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    cex.data <- 0.9
    text(  1:datagrid.x, rep(4, 11), 
           DP2(Wts[1:11]),  
	   cex = cex.data )
    text(  1:datagrid.x, rep(3, 11), 
           DP2(Wts[12:22]), 
	   cex = cex.data )
    text(  1:datagrid.x, rep(2, 11), 
           DP2(Wts[23:33]), 
	   cex = cex.data )
    text(  1:datagrid.x, rep(1, 11), 
           DP2(Wts[34:44]), 
	   cex = cex.data )
    
    # Now highlight the observation being used
    if (i <= 11) {
      points(  i, 4, 
               pch = 1, 
               col = plotDark, 
               lwd = 3, 
               cex = 4 )
    }

   
    # Produce the stemplot
    stemgrid.x <- 34
    stemgrid.y <- 4
    stemcounts <- array(1, dim = stemgrid.y)
    
    plot( expand.grid(1:stemgrid.x, 1:stemgrid.y), 
          type = "n",
          xlab = "2 | 6 means 2.6 kg",
          ylab = "",
          main = "Stemplot of birthweights",
          xlim = c(0, 34), 
          ylim = c(0.8, 4.2),
          axes = FALSE)
    
    if (Wts[i] < 2 ) {
      mtext("1", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 4, 
            col = plotDark, 
            font = 2, 
            cex = 1.5)
      mtext("2", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 3)
      mtext("3", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 2)
      mtext("4", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 1)
    }
    if ( (Wts[i] >= 2) & (Wts[i] < 3) ){
      mtext("1", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 4)
      mtext("2", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 3, 
            col = plotDark, 
            font = 2, 
            cex = 1.5)
      mtext("3", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 2)
      mtext("4", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 1)
    }
    if ( (Wts[i] >= 3) & (Wts[i] < 4) ){
      mtext("1", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 4)
      mtext("2", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 3)
      mtext("3", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 2, 
            col = plotDark, 
            font = 2, 
            cex = 1.5)
      mtext("4", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 1)
    }
    if ( Wts[i] >= 4 ){
      mtext("1", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 4)
      mtext("2", 
            side = 2,
            line = 1, 
            las = 1, 
            at = 3)
      mtext("3", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 2)
      mtext("4", 
            side = 2, 
            line = 1, 
            las = 1, 
            at = 1, 
            col = plotDark, 
            font = 2, 
            cex = 1.5)
    }
    
    # Divide the stems and leaves
    abline(v = 0,  
           lwd = 2)
    
    # Plot leaves
    stemcounts <- rep(1, 4)
    
    for (j in (1:4) ){
      
      if ( Wts[j] < 2) {
        if (j != i) points( stemcounts[1], 4, 
                            pch = substr( as.character(Wts[j] - 1), 3, 3) )
        if (j == i ) text( stemcounts[1], 4, 
                           substr( as.character(Wts[j] - 1), 3, 3),
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[1] <- stemcounts[1] + 1
      }
      
      if ( (Wts[j] >= 2) & (Wts[j] < 3) ) {
        if (j != i ) points( stemcounts[2], 3, 
                             pch = substr( as.character(Wts[j] - 2), 3, 3) )
        if (j == i ) text( stemcounts[2], 3, 
                           substr( as.character(Wts[j] - 2), 3, 3),
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[2] <- stemcounts[2] + 1
      }
      
      if ( (Wts[j] >= 3) & (Wts[j] < 4) ) {
        if (j != i ) points( stemcounts[3], 2, 
                             pch = substr( sprintf("%1.1f", Wts[j]), 3, 3) )
        if (j == i ) text( stemcounts[3], 2, 
                           substr( sprintf("%1.1f", Wts[j]), 3, 3),
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[3] <- stemcounts[3] + 1
      }
      if ( Wts[j] >= 4 ) {
        if (j != i ) points( stemcounts[4], 1, 
                             pch = substr( as.character(Wts[j] - 2), 3, 3) )
        if (j == i ) text( stemcounts[4], 1, 
                           substr( as.character(Wts[j] - 4), 3, 3),
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[4] <- stemcounts[4] + 1
      }
    }
}
```


```{r BBStem44,  fig.align="center", fig.cap="The final stemplot for the baby-weight data. The data are on the left; the final stemplot on the right.", fig.width=10, fig.height=3, out.width='100%' }
if (knitr::is_latex_output()){
  Wts <- sort(BB$Weight)
  
  DP2 <- function(x){format(round(x, 1), nsmall = 1) }
  
  par(mfrow = c(1, 2), 
      mar = c(3.85, 2, 3.5, 1) + 0.1)
  for (i in (length(Wts):length(Wts))){
    # Plot the data
    datagrid.x <- 11
    datagrid.y <- 4
    
    plot( expand.grid(1:datagrid.x, 1:datagrid.y),
          type = "n",
          ylim = c(0.7, 4.3),
          xlim = c(1, 13),
          main = "The birthweights (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    cex.data <- 0.9
    text(  1:datagrid.x, 
           rep(4, 11), 
           DP2(Wts[1:11]),  
           cex = cex.data )
    text(  1:datagrid.x, 
           rep(3, 11), 
           DP2(Wts[12:22]), 
           cex = cex.data )
    text(  1:datagrid.x, 
           rep(2, 11), 
           DP2(Wts[23:33]), 
           cex = cex.data )
    text(  1:datagrid.x, 
           rep(1, 11), 
           DP2(Wts[34:44]), 
           cex = cex.data )
    
    # Produce the stemplot
    stemgrid.x <- 34
    stemgrid.y <- 4
    stemcounts <- array(1, dim = stemgrid.y)
    
    plot( expand.grid(1:stemgrid.x, 1:stemgrid.y), 
          type = "n",
          xlab = "2 | 6 means 2.6 kg",
          ylab = "",
          main = "Stemplot of birthweights",
          xlim = c(0, 34),
          ylim = c(0.8, 4.2),
          axes = FALSE)
    
    mtext("1", 
          side = 2, 
          line = 1, 
          las = 1, 
          at = 4)
    mtext("2", 
          side = 2, 
          line = 1, 
          las = 1, 
          at = 3)
    mtext("3", 
          side = 2, 
          line = 1, 
          las = 1, 
          at = 2)
    mtext("4", 
          side = 2, 
          line = 1, 
          las = 1, 
          at = 1)
    
    # Divide the stems and leaves
    abline(v = 0, 
           lwd = 2)
    
    # Plot leaves
    for (j in (1:44) ){
      if ( Wts[j] < 2) {
        if (j != i) points( stemcounts[1], 4, 
                            pch = substr( as.character(Wts[j] - 1), 3, 3) )
        stemcounts[1] <- stemcounts[1] + 1
      } 
      
      if ( (Wts[j] >= 2) & (Wts[j] < 3) ) {
        if (j != i ) points( stemcounts[2], 3, 
                            pch = substr( as.character(Wts[j] - 2), 3, 3) )
        stemcounts[2] <- stemcounts[2] + 1
      }
      
      if ( (Wts[j] >= 3) & (Wts[j] < 4) ) { 
        if (j != i ) points( stemcounts[3], 2, 
                            pch = substr( sprintf("%1.1f", Wts[j]), 3, 3) )
        stemcounts[3] <- stemcounts[3] + 1
      }
      if ( Wts[j] >= 4 ) {
        points( stemcounts[4], 1, 
                            pch = substr( as.character(Wts[j] - 4), 3, 3) )
        stemcounts[4] <- stemcounts[4] + 1
      }
    }
  }
}
```

For stemplots:

* Place the larger unit (e.g., kilograms) on the left (stems).
* Place the next smallest unit (e.g., first decimal place of a kilogram) on the right (leaves).
* Some data do not work well with stemplots.
* Sometimes, data may need suitable rounding before creating the stemplot (the baby weights were originally given to three decimal places).
* The numbers in each row should be evenly spaced, so that the numbers in the columns are under each other.
  This means the length of each stem is proportional to the number of observations.
* Within each stem, the observations are *ordered* so patterns can be seen.
* Add an explanation for reading the stemplot.
  For example, the stemplot for the baby-birth data says '2 | 6 means 2.6kg'.
  For instance, '2 | 6' could mean 26kg, or 0.26kg, or 2lb 6oz.


::: {.example #StemLeafPlots name="Stemplots"}
A study of krill [@data:Greenacre2016:reporting] produced 15 measurements of the number of eggs.
The stemplot 
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:KrillStem44); the online version has an animation)'
} else {
   'in the animation below shows the stemplot being constructed, and'
}`
shows that the number of eggs is usually under 10, but occasionally a large number of eggs are seen.
For the purposes of the stemplot, single-digit observations are treated as `00`, `01`, `02` and so on, using `0` as the stem. 
:::


```{r animation.hook="gifski", interval=0.85, dev=if (is_latex_output()){"pdf"}else{"png"}}
if (knitr::is_html_output()){
#KRILL

  Eggs <- c(1, 0, 1, 2, 2, 8, 16, 31, 0, 0, 20, 0, 26, 1, 3)
  Eggs <- sort(Eggs)
  
  #  par(width=5, height=5, mar=c(0, 0, 2, 0)+0.1, mfrow=c(2, 1))
  par(mfrow = c(2, 1),
      mar = c(3.85, 3, 4, 1) + 0.1)
  for (i in (1:length(Eggs))){
    # Plot the data
    datagrid.x <- 5
    datagrid.y <- 3
    
    plot( expand.grid(1:datagrid.x, 
                      1:datagrid.y),
          type = "n",
          ylim = c(0, 3.3),
          main = "The egg counts (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    text(  1:datagrid.x, 
           rep(3, 5), 
	   as.character(Eggs[1:5]) )
    text(  1:datagrid.x, 
           rep(2, 5), 
	   as.character(Eggs[6:10]) )
    text(  1:datagrid.x, 
           rep(1, 5), 
	   as.character(Eggs[11:15]) )
    
    # Now highlight the observation being used
    if (i <= 5) {
      points(  i, 3, 
               pch = 1, 
	       col = plotDark, 
	       lwd = 3, 
	       cex = 4 )
    }
    if ( (i > 5) & (i <= 10)){
      points(  i - 5, 2, 
               pch = 1, 
               col = plotDark, 
               lwd = 3, 
               cex = 4 )
    }
    if ( (i > 10) ){
      points(  i - 10, 1, 
               pch = 1, 
               col = plotDark, 
               lwd = 3, 
               cex = 4 )
    }
    
    # Produce the stemplot
    stemgrid.x <- 12
    stemgrid.y <- 4
    stemcounts <- array(1, dim = stemgrid.y)
    
    plot( expand.grid(1:stemgrid.x, 
                      0:stemgrid.y), 
          type = "n",
          xlab = "1|6 means 16 eggs",
          ylab = "",
          main = "Stemplot of egg counts of krill",
          xlim = c(0, 11),
          ylim = c(-0.2, 3.2),
          axes = FALSE)
    
    if (Eggs[i] < 10 ) {
      mtext("0", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 0)
    }
    if ( (Eggs[i] >= 10) & ( Eggs[i] < 20) ){
      mtext("0", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 0)
    }
    if ( (Eggs[i] >= 20) & (Eggs[i] < 30) ){
      mtext("0", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
      mtext("3",  
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 0)
    }
    if ( Eggs[i] >= 30 ){
      mtext("0", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 0, 
	    col = plotDark, 
	    font = 2, 
	    cex = 1.5)
    }
    # Divide the stems and leaves
    abline(v = 0, 
           lwd = 2)
    
    # Plot leaves
    for (j in (1:i) ){
      if ( Eggs[j] < 10) {
        if ( j != i) points( stemcounts[1], 
                             3, 
                             pch = as.character(Eggs[j]))
        if (j == i ) text( stemcounts[1], 3, 
                           as.character(Eggs[j]), 
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[1] <- stemcounts[1] + 1
      }
      
      if ( (Eggs[j] >= 10) & (Eggs[j] < 20) ) {
        if (j != i ) points( stemcounts[2], 
                             2, 
                             pch = as.character(Eggs[j] - 10) )
        if (j == i ) text( stemcounts[2], 
                           2, 
                           as.character(Eggs[j] - 10), 
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[2] <- stemcounts[2] + 1
      }
      
      if ( (Eggs[j] >= 20) & (Eggs[j] < 30) ) {
        if (j != i ) points( stemcounts[3], 
                             1, 
                             pch = as.character(Eggs[j] - 20))
        if (j == i ) text( stemcounts[3], 
                           1, 
                           as.character(Eggs[j] - 20), 
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[3] <- stemcounts[3] + 1
      }
      if ( Eggs[j] >= 30 ) {
        if (j != i ) points( stemcounts[4], 
                             0, 
                             pch = as.character(Eggs[j] - 30))
        if (j == i ) text( stemcounts[4], 
                           0, 
                           as.character(Eggs[j] - 30), 
                           col = plotDark, 
                           cex = 1.5, 
                           font = 2)
        stemcounts[4] <- stemcounts[4] + 1
      }
    }
  }
}
```


```{r KrillStem44, fig.align="center", fig.width=5, fig.cap="The krill data. The data are on the left, and the stemplot on the right.", fig.width=10, fig.height=2.75, out.width='100%' }
if (knitr::is_latex_output()){
#KRILL

  Eggs <- c(1, 0, 1, 2, 2, 8, 16, 31, 0, 0, 20, 0, 26, 1, 3)
  Eggs <-sort(Eggs)
  
  #  par(width=5, height=5, mar=c(0, 0, 2, 0)+0.1, mfrow=c(2, 1))
  par(mfrow = c(1, 2),
      mar = c(3.85, 3, 4, 1) + 0.1)
  for (i in (length(Eggs):length(Eggs))){
    # Plot the data
    datagrid.x <- 5
    datagrid.y <- 3
    
    plot( expand.grid(1:datagrid.x, 
                      1:datagrid.y),
          type = "n",
          ylim = c(0, 3.3),
          xlim = c(1, 7),
          main = "The egg counts (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    text(  1:datagrid.x, 
           rep(3, 5), 
	   as.character(Eggs[1:5]) )
    text(  1:datagrid.x, 
           rep(2, 5), 
	   as.character(Eggs[6:10]) )
    text(  1:datagrid.x, 
           rep(1, 5), 
	   as.character(Eggs[11:15]) )
    
    # Now highlight the observation being used

    # Produce the stemplot
    stemgrid.x <- 12
    stemgrid.y <- 4
    stemcounts <- array(1, dim = stemgrid.y)
    
    plot( expand.grid(1:stemgrid.x, 0:stemgrid.y), 
          type = "n",
          xlab = "1|6 means 16 eggs",
          ylab = "",
          main = "Stemplot of egg counts of krill",
          xlim = c(0, 11),
          ylim = c(-0.2, 3.2),
          axes = FALSE)
    
      mtext("0", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 3)
      mtext("1", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 2)
      mtext("2", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 1)
      mtext("3", 
            side = 2, 
	    line = 1, 
	    las = 1, 
	    at = 0)
    }

    # Divide the stems and leaves
    abline(v = 0, 
           lwd = 2)
    
    # Plot leaves
    for (j in 1:length(Eggs) ){
      if ( Eggs[j] < 10 ) {
        points( stemcounts[1], 3, 
	        pch = as.character(Eggs[j]))
        stemcounts[1] <- stemcounts[1] + 1
      }
      
      if ( (Eggs[j] >= 10) & (Eggs[j] < 20) ) {
        points( stemcounts[2], 2, 
	        pch = as.character(Eggs[j] - 10) )
        stemcounts[2] <- stemcounts[2] + 1
      }
      
      if ( (Eggs[j] >= 20) & (Eggs[j] < 30) ) {
        points( stemcounts[3], 1, 
	        pch = as.character(Eggs[j] - 20))
        stemcounts[3] <- stemcounts[3] + 1
      }
      if ( Eggs[j] >= 30 ) {
        points( stemcounts[4], 0, 
	        pch = as.character(Eggs[j] - 30))
        stemcounts[4] <- stemcounts[4] + 1
      }
    }
}
```


`r if (knitr::is_html_output()){
  'The following short video may help explain some of these concepts:'
}`

<div style="text-align:center;">
<iframe width="560" height="315" src="https://www.youtube.com/embed/KfkvedIBPjE" frameborder="0" allow="accelerometer; encrypted-media; gyroscope; picture-in-picture"></iframe>
</div>



### Dot charts (quantitative data) {#DotChartsOneVar}

Dot charts show the original data on a single axis, with each observation represented by a dot.

   
<div style="float:right; width: 222x; border: 1px; padding:10px">
<img src="Illustrations/pexels-dzenina-lukac-1583884.jpg" width="200px"/>
</div>


::: {.example #DotChartsQuant name="Dot charts"}
A study examined the serving size of fries at McDonald's [@data:Wetzel2005:McDonalds].
The dot chart is shown in 
`r if (knitr::is_latex_output()) {
   'Fig.\\ \\@ref(fig:FriesBabiesPDF) (left panel)'
} else {
   'Fig.\\ \\@ref(fig:FriesHTML) (top panel)'
}`
(based on @data:Wetzel2005:McDonalds).
The mass of fries is almost always *under* the target, and often substantially so.

Notice that observations are *stacked* on top of each other when  multiple observations are the same; for example, two orders of fries were 133g.
*Jittering* is also an option (Example\ \@ref(exm:DotsChartsQuant2)).
:::


```{r FriesBabiesPDF, fig.cap="Large orders of French fries (left panel), where similar observations are stacked. A dot chart of the baby-weight data (right panel), where similar observations are jittered.", fig.align="center", fig.width=9, fig.height=3, out.width='95%'}
if (knitr::is_latex_output()){

  par( mfrow = c(1, 2) )

  Fries <- c(117, 126, 128, 132, 133, 133, 134, 137, 138, 139, 139, 140, 141, 142, 142.5, 143, 143.5,  
             145, 146, 146, 151, 152, 152, 154, 154, 154.5, 154.5, 155, 156, 156.5, 157, 176)
  stripchart(Fries, 
             method = "stack",
             main = "Dotplot for mass of large\norder of fries",
             xlab = "Serving weight (in g)",
             axes = FALSE,
             pch = 20,
  	   col = plotDark,
             ylim = c(-1, 10),
             xlim = c(112, 180) )
  axis(side = 1, 
       at = seq(115, 175, 
                by = 10))
  arrows(x0 = 145, 
         y0 = 5.5, 
         x1 = 171, 
         y1 = 2, 
         length = 0.10,
         angle = 15)
  lines( x = c(171, 171),
         y = c(-2, 4),
         col = "grey")
  text(145, 6.25, 
       "McDonald's target of\n171 grams per serving", 
       cex = 0.9,
       pos = 3)

########
  
  stripchart(BB$Weight, 
           method = "jitter",
           pch = 19,
	   col = plotDark,
	   xlim = c(1.5, 4.5),
           main = "Weight of babies",
           axes = FALSE,
           xlab = "Baby weight (in kg)")
axis(side = 1,
     at = seq(1.5, 4.5, by = 0.5))
}
  
```


```{r FriesHTML, fig.height=5.0, fig.width=8, fig.cap="Large orders of French fries: Mass measurements (top panel) and percentage of target mass (below panel)", fig.align="center", out.width='95%'}
if (knitr::is_html_output()){

  par( mfrow = c(2, 1) )
  
  Fries <- c(117, 126, 128, 132, 133, 133, 134, 137, 138, 139, 139, 140, 141, 142, 142.5, 143, 143.5,  
             145, 146, 146, 151, 152, 152, 154, 154, 154.5, 154.5, 155, 156, 156.5, 157, 176)
  stripchart(Fries, 
             method = "stack",
             main = "Dotplot for mass of large\norder of fries",
             xlab = "Serving weight (in g)",
             axes = FALSE,
             pch = 20,
  	   col = plot.colour,
             ylim = c(-1, 10),
             xlim = c(112, 180) )
  axis(side = 1, 
       at = seq(115, 175, 
                by = 10))
  arrows(x0 = 145, 
         y0 = 5.5, 
         x1 = 171, 
         y1 = 2, 
         length = 0.10,
         angle = 15)
  lines( x = c(171, 171),
         y = c(-2, 4),
         col = "grey")
  text(145, 6.5, 
       "McDonald's target of\n171 grams per serving", 
       cex = 0.8,
       pos = 3)
  
  
  
  stripchart( Fries/171 * 100, 
             method = "stack",
             main = "Dotplot for percentage of target\nmass of large order of Fries",
             axes = FALSE,
             xlab = "Percentage of target weight (171g)",
             pch = 20,
  	         col = plot.colour,
             ylim = c(-1, 10),
             xlim = c(65, 105) )
  axis(side = 1, 
       at = seq(65, 105, by = 5))
  arrows(x0 = 85, 
         y0 = 5.5, 
         x1 = 100, 
         y1 = 2, 
         length = 0.10,
         angle = 15)
  lines( x = c(100, 100),
         y = c(-2, 4),
         col = "grey")
  text(85, 6.5, 
       "McDonald's target of\n171 grams per serving", 
       cex = 0.8,
       pos = 3)
}
```


::: {.example #DotsChartsQuant2 name="Dot charts"}
Consider again the weights (in kg) of babies born in a Brisbane hospital
`r if( knitr::is_html_output() ) {
  '(Fig.\\ \\@ref(fig:BabyBoomDataHTML)).'
} else
{
  '(Table\\ \\@ref(tab:BabyBoomDataLATEX)).'
}`
A dot chart (Fig.\ \@ref(fig:BabyBoomDot)) 
`r if( knitr::is_html_output() ) {
  '(Fig.\\ \\@ref(fig:BabyBoomDotTML)).'
} else
{
  '(Fig\\ \\@ref(tab:FriesBabiesPDF), right panel).'
}`
shows that most babies were born between 3 and 4kg.
Notice that the observations have been *jittered* (i.e., placed with some added randomness in the vertical direction) to avoid overplotting.
:::


```{r BabyBoomDotHTML, fig.cap="A dot chart of the baby-weight data", fig.align="center", fig.height=2.5, fig.width=5.5}
if (knitr::is_html_output()){
  stripchart(BB$Weight, 
             method = "jitter",
             pch = 19,
             col = plot.colour,
             xlim = c(1.5, 4.5),
             main = "Weight of babies",
             axes = FALSE,
             xlab = "Baby weight (in kg)")
  axis(side = 1,
       at = seq(1.5, 4.5, by = 0.5))
}
```


### Histograms {#Histograms}

Histograms are a series of boxes, where the width of the box represents a range of *values* of the variable being graphed, and the height of the box represents the *number* (or *percentage*) of observations within that range of values^[Technically, the *area* of the box is proportional to the number of observations. 
Since we only consider histograms where the boxes are all the same width, this is equivalent.].
The vertical axis can be counts (labelled as 'Counts', 'Frequency', or similar) or percentages.


<div style="float:right; width: 222x; border: 1px; padding:10px">
<img src="Illustrations/pexels-nappy-3569597.jpg" width="200px"/>
</div>


::: {.example #Histograms name="Histograms"}
Consider again the weights (in kg) of babies born in a Brisbane hospital in one day
`r if( knitr::is_html_output() ) {
  '(Fig.\ \\@ref(fig:BabyBoomDataHTML)).'
} else {
  '(Table\ \\@ref(tab:BabyBoomDataLATEX)).'
}`
A histogram (below) can be constructed for these data.
When an observation occurs on a boundary between the boxes, software usually (but not universally) places it in the *higher* box (so 2.5kg would be counted in the '2.5 to 3.0kg' box, not the '2.0 to 2.5kg' box).
That is, the boxes *include* the lower limit, but *exclude* the upper limit.
The histogram shows, for example, that 17 babies weighed 3.0kg or more, but under 3.5kg.

`r if (knitr::is_latex_output()) {
   'Figure\\ \\@ref(fig:BBHist1) shows the histogram starting to be built, and Fig.\\ \\@ref(fig:BBHist44) shows the final histogram. (The online version has an animation.)'
} else {
   'The animation below shows how the histogram is constructed.'
}`
:::


```{r animation.hook="gifski", interval=0.55, fig.height=7}
### BUILD A HISTOGRAM
if (knitr::is_html_output()) {
  
  Wts <- sort(BB$Weight)
  DP2 <- function(x){format(round(x, 1), nsmall = 1) }
  
  #  par(width=5, height=5, mar=c(0, 0, 2, 0) + 0.1, mfrow=c(2, 1))
  par(mfrow = c(2, 1))
  for (i in (1:length(Wts))){
    # Plot the data
    datagrid.x <- 11
    datagrid.y <- 4
    
    plot( expand.grid(1:datagrid.x, 
                      1:datagrid.y),
          type = "n",
          ylim = c(0.7, 4.3),
          main = "The data (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    text(  1:datagrid.x, 
           rep(4, 11), 
	   DP2(Wts[1:11]) )
    text(  1:datagrid.x, 
           rep(3, 11), 
	   DP2(Wts[12:22]) )
    text(  1:datagrid.x, 
           rep(2, 11), 
	   DP2(Wts[23:33]) )
    text(  1:datagrid.x, 
           rep(1, 11), 
	   DP2(Wts[34:44]) )
    
    # Now highlight the observation being used
    if (i <= 11) {
      points(  i, 4, 
               pch = 1, 
	             col = plotDark, 
	             lwd = 3, 
	             cex = 4 )
    }
    if ( (i > 11) & (i <= 22)){
      points(  i - 11, 3, 
               pch = 1, 
	             col = plotDark, 
	             lwd = 3, 
	             cex = 4 )
    }
    if ( (i > 22) & (i <= 233)){
      points(  i - 22, 2, 
               pch = 1, 
	             col = plotDark, 
	             lwd = 3, 
	             cex = 4 )
    }
    if ( (i > 33) ){
      points(  i - 33, 1, 
               pch = 1, 
	             col = plotDark, 
	             lwd = 3, 
	             cex = 4 )
    }
    
    # Produce the histogram
    
    hist( Wts[1:i],
          right = FALSE,
          xlab = "Birthweight (in kg)",
          ylab = "Number of babies",
          main = "Histogram of birthweights of babies born\non one day in a Brisbane hospital",
          las = 1,
          xlim = c(1.5, 4.5),
          ylim = c(0, 20),
          breaks = seq(1.5, 4.5, by = 0.5),
          col = plot.colour)
       points(Wts[i],
          0,
          pch = 19,
          col = plot.colour,
          cex = 2)
  }
} 
```


```{r BBHist1, fig.align="center", fig.width=10, fig.height=3, out.width='95%', fig.cap="Starting to make the histogram for the baby-birth data: the first 6 observations added" }
if (knitr::is_latex_output()){
  Wts <- sort(BB$Weight)

  par(mfrow = c(1, 2))
  
  i <- 6
    # Plot the data
    datagrid.x <- 11
    datagrid.y <- 4
    
    plot( expand.grid(1:datagrid.x, 
                      1:datagrid.y),
          type = "n",
          ylim = c(0.7, 4.5),
          xlim = c(1, 13),
          main = "The birthweights (in order)",
          xlab = "",
          ylab = "",
          axes = FALSE)
    text(  1:datagrid.x, 
           rep(4, 11), 
	   DP2(Wts[1:11]) )
    text(  1:datagrid.x, 
           rep(3, 11), 
	   DP2(Wts[12:22]) )
    text(  1:datagrid.x, 
           rep(2, 11), 
	   DP2(Wts[23:33]) )
    text(  1:datagrid.x, 
           rep(1, 11), 
	   DP2(Wts[34:44]) )
    
    # Now highlight the observation being used
    if (i <= 11) {
      points(  i, 4, 
               pch = 1, 
	             col = plotDark, 
	             lwd = 3, 
	             cex = 4 )
    }
    
    # Produce the histogram
    hist( Wts[1:i],
          right = FALSE,
          xlab = "Birthweight (in kg)",
          ylab = "Number of babies",
          main = "Histogram of birthweights of babies born\non one day in a Brisbane hospital",
          las = 1,
          xlim = c(1.5, 4.5),
          ylim = c(0, 20),
          breaks = seq(1.5, 4.5, 
                       by = 0.5),
          col = plot.colour)
    
    # Annotate histogram
    arrows(x0 = 1.75,
           x1 = 1.75,
           y0 = 19,
           y1 = 3,
           length = 0.15,
           angle = 15)
    text(1.75, 19,
         pos = 4,
         labels = "One observation: 1.7")
    
    arrows(x0 = 2.25,
           x1 = 2.25,
           y0 = 14,
           y1 = 5,
           length = 0.15,
           angle = 15)
    text(2.25, 14,
         pos = 4,
         labels = "Four observations: 2.1, 2.1, 2.2, 2.4")  
    
    arrows(x0 = 2.75,
           x1 = 2.75,
           y0 = 9,
           y1 = 2,
           length = 0.15,
           angle = 15)
    text(2.75, 9,
         pos = 4,
         labels = "One observation so far: 2.6")  
}
```


```{r BBHist44, fig.align="center", fig.width=5, fig.height=3, out.width='45%', fig.cap="The final histogram for the baby-birth data" }
if (knitr::is_latex_output()){
  Wts <- sort( BB$Weight )
  
  # Produce the histogram
  
  hist( Wts,
        right = FALSE, # Right-closed and left-open (i.e., like 2.5 to under 3 kg)
        xlab = "Birthweight (in kg)",
        ylab = "Number of babies",
        main = "Histogram of birthweights of babies\non one day in Brisbane",
        las = 1,
        xlim = c(1.5, 4.5),
        ylim = c(0, 20),
        breaks = seq(1.5, 4.5, by = 0.5),
        col = plot.colour)
  
}
```


::: {.example #Histograms2 name="Histograms"}
A study of *brain freezes* after eating cold food or drinking a cold drink measured the duration of the brain freeze [@data:Mages2017:BrainFreeze].

A histogram of the data (Fig.\ \@ref(fig:BrainFreeze), based on @data:Mages2017:BrainFreeze, Fig. 2b), shows 11 people experience HICS symptoms less than 5 seconds in length; nine people experienced symptoms for at least 5 but less than 10 seconds; and 1 person experienced symptoms for at least 35 seconds but under 40 seconds.
:::

```{r BrainFreeze, fig.cap="Duration of HICS (brain freeze) after drinking ice water", fig.width=5.5, fig.height=3.5, fig.align="center"}
# Read from Mages 2017, Figure 2b

Hts <- c(11, 9, 4, 4, 2, 5, 0, 1)
Borders <- seq(0, 40, by=5)

old.par <- par()

barplot(Hts,
        space = 0,
        ylim = c(0, 12),
        col = plot.colour,
        main = "Duration of brain freeze after drinking ice water",
        ylab = "Number of volunteers",
        xlab = "Duration (s)",
        axes = FALSE)
axis(side = 2, 
     las = 1)
axis(side = 1, 
     at = 0:8, 
     labels = Borders)
box()
```
   

### Describing the distribution {#SummariseData}

Graphs are constructed to help readers understand the data.
After producing a graph for one quantitative variable, describe these features:

1. The *average*:
   What is an "average" value, or central value?
2. The *variation*: 
   How much variation is present in the bulk of the data?
3. The *shape*:
   What is the shape of the distribution?
   That is, are most of the values smaller or larger, or about even distributed between smaller and larger values?
4. Mention any *outliers* (observations unusually large or small) or unusual features.

::: {.example #DescribeQuantData name="Describing quantitative data"}
The weights of babies (displayed in Example\ \@ref(exm:Histograms)) are typically between about 2.5kg and 3kg. (the *average*), with most between 1.5kg and 4.5 kg (*variation*).
A few babies have very lower weights, probably premature births (*shape*).
No unusual values are present.
:::

Describing the *shape* can be difficult, but sometimes terminology helps:

* *Right* (or *positively*) skewed: most data are smaller, with some larger values.
* *Left* (or *negatively*) skewed: most the data are larger, with some smaller values.
* Symmetric data (and perhaps bell-shaped):
  Approximately equal numbers of values are smaller and larger.
* Bimodal data: The distribution has two peaks.

`r if (knitr::is_html_output()) {
  'The carousel below (click the left and right arrows to move through the example plots) shows typical shapes.'
} else {
  'Figure\\ \\@ref(fig:ShapeDescriptionExamples) shows typical shapes.'
}`
Sometimes, no short descriptions are suitable.


```{r, child = if (knitr::is_html_output())  './children/ShapeExampleCarousel.Rmd'}
```

```{r, child = if (knitr::is_latex_output()) './children/ShapeExampleImages.Rmd' }
```


```{r}
if (FALSE){
par(mfrow = c(2, 2))

set.seed(110110)

xx1 <- scales::rescale( c( rchisq(300, df = 2), 
                           7, 7, 8, 7, 7, 7, 8, 7, 7, 6), 
                        to = c(0.01, 4.99) )

xx2 <- scales::rescale( c( rchisq(300, df = 2), 
                           7, 7, 8, 7, 7, 7, 8, 7, 7, 6), 
                        to = c(0.01, 4.99) )
xx2 <- max(xx1) - xx1

xx3 <- scales::rescale( rnorm(500, mean = 10), 
                        to = c(0.01, 4.99) )

xxA <- rnorm(500,
             mean = 0, sd = 1)
xxB <- rnorm(250,
             mean = 4, 
             sd = 1)
xx4 <- scales::rescale( c( xxA[ xxA > 0 ], xxB), 
                        to = c(0.01, 4.99) )

plotA <- ggplot( data.frame(x = xx1), 
              aes( x = x) ) +
              geom_histogram( binwidth = 0.5, 
                              color = "black", 
                              fill = "blue", 
                              alpha = 0.2, 
                              breaks = seq(0, 5, by = 0.5)) +
              scale_y_continuous(breaks = seq(0, 140, by = 40) ) +
              labs( x = "Observations", 
                    y = "Number of obs.", 
                    title = "Skewed right") +
              theme_light() 
plotB <- ggplot( data.frame(x = xx2), 
              aes(x = x) ) +
              geom_histogram(binwidth = 0.5, 
                             color = "black", 
                             fill = "blue", 
                             alpha = 0.2, 
                             breaks = seq(0, 5, by = 0.5)) +
              scale_y_continuous(breaks = seq(0, 140, by = 40) ) +
              labs( x = "Observations", 
                    y = "Number of obs.", 
                    title = "Skewed left") +
              theme_light() 
plotC <- ggplot( data.frame(x = xx3), 
              aes(x = x) ) +
              geom_histogram(binwidth = 0.5, 
                             color = "black", 
                             fill = "blue", 
                             alpha = 0.2, 
                             breaks = seq(0, 5, by = 0.5)) +
              scale_y_continuous(breaks = seq(0, 120, by = 40) ) +
              labs( x = "Observations", 
                    y = "Number of obs.", 
                    title = "Bell-shaped") +
              theme_light() 
plotD <- ggplot( data.frame(x = xx4), 
              aes(x = x) ) +
              geom_histogram(binwidth = 0.5, 
                             color = "black", 
                             fill = "blue", 
                             alpha = 0.2, 
                             breaks = seq(0, 5, by = 0.5)) +
              scale_y_continuous(breaks = seq(0, 120, by = 40) ) +
              labs( x = "Observations", 
                    y = "Number of obs.", 
                    title = "Bimodal") +
              theme_light() 
cowplot::plot_grid(plotA, plotB, plotC, plotD)

###################

hist( scales::rescale(xx, to = c(0, 5)), 
      xlab = "Values", 
      ylab = "Number of obs.",
      main = "Skewed right",
      axes = FALSE,
      col = plot.colour, 
      las = 1)
axis(side = 2, 
     las = 1)
box()


xx <- max(xx) - xx
hist( scales::rescale(xx, to = c(0, 5)), 
      xlab = "Values", 
      ylab = "Number of obs.",
      main = "Skewed left",
	    axes = FALSE,
      col = plot.colour, 
	    las = 1)
axis(side = 2, 
     las = 1)
box()

hist( scales::rescale(xx, to = c(0, 10)), 
      xlab = "Values", 
      ylab = "Number of obs.",
      main = "Approx. symmetric",
	    axes = FALSE,
      col = plot.colour, 
	    las = 1)
axis(side = 2, 
     las = 1)
box()



hist( scales::rescale(xx, to = c(0, 10)),  
	    right = TRUE, 
	    axes = FALSE,
      col = plot.colour, 
	    las = 1,
      xlab = "Values",
      ylab = "Number of obs.", 
      main = "Bimodal distribution")
axis(side = 2, 
     las = 1)
box()

par(mfrow = c(1, 1))
}
```


::: {.example #BimodalFaithful name="Bimodal data"}
The *Old Faithful* geyser in Yellowstone National Park (USA) erupts regularly [@hardle1991smoothing].
The time between eruptions (Fig.\ \@ref(fig:BimodalFaithfulHistoChangeBins), left panel) is bimodal, with peaks near 55 minutes and 80 minutes.
:::


```{r}
data(faithful) # Built-in data set
```

The number of boxes used in the histogram can change the impression of the distribution.
Software makes it easy to try different bin sizes to find one that suitably displays the overall distribution.
```{r, child = if (knitr::is_html_output())  './children/HistogramBins/histogramChangeBinsHTML.Rmd'}
```
```{r, child = if (knitr::is_latex_output()) './children/HistogramBins/histogramChangeBinsLATEX.Rmd'}
```



::: {.thinkBox .think data-latex="{iconmonstr-light-bulb-2-240.png}"}
Describe the histogram in Fig.\ \@ref(fig:BrainFreeze), the brain freeze durations.\label{thinkBox:DescribeBrainFreeze}

`r if (knitr::is_latex_output()) '<!--'`
`r webexercises::hide()`
* *Average*: Hard to be sure... maybe between 10 or 15. (More observations appear at the smaller values (as the bars are higher).) 
* *Variation*: From about 0 to about 40. 
* *Shape*: Slightly skewed right. 
* *Outliers*: No outliers or unusual observations. The observation between 35 and 40 *may* be an outlier.
  I suspect it is *not* an outlier, as a larger sample may very well have observations between 30 and 35.
  Of course, I could be wrong.
`r webexercises::unhide()`
`r if (knitr::is_latex_output()) '-->'`
:::


## One qualitative variable {#GraphsOneQual}

For qualitative data, like quantitative data, a graph shows the *distribution* of the data: the values taken by the variable, and how often they appear.

The three options for graphing qualitative data are:

* [*Dot chart*](#DotChartsOneQual):
  Usually a good choice.
* [*Bar chart*](#BarCharts):
  Usually a good choice.
* [*Pie chart*](#PieCharts):
  Only useful in special circumstances, and can be harder to interpret.

For *nominal* data, the levels of the variables could be displayed alphabetically, by size, by personal preference, or any other way: use the order most likely to be useful to readers.
For *ordinal* data, the natural order of the levels should almost always be used.

Sometimes these graphs are also used for discrete quantitative data with a small number of possible options.


### Dot charts  (qualitative data) {#DotChartsOneQual}

Dot charts indicate the counts (or the corresponding percentages) in each level, using dots on a line *starting at zero*.
The levels can be on the horizontal or vertical axis, and the counts or percentages in the other direction.
Placing the level names on the vertical axis often makes for easier reading, and space for long labels. 


::: {.example #DotPlotsQual name="Dot plots"}
A study of spider monkeys [@data:Chapman1990:SpiderMonkeys] examined the types of social groups present (Fig.\ \@ref(fig:SpiderMonkeys)).
The most common social group has many females plus offspring  (with almost 50 social groups).
Also, none of the commonly-observed social group include males.
:::


```{r SpiderMonkeys, fig.cap="Dot chart of spider monkey family groups", fig.align="center", fig.width=5.5, fig.height=3.5}
counts <- c(8, 3, 2, 15, 1, 23, 48)
names <- c("Solitary", 
           "All males", 
           "Female + no young", 
           "Mixed young", 
           "Mixed + no young", 
           "One female + offspring", 
           "Many females + offspring")

sort.order <- sort( counts, 
                    index.return = TRUE)

old.mar = par()$mar
new.mar <- old.mar
new.mar[2] <- old.mar[2] * 2.75
par(mar = new.mar )

dotchart( counts[sort.order$ix], 
          pch = 19,
          las = 1, 
          xlim = c(0, 50),
          xlab = "Percentage of observations",
          ylab = "",
          main = "Spider monkey sub-groups")
axis(side = 2, 
     seq_along(counts[ sort.order$ix ]), 
     names[sort.order$ix], 
     las = 1)
abline(v = 0,
	     col = "grey")   # Add zero line for clarity

```

For dot charts:

* Place the qualitative variable on the horizontal or vertical axis (and label with the levels of the variable).
* Use counts or percentages on the other axis.
* For nominal data, the levels can be ordered any way: *Think about the most helpful order*.


### Bar charts {#BarCharts}

Bar charts indicate the counts in each category using a bar *starting from zero*.
As with dot charts, the levels can be on the horizontal or vertical axis, but placing the level names on the vertical axis often makes for easier reading, and room for long labels. 


::: {.example #IndependenceBarText name="Bar charts"}
In a study of functional independence [@data:Ocepek2013:SmartHome], the type of diagnoses were graphed using a bar chart (Fig.\ \@ref(fig:IndependenceBar)).
For example, two people in the sample have cerebral palsy.
The reason for the different coloured bars is explained in Fig.\ \@ref(fig:IndependencePie).
:::


```{r IndependenceBar, fig.cap="Diagnoses of participants: numbers (left panel) and percentges (right panel)", fig.align="center", out.width='100%', fig.height=3.5, fig.width=10}
Number <- rev( c(11, 10, 9, 8, 5, 4, 3, 3, 2, 4) )
Diag <- rev( c("Amputation-one leg",
               "Neuromuscular disease",
              "SCI-paraplegia",
              "Rheum. dis. or multiple skeletal injury",
              "SCI-tetraplegia",
              "Cerebral vascular insult",
              "Amputation-both legs",
              "Amputation-arm",
              "Cerebral palsy",
              "Other") )
Diag.Cols <- rev( c( 
                grey(0.8), 
                grey(0.2), 
                grey(0.8),
                grey(0.8),
                grey(0.2),
                grey(0.8),
                grey(0.2),
                rep( grey(0.8), 3) ) )

par( mfrow = c(1, 2),
     mar = c(4, 16, 4, 2) + 0.1)
barplot( Number, 
         names.arg = Diag,
         horiz = TRUE,
         las = 1,
         xlim = c(0, 12),
         main = "Diagnoses of\nparticipants",
         xlab = "Number of participants",
         col = Diag.Cols)


barplot( Number/sum(Number) * 100 , 
         names.arg = Diag,
         horiz = TRUE,
         las = 1,
         xlim = c(0, 20),
         main = "Diagnoses of\nparticipants",
         xlab = "Percentage of participants",
         col = Diag.Cols)
```


For bar charts:

* Place the qualitative variable on the horizontal or vertical axis (and label with the levels of the variable).
* Use counts or percentages on the other axis.
* For nominal data, levels can be ordered any way: *Think about the most helpful order*.
* Bars have gaps between bars, as the bars represent distinct categories.
  In contrast, the bars in histograms are butted together (except when an interval has a count of zero), as the bars represent a numerical scale.



### Pie charts {#PieCharts}

In pie charts, a circle is divided into segments proportional to the number in each level of the qualitative variable.


::: {.example #PieCharts name="Pie charts"}
In a study of functional independence [@data:Ocepek2013:SmartHome], the *severity* of the diagnoses were shown in a pie chart (Fig.\ \@ref(fig:IndependencePie)).
This picture conveys only one piece of information ("69% of patients had a less severe injury"), so a graph of any kind is probably unnecessary.

The pie chart colours are used in the bar chart in Example\ \@ref(exm:IndependenceBarText), indicating the severity of the diagnoses.
This is called *encoding* extra information into the bar chart.
:::

   
```{r IndependencePie, fig.cap="Severity of diagnoses of participants", fig.align="center", fig.height=2.75, out.width='40%'}
par( mar = c(0.5, 0.5, 3.5, 0.5) )
pie( c(31, 69),
     labels = c("More severe (31)", 
                "Less severe (69%)"),
     col = c( grey(0.2), 
              grey(0.8)),
     main = "Severity of diagnoses of participants",
     init.angle = 90)
```

Pie charts present challenges:

* Pie charts only work when graphing parts of a whole.
* Pie charts only work when *all* options are present ('exhaustive').
* Pie charts are difficult to use with levels having zero counts, or small counts.
* Pie charts are difficult to read with many categories present.
* Pie charts are hard to read: Humans comparing *lengths* (as in bar and dot charts) better than *angles* (as in pie charts) [@data:Friel:Graphs].


::: {.thinkBox .think data-latex="{iconmonstr-light-bulb-2-240.png}"}
In which of these situations is a pie chart *appropriate*?\label{thinkBox:WhenPiechart}
  
1. The percentage of people who use these web browsers: Firefox, Chrome, and Safari.
2. For each country in Asia, the percentage of people who own an iPhone.

`r if (knitr::is_latex_output()) '<!--'`
`r webexercises::hide()`
1. A pie chart is *not suitable*.
   The three browsers are not *mutually exclusive* (people can use more than one of these browsers) nor *exhaustive* (some people may use browsers not listed). 
   For example, the percentages *could* be that 65% use Firefox, 84% use Chrome, and 20% use Safari. 
   These add to more than 100%.
2. A pie chart is *not suitable*, as the percentages are not parts of a whole.
   Each individual (person) has information recorded on *two* qualitative variables: (a) which state the person lives in (many levels); and (b) whether or not they own an iPhone ('yes' or 'no').  
   For example, the percentages *could* be 53% in Queensland, 61% in NSW, 41% in Victoria, and so on.
   They could possibly add to more than 100%.
`r webexercises::unhide()`
`r if (knitr::is_latex_output()) '-->'`
:::



### Comparing pie charts and bar charts {#CompareBarPie}

Consider the pie chart in 
`r if (knitr::is_latex_output()) {
   'Fig.\\ \\@ref(fig:LCgraphsPDF) (left panel),'
} else {
   'the top panel of Fig.\\ \\@ref(fig:LCgraphsHTML),'
}`
based on data in @data:andersen:1977.
The pie chart displays the number of lung cancer deaths in Fredericia (Denmark) from 1968 and 1971 inclusive, for various age groups (qualitative).
A pie chart *is* appropriate: only one variable is recorded on each individual (the age of each individual person), and the counts are parts of a whole.
However, determining *which* age groups have the most lung cancer deaths is hard.

The equivalent bar chart 
`r if (knitr::is_latex_output()) {
   '(right panel)'
} else {
   '(bottom panel)'
}`
makes the comparison easy: clearly the age groups '65 to 69' and 'Over 74' have slightly fewer deaths than the other age groups.


```{r LCgraphsHTML, fig.cap="Graphs of deaths from lung cancer in Fredericia (Denmark) from 1968 and 1971", fig.align="center", fig.width=5, fig.height=8}
if (knitr::is_html_output()) {
  
  data(DanishLC)
  
  par( mfrow = c(2, 1), 
       mar = c(5.5, 5, 4.1, 1))
  
  pie( DanishLC$Cases[DanishLC$City == "Fredericia"],
       labels = c("40 to 54", 
                  "55 to 59", 
                  "60 to 64", 
                  "65 to 69", 
                  "70 to 74", 
                  "Over 74"),
       col = viridis(6),
       mar = c(0.5, 0.5, 4.1, 0.5),
       main = "Pie chart of lung cancer deaths\nby age group")
  
  barplot(DanishLC$Cases[DanishLC$City == "Fredericia"], 
          las = 2, 
          names.arg = c("40 to 54", 
                        "55 to 59", 
                        "60 to 64", 
                        "65 to 69", 
                        "70 to 74", 
                        "Over 74"),
          col = viridis(6),
          ylim = c(0, 12),
          sub = "Age group", 
          ylab = "Number of deaths", 
          main = "Bar chart of lung cancer deaths\nby age group")
}
```

```{r LCgraphsPDF, fig.cap="Graphs of deaths from lung cancer in Fredericia (Denmark) from 1968 and 1971 (inclusive)", fig.align="center", fig.width=10, fig.height=3.5, out.width='95%'}
if (knitr::is_latex_output()) {
  data(DanishLC)
  
  par( mar = c(0.5, 0.5, 0.5, 0.5) )
  
  nf <- layout( matrix(c(1, 2), 
                       ncol = 1),
                widths = c(3, 1) )
  
  par( mfrow = c(1, 2), 
       mar = c(5.5, 5, 4.1, 1))
  
  pie( DanishLC$Cases[DanishLC$City == "Fredericia"],
       labels = c("40 to 54", 
                  "55 to 59", 
                  "60 to 64", 
                  "65 to 69", 
                  "70 to 74", 
                  "Over 74"),
       col = viridis(6),
       main = "Pie chart of lung cancer deaths\nby age group")
  
  barplot(DanishLC$Cases[DanishLC$City == "Fredericia"], 
          las = 2, 
          names.arg = c("40 to 54", 
                        "55 to 59", 
                        "60 to 64", 
                        "65 to 69", 
                        "70 to 74", 
                        "Over 74"),
          col = viridis(6),
          ylim = c(0, 12),
          sub = "Age group", 
          ylab = "Number of deaths", 
          main = "Bar chart of lung cancer deaths\nby age group")
}
```


Recall that the *purpose of a graph is to display information in the clearest, simplest possible way, to help the reader understand the message(s) in the data*.
A pie chart often makes the message hard to see [@siegrist1996use].


::: {.example #ComparingGraphs name="Comparing graphs"}
In the NHANES study [@data:NHANES3:Data], the age of each participant was recorded.
Rank the age groups from largest group to smallest group using each graph in Fig.\ \@ref(fig:3Dcharts), all constructed from the same data.
Which graph makes comparing the category sizes easier?
:::

```{r 3Dcharts, fig.cap="Four different graphs for the same data", fig.height=6, out.width='85%',fig.align="center"}
data(NHANES) # NHANES package
AgeD <- NHANES$AgeDecade
AgeD.levels <- levels(NHANES$AgeDecade)
num.levels <- length(AgeD.levels)

levels(AgeD) <- c( AgeD.levels[1:(num.levels - 2)],
                   rep(" 60+", 2))
cols <- viridis(9)[3:8]

# New data
par(mfrow = c(2, 2)) #, mar=c(3, 1, 1, 1))

plotrix::pie3D( table(AgeD),
                labels = levels(AgeD),
                main = "Number in each age group",
                theta = pi/4,
                labelcex = 0.9,
                col = cols)

pie(table(AgeD),
    las = 1,
    col = cols,
    main = "Number in each age group")

barplot(table(AgeD),
        las = 2,
        ylim = c(0, 1600),
        col = cols,
        xlab = "Age group",
        cex.names = 0.9,
        ylab = "Number in sample",
        main = "Number in each age group")

abline(h = seq(0, 1600, by = 250), 
       col = "white")


dotchart( as.vector(table(AgeD)),
          labels = names(table(AgeD)),
          las = 2,
          cex = 1,
          pch = 19,
          xlim = c(0, 1600),
          col =  cols,
          ylab = "Age group",
          xlab = "Number in sample",
          main = "Number in each age group")
```

<iframe src="https://learningapps.org/watch?v=pf4om4k5t22" style="border:0px;width:100%;height:500px" allowfullscreen="true" webkitallowfullscreen="true" mozallowfullscreen="true"></iframe>


## Tabular summaries {#GraphNeeded}
 
Graphs are generally excellent for summarising data.
However, sometimes data can be collated in a table, or a *tabular summary*.
For *qualitative* data, the tabular summary records the numbers (or percentages, or both) for each level of the variable (Example\ \@ref(exm:NHANESdata)).
For *quantitative* data, the tabular summary effectively produces a histogram as a table (Example\ \@ref(exm:BabyBoomTable)).

With small amounts of information, sometimes just writing the information is better ('69% of diagnoses were less severe'; Example\ \@ref(exm:PieCharts)).



::: {.example #NHANESdata name="Tables for data"}
The NHANES age data in Example\ \@ref(exm:ComparingGraphs) may be displayed using a tabular summary (Table\ \@ref(tab:SimpleTable)) rather than a graphical summary.
:::

```{r SimpleTable}
simpletab <- t(
  rbind("Age group" = levels(AgeD),
        Number = table(AgeD), 
        Percentage = format(table(AgeD)/sum(table(AgeD)) * 100, 
                            nsmall = 1, # Number of digits to right of decimal point
                            digits = 3))
)

if( knitr::is_latex_output() ) {
  kable(simpletab,
        format = "latex",
        longtable = FALSE,
        align = c("r", "r", "r"),
        booktabs = TRUE,
        col.names = c("Age group",
                      "Number of people",
                      "Percentage of people"),
        linesep = c("","","", "\\addlinespace"),
        row.names = FALSE,
        caption = "The NHANES age distribution, displayed in a tabular summary") %>%
    kable_styling(font_size = 10) %>%
    row_spec(0, bold = TRUE)
}
if( knitr::is_html_output() ) {
  kable(simpletab,
        format = "html",
        longtable = FALSE,
        align = c("r", "r", "r"),
        booktabs = TRUE,
        col.names = c("Age group",
                      "Number of people",
                      "Percentage of people"),
        row.names = FALSE,
        caption = "The NHANES age distribution, displayed in a tabular summary") %>%
    kable_styling(full_width = FALSE)
}
```



:::{.example #BabyBoomTable}
In Example\ \@ref(exm:Histograms), a histogram was constructed for the weight of babies born in a Brisbane hospital in one day.
The data can also be presented as the number of babies within each weight category, just as with the histogram (Table\ \@ref(tab:BabyBoomTable)).
The histogram is essentially a picture of this table.
:::



```{r BabyBoomTable}
Wts <- BB$Weight
out <- hist(Wts,
            plot = FALSE,
            right = FALSE)
BBpercentage <- out$counts / sum(out$counts) * 100
  
BBtab <- array(dim = c(6, 3) )
BBtab[, 2] <- out$counts 
BBtab[, 3] <- round(BBpercentage) 
BBtab[, 1] <- c("1.5kg to under 2.0kg",
                "2.0kg to under 2.5kg",
                "2.5kg to under 3.0kg",
                "3.0kg to under 3.5kg",
                "3.5kg to under 4.0kg",
                "4.0kg to under 4.5kg") 
colnames(BBtab) <- c("Weight group", 
                     "Number of babies born",
                     "Percentage of babies born")

if( knitr::is_latex_output() ) {
  kable(BBtab,
        format = "latex",
        longtable = FALSE,
        align = c("c", "r", "r"),
        booktabs = TRUE,
        linesep = c("","","\\addlinespace"),
        row.names = FALSE,
        caption = "The baby-weights data, displayed in a tabular summary") %>%
	kable_styling(font_size = 10) %>%
	row_spec(0, bold = TRUE)
}
if( knitr::is_html_output() ) {
 kable(BBtab,
               format = "html",
               longtable = FALSE,
               align = c("c", "r", "r"),
               booktabs = TRUE,
               row.names = FALSE,
               caption = "The baby-weights data, displayed in a tabular summary") %>%
      kable_styling(full_width = FALSE)
}
```


## One qualitative variable and one quantitative variable {#GraphOneQualOneQuant}

Relationships between *one qualitative variable* and *one quantitative variable* can be displayed using: 

* [*Back-to-back stemplots*](#BackToBackStem):
  Best for small amounts of data when the qualitative variable has *two levels*; 
* [*2-D dot charts*](#TwoDDot): 
  Best choice for small to moderate amounts of data;
* [*Boxplots*](#Boxplot):
  Best choice, except for small amounts of data.        


### Back-to-back stemplot {#BackToBackStem}

Back-to-back stemplots are two stemplots (Sect.\ \@ref(StemAndLeafPlots)) sharing the same stems; one group has the leaves emerging left-to-right from the stem, and the second group has the leaves emerging right-to-left from the stem.
Back-to-back stemplots can only be used when *two* groups are being compared.


::: {.example #KrillData name="Back-to-back stemplots"}
A study of krill [@data:Greenacre2016:reporting] produced the data in 
`r if (knitr::is_latex_output()) {
   'Table\\ \\@ref(tab:KrillDataTable).'
} else {
   'Table\\ \\@ref(tab:KrillDataTableHTML).'
}`
A back-to-back stemplot allows the two groups to be compared visually 
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:KrillEggsDoubleStem)).'
} else {
   '(Fig.\\ \\@ref(fig:KrillEggsDoubleStemHTML)).'
}`
The *Treatment* leaves go from right-to-left, and the leaves for the *Control* group left-to-right, sharing the same stems.
The control group tends to produce more eggs, in general.
:::


```{r}
Eggs.T <- c(0, 0, 1, 1, 3, 8, 8, 12, 18, 21, 26, 30, 35, 48, 50)
Eggs.C <- c(0, 0, 0, 0, 1, 1, 1,  2,  2,  3,  8, 16, 20, 26, 31)

KrillEggs2 <- cbind( Eggs.T[1:8], 
                     c(Eggs.T[9:15], ""), 
                     Eggs.C[1:8], 
		     c(Eggs.C[9:15], "") )

```

`r if (knitr::is_html_output()) '<!--'`
\begin{figure}
\begin{minipage}{0.5\textwidth}
\captionof{table}{The number of eggs laid by krill, for those in a treatment group and for those in a control group\label{tab:KrillDataTable}}
\fontsize{10}{12}\selectfont
```{r KrillDataTable}

if( knitr::is_latex_output() ) {
  kable( KrillEggs2,
         format = "latex",
         longtable = FALSE,
         table.env = "@empty",
         booktabs = TRUE,
         linesep = c("", "", "\\addlinespace"),
         align = c("c","c")) %>%
         #col.names = c("Treatment group", "Control group"),
    add_header_above(header = c("Treatment group " = 2, 
                                " Control group" = 2), 
                     bold = TRUE, 
                     align = "c")
}
```
\end{minipage}
\hspace{0.05\textwidth}
\begin{minipage}{0.50\textwidth}%
\centering
```{r KrillEggsDoubleStem, fig.width=4.0, fig.height=3, out.width='100%'}

  stemgrid.x <- 1
  stemgrid.y <- 6
  stemcounts <- array(1, dim=stemgrid.y)

  par( mar = c(0.25, 0.25, 4, 0.25) + 0.1)

  plot( expand.grid(1:stemgrid.x, 0:stemgrid.y), 
      type = "n",
      xlab = "",
      ylab = "",
      main = "Stemplot of egg counts of krill\n(1|6 means 16 eggs)",
      xlim = c(-11, 7),
      ylim = c(-0.2, 6.2),
      axes = FALSE)
  
### Treatment data in cols 1  to 11 (x = -11 to x = -1)
### Stems in col 12  (x=0)
### Control data in cols 13 to 19 (x= 1 to x = 7)
  
  # Divide the stems and leaves
  lines(c(-0.5, -0.5),
        c(0.5, 6.5),
         lwd = 2,
         col = "grey")
  lines(c(0.5, 0.5),
        c(0.5, 6.5),
         lwd = 2,
         col = "grey")
  text(-5, 0, 
       "Treatment data", 
       font = 2)
  text( 4, 0, 
        "Control data", 
	font = 2)
  
  # Plot stems
  points( c(0, 0, 0, 0, 0, 0), 
          1:6,
          pch = as.character(5:0))
  # Plot leaves
  for (j in (1:2) ){

    if (j == 1) { # Right-hand side
      Eggs <- Eggs.C
      step <- 1
      stemcounts <- array(1, dim=stemgrid.y)
    }
    if (j == 2) { # Left-hand side
      Eggs <- Eggs.T
      step <- -1
      stemcounts <- array(-1, dim = stemgrid.y)
    }
    
    for (i in 1:length(Eggs)){
    
      if ( Eggs[i] < 10) {
        points( stemcounts[1], 6, 
	        pch = as.character(Eggs[i]))
        stemcounts[1] <- stemcounts[1] + step
      }

    if ( (Eggs[i] >= 10) & (Eggs[i] < 20) ) {
      points( stemcounts[2], 5, 
              pch = as.character(Eggs[i] - 10) )
      stemcounts[2] <- stemcounts[2] + step
    }

    if ( (Eggs[i] >= 20) & (Eggs[i] < 30) ) {
      points( stemcounts[3], 4, 
              pch = as.character(Eggs[i] - 20))
      stemcounts[3] <- stemcounts[3] + step
    }
    if ( ( Eggs[i] >= 30 ) & (Eggs[i] < 40) ) {
      points( stemcounts[4], 3, 
              pch = as.character(Eggs[i] - 30))
      stemcounts[4] <- stemcounts[4] + step
    }
    if ( ( Eggs[i] >= 40 ) & (Eggs[i] < 50) ) {
      points( stemcounts[5], 2, 
              pch = as.character(Eggs[i] - 40))
      stemcounts[5] <- stemcounts[5] + step
    }
    if ( ( Eggs[i] >= 50 ) & (Eggs[i] < 60) ) {
      points( stemcounts[6], 1, 
              pch = as.character(Eggs[i] - 50))
      stemcounts[6] <- stemcounts[6] + step
    }
  }
}
```
\caption{The number of eggs from krill, for treatment and control groups}\label{fig:KrillEggsDoubleStem}
\end{minipage}
\end{figure}
`r if (knitr::is_html_output()) '-->'`


`r if (knitr::is_latex_output()) '<!--'`
```{r KrillEggsDoubleStemHTML, fig.width=6.0, out.width='100%'}

  stemgrid.x <- 1
  stemgrid.y <- 6
  stemcounts <- array(1, dim=stemgrid.y)
  
  par( mar = c(0.5, 0.5, 4, 0.5) + 0.1)

  plot( expand.grid(1:stemgrid.x, 0:stemgrid.y), 
      type = "n",
      xlab = "",
      ylab = "",
      main = "Stemplot of egg counts of krill\n(1|6 means 16 eggs)",
      xlim = c(-11, 7),
      ylim = c(-0.2, 6.2),
      axes = FALSE)
  
### Treatment data in cols 1  to 11 (x = -11 to x = -1)
### Stems in col 12  (x=0)
### Control data in cols 13 to 19 (x= 1 to x = 7)
  
  # Divide the stems and leaves
  lines(c(-0.5, -0.5),
        c(0.5, 6.5),
         lwd = 2,
         col = "grey")
  lines(c(0.5, 0.5),
        c(0.5, 6.5),
         lwd = 2,
         col = "grey")
  text(-5, 0, 
       "Treatment data", 
       font = 2)
  text( 4, 0, 
        "Control data", 
	font = 2)
  
  # Plot stems
  points( c(0, 0, 0, 0, 0, 0), 
          1:6,
          pch = as.character(5:0))
  # Plot leaves
  for (j in (1:2) ){

    if (j == 1) { # Right-hand side
      Eggs <- Eggs.C
      step <- 1
      stemcounts <- array(1, dim=stemgrid.y)
    }
    if (j == 2) { # Left-hand side
      Eggs <- Eggs.T
      step <- -1
      stemcounts <- array(-1, dim = stemgrid.y)
    }
    
    for (i in 1:length(Eggs)){
    
      if ( Eggs[i] < 10) {
        points( stemcounts[1], 6, 
	        pch = as.character(Eggs[i]))
        stemcounts[1] <- stemcounts[1] + step
      }

    if ( (Eggs[i] >= 10) & (Eggs[i] < 20) ) {
      points( stemcounts[2], 5, 
              pch = as.character(Eggs[i] - 10) )
      stemcounts[2] <- stemcounts[2] + step
    }

    if ( (Eggs[i] >= 20) & (Eggs[i] < 30) ) {
      points( stemcounts[3], 4, 
              pch = as.character(Eggs[i] - 20))
      stemcounts[3] <- stemcounts[3] + step
    }
    if ( ( Eggs[i] >= 30 ) & (Eggs[i] < 40) ) {
      points( stemcounts[4], 3, 
              pch = as.character(Eggs[i] - 30))
      stemcounts[4] <- stemcounts[4] + step
    }
    if ( ( Eggs[i] >= 40 ) & (Eggs[i] < 50) ) {
      points( stemcounts[5], 2, 
              pch = as.character(Eggs[i] - 40))
      stemcounts[5] <- stemcounts[5] + step
    }
    if ( ( Eggs[i] >= 50 ) & (Eggs[i] < 60) ) {
      points( stemcounts[6], 1, 
              pch = as.character(Eggs[i] - 50))
      stemcounts[6] <- stemcounts[6] + step
    }
  }
}
```
`r if (knitr::is_latex_output()) '-->'`



### 2-D dot charts {#TwoDDot}

A 2-D dot chart places a dot for each observation, separated for each level of the qualitative variable (also see Sect.\ \@ref(DotChartsOneQual)).
For the same krill data in Example\ \@ref(exm:KrillData), a dot chart is shown in Fig.\ \@ref(fig:TwoDDotchart).

Many observations are the same, so some points would be *over-plotted* if points were not *stacked* (left panel), or some randomness added (a 'jitter') in the vertical direction, before plotting (right panel).


```{r TwoDDotchart, fig.height=6, out.width='100%',fig.cap="Two variations of a 2-D dot chart for the krill-egg data to avoid overplotting: stacking (left) and jittering (right)", fig.align="center", fig.width=8.75, fig.height=3}
old.par <- par()
par(mar = c(5, 6, 4, 2) + 0.1)
par(mfrow = c(1,2))

stripchart( list(Treatment = Eggs.T, 
                 Control = Eggs.C), 
            main = "A 2-D dot chart: the number of eggs\nlaid by krill (points stacked)",
            method = "stack", 
            las = 1,
            pch = 19,
            xlim = c(0, 50),
            cex = 0.7,
            ylim = c(0.6, 2.4),
            xlab = "Number of eggs")

stripchart( list(Treatment = Eggs.T, 
                 Control = Eggs.C), 
            main = "A 2-D dot chart; the number of eggs\nlaid by krill (points jittered)",
            method = "jitter",
            jitter = 0.18, 
            las = 1,
            pch = 19,
            xlim = c(0, 50),
            cex = 0.7,
            ylim = c(0.6, 2.4),
            xlab = "Number of eggs")

#par( old.par )
```


### Boxplots {#Boxplot}

For the krill data used in Example\ \@ref(exm:KrillData), a boxplot is shown in Fig.\ \@ref(fig:BoxplotKrill).
Understanding boxplots is challenging, so boxplots will be discussed again in Sect.\ \@ref(VariationIQR).


```{r BoxplotKrill, fig.cap="A boxplot  for the krill-egg data", fig.align="center",  fig.height=3.5, , fig.width=5}
boxplot( list(Treatment = Eggs.T, 
              Control = Eggs.C), 
            main = "A boxplot of the\nnumber of eggs laid by krill",
            las = 1,
            pch = 19,
            col = plot.colour,
            ylim = c(0, 50),
            ylab = "Number of eggs")
```


```{r}
KrillT.quantiles <- quantile(Eggs.T)
Eggs.C2 <- Eggs.C[ -which.max(Eggs.C)]
KrillC2.quantiles <- quantile(Eggs.C2)
```

To explain boxplots, first focus on just *one* boxplot from Fig.\ \@ref(fig:BoxplotKrill): the boxplot for the *Treatment* group (the left box).
Boxplots have five horizontal lines; from the top to the bottom of the plot (Fig.\ \@ref(fig:BoxplotKrillTreatment)):

1. *Top line*:
   The *largest* number of eggs is `r KrillT.quantiles[5]`.
2. *Second line from top*:
   75% of the observations are smaller than about `r KrillT.quantiles[4]`, represented by the line at the top of the central box. 
   This is called the *third quartile* ($Q_3$).
3. *Middle line*:
   50% of the observations are smaller than about `r KrillT.quantiles[3]`, represented by the line in the centre of the central box. 
   This is an 'average' value for the data, the *second quartile* ($Q_2$).
4. *Second line from bottom*:
   25% of the observations are smaller than about `r KrillT.quantiles[2]`, represented by the line at the bottom of the central box. 
   This is called the *first quartile* ($Q_1$).
5. *Bottom line*:
   The *smallest* number of eggs is `r KrillT.quantiles[1]`. 


```{r BoxplotKrillTreatment, fig.cap="A boxplot for the krill-egg data; the boxplot and dotplot just for the treatment group. Thwe circled numbers on the left refer to the explantions in-text.", fig.width=6, fig.height=6.5, out.width='60%', fig.align="center" }
Eggs.quantiles <- quantile(Eggs.T)
out <- boxplot( list(Treatment = Eggs.T), 
            main = "A boxplot of the number of\neggs laid by krill: Treatment group",
            las = 1,
            pch = 19,
            col = "white",
            ylim = c(-4, 50),
            xlim = c(0.75, 2.75),
            ylab = "Number of eggs")
abline( h = Eggs.quantiles, 
        col = "grey",
        lty = 2,
        lwd = 1)
bxp(out, # Replot to put over the dashed lines
    las = 1, 
    add = TRUE, 
    boxfill = plot.colour)

text( 1.5, Eggs.quantiles[5], 
      "Maximum value", 
      pos = 1, 
      cex = 0.9)
text( 1.5, Eggs.quantiles[4], 
      "75% smaller than this", 
      pos = 3, 
      cex = 0.9)
text( 1.5, Eggs.quantiles[3], 
      "Middle value\n('average')", 
      pos = 3, 
      cex = 0.9)
text( 1.5, Eggs.quantiles[2], 
      "25% smaller\nthan this", 
      pos = 3, 
      cex = 0.9)
text( 1.5, Eggs.quantiles[1], 
      "Minimum value", 
      pos = 1, 
      cex = 0.9)
stripchart( list(Treatment = Eggs.T), 
            vertical = TRUE, 
            at = 2.0, 
            add = TRUE, 
            method = "stack",
            pch = 19,
            col = plot.colour0)
text( x = 2.5,
      y = c( mean(Eggs.quantiles[2:3]),
             mean(Eggs.quantiles[3:4]),
             mean(Eggs.quantiles[4:5])
             ),
      "25% of data\n in here")
text( x = 2.5,
      y = mean(Eggs.quantiles[1:2]),
      "25% in here")

arrows(x0 = 2.25,
       y0 = Eggs.quantiles[1],
       x1 = 2.25,
       y1 = Eggs.quantiles[2],
       code = 3, # Arrow at both ends
       length = 0.1,
       angle = 15)
arrows(x0 = 2.2,
       y0 = Eggs.quantiles[2],
       x1 = 2.2,
       y1 = Eggs.quantiles[3],
       code = 3, # Arrow at both ends
       length = 0.1,
       angle = 15)
arrows(x0 = 2.25,
       y0 = Eggs.quantiles[3],
       x1 = 2.25,
       y1 = Eggs.quantiles[4],
       code = 3, # Arrow at both ends
       length = 0.1,
       angle = 15)
arrows(x0 = 2.2,
       y0 = Eggs.quantiles[4],
       x1 = 2.2,
       y1 = Eggs.quantiles[5],
       code = 3, # Arrow at both ends
       length = 0.1,
       angle = 15)

# Add values to right axis
Eggs.quantiles.labs <- Eggs.quantiles # So we get all the labels appearing; these two are very close, so only one gets printed
Eggs.quantiles.labs[1] <- Eggs.quantiles.labs[1] - 1
Eggs.quantiles.labs[2] <- Eggs.quantiles.labs[2] + 1
axis(side = 4, # This places the tick marks
     at = Eggs.quantiles,
     labels = rep("", 5) )
mtext(Eggs.quantiles,
      side = 4,
      at = Eggs.quantiles.labs,
      las = 1,
      line = 1,
      cex = 1)

axis(side=1)

points(x = rep(0.75, length(Eggs.quantiles) ),
      y = Eggs.quantiles,
      pch = 1,
      cex = 2.5)
text( x = rep(0.75, length(Eggs.quantiles) ),
      y = Eggs.quantiles,
      label = rev(c("1", "2", "3", "4", "5")),
      cex = 0.9
      )
```


The box for the *Control* group krill is slightly different (Fig\. \@ref(fig:BoxplotKrill)): One observation is identified with a point, *above* the top line.
Computer software has identified this observation as potentially unusual (in this case, unusually *large*), and has plotted this point separately.
(Unusually large or small observations are called *outliers*.)
The values of the quantiles ($Q_1$, $Q_2$ and $Q_3$) are computed as usual (with outliers included).

For the *Control* data, the largest observation (31 eggs) is deemed unusually large (explained in Sect.\ \@ref(OutliersIQRrule)), and then the boxplot is constructed as follows (Fig.\ \@ref(fig:BoxplotKrillControl)):

```{r}
outC <- boxplot( list(Treatment = Eggs.C),
                 plot = FALSE)
Eggs.quantiles <- outC$stats
```


* *Top line*:
  The *largest* number of eggs (*excluding* the outlier of 31) is about `r Eggs.quantiles[5]`: 
  This is the line at the top of the boxplot.
* *Second line from top*:
  75% of the observations (*including* the 31 eggs) are smaller than about `r Eggs.quantiles[4]`, represented by the line at the top of the central box. 
  This is called the *third quartile*, or $Q_3$.
* *Middle line*:
  50% of the observations (*including* the 31 eggs) are smaller than about `r Eggs.quantiles[3]`,  represented by the line in the centre of the central box.
  This is an 'average' value for the data, the *second quartile*, or $Q_2$.
* *Second line from bottom*:
  25% of the observations (*including* the 31 eggs) are smaller than about `r Eggs.quantiles[2]`,  represented by the line at the bottom of the central box.
  This is called the *first quartile*, or $Q_1$.
* *Bottom line*:
  The smallest number of eggs is `r Eggs.quantiles[1]`. 
  This is the line at the bottom of the boxplot.

Clearly, decimal numbers of eggs (such as $Q_2 = `r Eggs.quantiles[2]`$) are not observable; however, with `r length(Eggs.C)` observations, finding exactly the value for which 25% of observations are smaller is not possible.
Software uses approximations to compute these values (and different software may use different rules, and hence may produce slightly different values for quartiles).


```{r BoxplotKrillControl, fig.cap="A boxplot  for the krill-egg data; the boxplot just for the control group", fig.width=5, fig.align="center"}
outC <- boxplot( list(Treatment = Eggs.C), 
            main = "A boxplot of the number of\neggs laid by krill: Control group",
            las = 1,
            pch = 19,
            col = plot.colour,
            ylim = c(-4, 35),
            xlim = c(0.7, 1.75),
            ylab = "Number of eggs")
Eggs.quantiles <- outC$stats

# Now, repeat without the outlier:
outC2 <- boxplot( list(Treatment = Eggs.C[ -length(Eggs.C)]),
                  plot = FALSE)
Eggs2.quantiles <- outC2$stats

abline( h = outC$stats[-3, 1], 
        col = "grey",
        lty = 2,
        lwd = 1)
lines( c(0, 1), 
       c(outC$stats[3, 1], outC$stats[3, 1]), 
       lty = 2, 
       col = "grey")

bxp(outC, # Replot to put over the dashed lines
    las = 1, 
    add = TRUE, 
    boxfill = plot.colour)


text( 1.4, 
      Eggs.quantiles[5, 1], 
      "Max. value (apart from outlier)", 
      pos = 1, 
      cex = 0.85)
text( 1.5, 
      Eggs.quantiles[4, 1], 
      "75% smaller than this", 
      pos = 3, 
      cex = 0.85)
text( 1.0, 
      Eggs.quantiles[3, 1], 
      "Middle value", 
      pos = 3, 
      cex = 0.85)
text( 1.5, 
      Eggs.quantiles[2, 1], 
      "25% smaller than this", 
      pos = 3, 
      cex = 0.85)
text( 1.5, 
      Eggs.quantiles[1, 1], 
      "Min. value", 
      pos = 1, 
      cex = 0.85)

arrows(1.2, max(Eggs.C), 1.03, 
       max(Eggs.C), 
       col = "grey", 
       lty = 1, 
       angle = 15)
text(1.4, 
     max(Eggs.C), 
     "Unusually large value", 
     cex = 0.85)

# Add values to right axis
Eggs.quantiles.labs <- Eggs.quantiles # So we get all the labels appearing; these two are very close, so only one gets printed
Eggs.quantiles.labs[1] <- Eggs.quantiles.labs[1] - 1.5
Eggs.quantiles.labs[2] <- Eggs.quantiles.labs[2] + 1.5
Eggs.quantiles.labs[3] <- Eggs.quantiles.labs[2] + 1.5

Eggs.quantiles[3] <- NA  # Remove median, as too cluttered
Eggs.quantiles.labs[3] <- NA  # Remove median, as too cluttered


axis(side = 4, # This places the tick marks
     at = c(Eggs.quantiles, 31), # Adding the outlier, too
     labels = rep("", 6) )
mtext(c(Eggs.quantiles, "31"),
      side = 4,
      at = c(Eggs.quantiles.labs, 31),
      las = 1,
      line = 1,
      cex = 1)


axis(side=1)

```


::: {.importantBox .important data-latex="{iconmonstr-warning-8-240.png}"}
Boxplots summarise data with only five numbers, so a lot of detail of the distribution is lost.
For this reason, boxplots are excellent for *comparing* distributions, but histograms are better for displaying a single quantitative distribution.
:::


::: {.example #BoxplotsExplained name="Boxplots explained"}
The NHANES study collects large amounts of information from about 10\,000 Americans each year (Sect.\ \@ref(NHANESGraphs)).
Consider the boxplot of ages
`r if (knitr::is_latex_output()) {
   ' from the 1999--2004 participants (Fig.\\ \\@ref(fig:BoxplotAge; the online version has an animation).'
} else {
   ' from the 1999--2004 participants. The animation below shows how the boxplot of the age of the Americans in the sample is constructed.'
}`
The "average" age of the subjects is about 38 years, and the ages range from almost zero to about 80 years of age.
:::


```{r animation.hook="gifski", interval=1, dev=if (is_latex_output()){"pdf"}else{"png"}}
data(NHANES)

# Grab some of the necessary stats
Qs <- quantile(NHANES$Age, 
               seq(0, 1, by = 0.25)) # See below for better way when there are outliers

left <- 0.75
right <- 1.25
axis.edge <- axis.edge.left <- 0.5
axis.edge.right <- 1.5
angle <- 15
mid <- mean( c(left, axis.edge.left) )
mid.r <- mean( c(right, axis.edge.right) )

# The box is at x = 1, and goes from x = 0.5 to 1.5 on the x-axis.
# See this by following the boxplot()) call with axis(side-1)

if (knitr::is_html_output()) {
  for (i in 1:6){
    boxplot( NHANES$Age, 
             col = plot.colour,
             las = 1,
             lwd = 2,
             ylim = c(0, 85),
             main = "Boxplot of Age (NHANES data)",
             ylab = "Age (in years)")
    
    if ( i == 2){
      arrows(left, Qs[5],
             axis.edge, Qs[5],
             angle = angle,
             lwd = 2)
      text(mid, Qs[5], 
           label = "Maximum", 
           pos = 1)
    }
    if ( i == 3 ){
      
      arrows(left, Qs[4],
             axis.edge, Qs[4],
             angle = angle,
             lwd = 2)
      text(mid, Qs[4], 
           label = expression( italic(Q)[3]), 
           pos = 1)
    }
    if (i == 4){
      arrows(left, Qs[3],
             axis.edge, Qs[3],
             angle = angle,
             lwd = 2)
      text(mid, Qs[3], 
           label = expression( italic(Q)[2]), 
           pos = 1)
    }
    if (i == 5){  
      arrows(left, Qs[2],
             axis.edge, Qs[2],
             angle = angle,
             lwd = 2)
      text(mid, Qs[2], label = expression( italic(Q)[1]), 
           pos = 1)
    }
    if ( i == 6){    
      arrows(left, Qs[1],
             axis.edge, Qs[1],
             angle = angle,
             lwd = 2)
      text(mid, Qs[1], 
           label = "Minimum", 
           pos = 3)
    }
    
    xL <- -10
    xR <- 10
    polygon( x = c(xL, xL, xR, xR),
             y = c(Qs[1], Qs[2], Qs[2], Qs[1]),
             col = makeTransparent(plot.colour2),
             border = NA)
    polygon( x = c(xL, xL, xR, xR),
             y = c(Qs[3], Qs[4], Qs[4], Qs[3]),
             col = makeTransparent(plot.colour2),
             border = NA)
    
    text(x = 1.35,
         y = c( mean(Qs[1:2]),
                mean(Qs[2:3]),
                mean(Qs[3:4]),
                mean(Qs[4:5])),
         c("25% of the data\nin this region"))
  }
}
```


```{r BoxplotAge, fig.align="center", fig.width=5,  fig.height=3.5, fig.cap="The boxplot for the age of people in the NHANES data" }
if (knitr::is_latex_output()){
    boxplot( NHANES$Age, 
             col = plot.colour,
             las = 1,
             ylim = c(-3, 85),
             main = "Boxplot of Age (NHANES data)",
             ylab = "Age (in years)")
}
```


::: {.example #BoxplotsHorizontal name="Boxplots"}
Boxplots can be plotted horizontally too, which leaves space for long labels.
In Fig.\ \@ref(fig:CementBuildBoxplot) (left panel; based on @data:Silva2016:rootcanal), the three dental cements are very different regarding their push-out forces.
:::


```{r CementBuildBoxplot, fig.cap="Comparing three push-out values for three dental cements", fig.align="center", fig.height=2.75, fig.width=4.0, fig.show='hold', out.width=c('49%', '49%')}
MTA    <- c(0.05, 0.1, 0.2, 0.25, 0.45)
Endo   <- c(0.1, 0.2, 0.4,   1.0, 2.05)
AHplus <- c(0.5, 1.1, 2.25, 2.9, 4.4)

CementStats <- list( n = c(5, 5, 5),
                     stats = matrix( cbind(MTA, 
                                           Endo, 
                                           AHplus), 
                                     ncol = 3, 
                                     byrow = FALSE),
                     conf = matrix( c(-1.5, 7, -2, 10, 3, 3), 
		                 ncol = 2, 
				             byrow = FALSE),  ### NOT USEFUL, only for notches
                  out = NULL,
                  group = NULL,
                  names = c("MTA Fillapex", 
		                        "Endo Seal", 
			                      "AH Plus"))
par( mar = c(4, 6, 3, 2) + 0.1)
plot( c(0, 5), 
      c(0.4, 3.6),
      type = "n",
      axes = FALSE,
      main = "Boxplot of push-out force",
      ylab = "",
      xlab = "Push-out force (MPa)")

axis(side = 1,
     las = 1)
abline(h = 0, 
       col = "grey", 
       lwd = 2)
bxp(CementStats,
    las = 1,
    horizontal = TRUE,
    add = TRUE,
    boxfill = plot.colour)

######

DB  <-  c( 7,   3, 1.5, 1, -1.5)
DBB <- c(10, 5.5, 3.5, 2, -2)

EngStats <- list( n = c(5, 5),
                  stats = matrix( cbind(DB, DBB), 
                                  ncol = 2, 
                                  byrow = FALSE),
                  conf = matrix( c(-1.5, 7, -2, 10), 
                                 ncol = 2, 
                                 byrow = FALSE),
                  out = c(8, 20, 12.5, 11.5, -7),
                  group = c(1, 2, 2, 2, 2),
                  names = c("DB", 
                            "DBB"))
plot( c(0.4, 2.6), 
      c(-10, 20),
      type = "n",
      axes = FALSE,
      main = "Boxplot of project costs",
      ylab = "Cost increase (%)",
      xlab = "Project delivery method")

axis(side = 2,
     las = 1)
abline(h = 0, 
       col = "grey", 
       lwd = 2)
bxp(EngStats,
    las = 1,
    add = TRUE,
    boxfill = plot.colour)
```


::: {.example #Boxplots2 name="Boxplots"}
A study of different engineering project delivery methods [@data:Hale2009:ProjectDelivery] produced the boxplot in Fig.\ \@ref(fig:CementBuildBoxplot) (right panel): the increase in the costs of projects seem to differ between the two methods.
The quantitative variable is the percentage increase in costs, and the qualitative variable is the method used.

The Design/Build (DB) method produces a smaller project cost growth on average (the centre line of the boxplot), but the Design/Bid/Build (DBB) method produces more variation in project cost growth.
Outliers are present for both methods, as indicated by the dots.
:::



<div style="text-align:center;">
```{r}
htmltools::tags$video(src = "./videos/Boxplots.mp4", 
                      width = "550", 
                      controls = "controls", 
                      loop = "loop", 
                      style = "padding:5px; border: 2px solid gray;")
```
</div>



## Two quantitative variables {#Scatterplots}

Scatterplots display the relationship between *two quantitative variables*.
Conventionally, when appropriate, the response variable is on the vertical axis, and the explanatory variable is on the horizontal axis.
As with any graph, describing the message in the graph is important, because the purpose of a graph is to display the information in the clearest, simplest possible way.
*Scatterplots* can be described by briefly explaining how the variables are related to each other.
Scatterplots are studied further in Sect.\ \@ref(Scatterplots).



<div style="float:right; width: 222x; border: 1px; padding:10px">
<img src="Illustrations/pexels-jim-de-ramos-1263349.jpg" width="200px"/>
</div>


::: {.example #Scatterplots2 name="Scatterplots"}
A study of athletes at the *Australian Institute of Sport* (AIS) measured numerous physical and blood measurements from high performance athletes [@data:Telford1991:sexsportsize].
Many relationships were of interest; Fig.\ \@ref(fig:ScatterAIS) shows the relationship between the sum of skin folds (SSF) and percentage body fat.
Each point represents the percentage body fat and the SSF for one athlete.
Clearly, the greater the percentage body fat, the greater the sum of skin folds, in general.
:::

	
```{r ScatterAIS, results ='hide', fig.cap="Scatterplot of SSF against percentage body fat for male AIS athletes", fig.width=4.5, fig.height=3.75, fig.align="center", out.width='50%'}
data(AISsub)

plot( SSF ~ PBF,
      data = AISsub,
      las = 1,
      xlim = c(5, 20),
      ylim = c(20, 120),
      main = "SSF vs PBF: Male AIS athletes",
      xlab = "Percentage body fat",
      ylab = "Sum of skin fold",
      pch = 19)
title(sub = "(Note: Units of measurement for SSF not given)", 
      cex.sub = 0.85)
```


## Two qualitative variables {#TwoQualVars}

The relationship between two qualitative variables can be explored using:

* *Stacked bar charts* (Sect.\ \@ref(StackedBarCharts));
* *Side-by-side bar charts* (Sect.\ \@ref(SideBySideBarCharts)); or
* *Dot charts* (Sect.\ \@ref(TwoWayCountsDotCharts)).

As an example, a study of road kill [@data:Russell2009:Possums] produced the data in Table\ \@ref(tab:PossumRoadkillNumber).
The two qualitative variables are the *season* (ordinal, with four levels) and the *sex* (nominal, with three levels including 'Unknown').


```{r PossumRoadkillNumber}
## Russell2009

Counts <- matrix( c(75, 25, 21,   
                    74, 27, 22,   
                    71, 10, 18,   
                    58, 10, 12), 
                  nrow = 4, 
		  ncol = 3, 
		  byrow = TRUE) 
rownames(Counts) <- c("Autumn", 
                      "Winter", 
                      "Spring", 
                      "Summer")
colnames(Counts) <- c("Unknown", 
                      "M", 
                      "F")

Counts2 <- Counts
rownames(Counts2) <- c("Aut.", 
                       "Wint.", 
                       "Spr.", 
                       "Sum.")

Percentages <- t( apply(Counts2, 1,
                  function(x){ round(x / sum(x) * 100, 1)}))
colnames(Percentages) <- c("Unknown", "Males", "Females")


if( knitr::is_latex_output() ) {
  T1 <- knitr::kable(Counts,
               format = "latex",
               longtable = FALSE,
               align = c("r", "c", "c", "c"),
               booktabs = TRUE) %>%
    #kable_styling(font_size = 10) %>%
    row_spec(0, bold = TRUE)
  
  rownames(Percentages) <- c("Autumn", 
                             "Winter",
                             "Spring",
                             "Summer")
    T2 <- knitr::kable(Percentages,
               format = "latex",
               longtable = FALSE,
               align = c("r", "r", "r"),
               booktabs = TRUE) %>%
    #kable_styling(font_size = 10) %>%
    row_spec(0, bold = TRUE)
    
   out <- knitr::kables(list(T1, T2),
                       format = "latex",
                       label = "PossumRoadkillNumber",
                       caption = "The number of possums found as road kill, by sex and season (left). The percentages of possums found as road kill by sex, within each season (i.e., rows sum to 100\\%) (right).") %>% 
    kable_styling(font_size = 10)
  out2 <- prepareSideBySideTable(out) 
  out2  
}
if( knitr::is_html_output() ) {
  knitr::kable(Counts,
               format = "html",
               longtable = FALSE,
               align = c("r", "c", "c", "c"),
               booktabs = TRUE,
               caption = "The number of possums found as road kill, by sex and season") %>%
     kable_styling(full_width = FALSE)
}
```


### Stacked bar charts {#StackedBarCharts}

The data can be graphed by using a bar for each season, *stacking* the bars by sex on top of each other, within each season (Fig.\ \@ref(fig:PossumStackedBarChart), left panel).


```{r PossumStackedBarChart, fig.height=3.75, fig.width=9, out.width='100%', fig.align="center", fig.cap="The number of possums found as road kill, by sex and season. Left panel: A stacked barchart; right panel: A side-by-side barchart"}

par(mfrow = c(1, 2))

colnames(Counts2) <- c("?", "M", "F")

out <- barplot( t(Counts2), 
   beside = FALSE,
   ylab = "Number of possums",
   xlab = "Season",
   main = "Stacked barchart of\nroadkill possums",
   col = viridis(3),
   #space = c(0.15, 1),
   ylim = c(0, 125),
   xlim = c(0, 5),
   las = 1)
grid(nx = NA, 
     ny = NULL)

text( out[4] + 0.6, 70, 
      "F", 
      adj = 0)
text( out[4] + 0.6, 60, 
      "M",
      adj = 0)
text( out[4] + 0.6, 30, 
      "?",
      adj = 0)

text( out[1] - 0.8, 110, 
      "F", 
      adj = 0)
text( out[1] - 0.8,  90, 
      "M", 
      adj = 0)
text( out[1] - 0.8,  40,
      "?", 
      adj = 0)


######

out <- barplot( t(Counts2), 
   beside = TRUE,
   ylab = "Number of possums",
   xlab = "Season",
   main = "Side-by-side barchart of\nroadkill possums",
   col = viridis(3),
   space = c(0.15, 1),
   ylim = c(0, 80),
   las = 1)
grid(nx = NA, 
     ny = 6)

mtext( "?", 
       at = out[1, ], 
       side = 1, 
       line = 0, 
       cex = 0.8, 
       padj = 0.5)
mtext( "M", 
       at = out[2, ], 
       side = 1, 
       line = 0, 
       cex = 0.8, 
       padj = 0.5)
mtext( "F", 
       at = out[3, ], 
       side = 1, 
       line = 0, 
       cex = 0.8, 
       padj = 0.5)
```



### Side-by-side bar charts {#SideBySideBarCharts}

Instead of stacking the bars within each season on top of each other, the bars can be placed *side-by-side* within each season (Fig.\ \@ref(fig:PossumStackedBarChart), right panel).



### Dot charts {#TwoWayCountsDotCharts}

Instead of bars, dots (or other symbols) can be used in place of a side-by-side bar chart (Fig.\ \@ref(fig:PossumDotCharts)).


```{r PossumDotCharts, fig.height=4.5, fig.width=5.0, fig.cap="A dot chart of the number of possums found as road kill, by sex and season", fig.align="center"}
colnames(Counts) <- c("Unknown", 
                      "Males", 
                      "Females")

dotchart(Counts,
         main = "Dot chart of possum road deaths\nby season and sex of possum",
         xlim = c(0, 80),
 	       xlab = "Number of possums",
	       pch = c(1, 2, 3, 4))
```


### Other variations {#OtherVariations}

Many variations of these charts are possible, such as using:

* a stacked bar chart, side-by-side bar chart, or dot chart can be used;
* percentages or counts can be used;
* the counts (or percentage) can appear on either the horizontal or vertical axis;
* either the sex of the possum or the season can be used as the first division of the data.

Some of these variations are shown in 
`r if (knitr::is_latex_output()) {
   'Fig.\\ \\@ref(fig:PossumBarChartsLATEX).'
} else {
   'Fig.\\ \\@ref(fig:PossumBarChartsHTML).'
}`
A table is also possible, of either counts
`r if (knitr::is_latex_output()) {
   '(Table\\ \\@ref(tab:PossumRoadkillNumber), left table) or percentages (right  table).'
} else {
   'Table\\ \\@ref(tab:PossumRoadkillPercentages) or percentages (Table\\ \\@ref(PossumRoadkillPercentagesHTML)).'
}`


```{r PossumBarChartsLATEX, fig.height=8, fig.width=8.0, out.width='100%', fig.cap="Four graphs: The number or percentage of possums found as road kill, by sex and season, and a dot chart", fig.align="center"}

if (knitr::is_latex_output()) {
  par( mfrow = c(2, 2))
  
  colnames(Counts2) <- c("?", 
                         "M", 
                         "F")
  
  barplot( Counts2, 
     beside = TRUE,
     xlab = "Number of possums",
     main = "Roadkill possums",
     #ylab = "Season",
     space = c(0.15, 1),
     ylim = c(0, 16),
     xlim = c(0, 80),
     horiz = TRUE,
     axes = FALSE,
     names.arg = c("", "", ""),
     las = 1)
  grid(ny = NA, 
       nx = 8)
  
  out <- barplot( Counts2, # Plot again to cover the lines added by  grid()
     beside =TRUE,
     xlab = "Number of possums",
     #ylab = "Season",
     main = "Bar chart of possum road deaths\nby season and sex of possum",
     col = viridis(4),
     space = c(0.15, 1),
     ylim = c(0, 16),
     xlim = c(0, 80),
     las = 1,
     horiz = TRUE,
     names.arg = c("", "", ""),
     add = TRUE)
  
  mtext( "Au", 
         at = out[1, ], 
         side = 2, 
         line = 1, 
         las = 1)
  mtext( "Wn",
          at = out[2, ], 
  	side = 2, 
  	line = 1, 
  	las = 1)
  mtext( "Sp", 
         at = out[3, ], 
         side = 2, 
         line = 1, 
         las = 1)
  mtext( "Su", 
          at = out[4, ], 
  	side = 2, 
  	line = 1, 
  	las = 1)
  
  mtext( "Females", 
          at = 14, 
  	side = 2, 
  	line = 3, 
  	las = 3)
  mtext( "Males",   
         at = 8, 
         side = 2, 
         line = 3, 
         las = 3)
  mtext( "Unknown", 
         at = 3, 
         side = 2, 
         line = 3, 
         las = 3)
  
  
  Percentages <- t( apply(Counts2, 1,
                          function(x){ round(x / sum(x) * 100, 1)}))
  colnames(Percentages) <- c("Unknown", "Males", "Females")
  
  barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     main = "Roadkill possums",
     #xlab = "Season",
     space = c(0.15,1),
     names.arg = rep("", 3),
     ylim = c(0,80),
     las = 1)
  grid(nx = NA, 
       ny = 8)
  
  out <- barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     #xlab =" Season",
     main = "Bar chart of possum road deaths\nby season and sex of possum",
     space = c(0.15, 1),
     col = viridis(4),
     names.arg = rep("", 3),
     ylim=c(0, 80),
     las = 1,
     add = TRUE)
  
  mtext( "Au", 
         at = out[1, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Wn", 
         at = out[2, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Sp", 
         at = out[3, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Su", at=out[4,], side=1, line=1, cex=0.8, padj=0.5)
  
  mtext( "Unknown", 
         at = mean(out[, 1]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Males",   
         at = mean(out[, 2]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Females", 
         at = mean(out[, 3]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  
  barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     main = "Roadkill possums",
     #xlab = "Season",
     space = c(0.15, 1),
     names.arg = rep("", 3),
     ylim = c(0, 80),
     las = 1)
  grid(nx = NA, 
       ny = 8)
  
  out <- barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     #xlab = "Season",
     main = "Bar chart of possum road deaths\nby season and sex of possum",
     space = c(0.15, 1),
     names.arg = rep("", 3),
     col = viridis(4),
     #horiz = TRUE,
     ylim = c(0, 80),
     las = 1,
     add = TRUE)
  
  mtext( "Au", 
         at = out[1, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Wn", 
         at = out[2, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Sp", 
         at = out[3, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Su", 
         at = out[4, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  
  mtext( "Unknown", 
         at = mean(out[, 1]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Males",   
         at = mean(out[, 2]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Females", 
         at = mean(out[, 3]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  
  
  
colnames(Counts) <- c("Unknown", 
                      "Males", 
                      "Females")

dotchart( t(Counts),
         main = "Dot chart of possum road deaths\nby season and sex of possum",
         xlim = c(0, 80),
 	       xlab = "Number of possums",
	       pch = c(1, 2, 3, 4))
}
```

```{r PossumBarChartsHTML, fig.height=10, fig.width=4.0, fig.cap="Four graphs: The number or percentage of possums found as road kill, by sex and season, and a dot chart", fig.align="center"}

if (knitr::is_html_output()) {
  par( mfrow = c(3, 1))

  colnames(Counts2) <- c("?", 
                         "M", 
                         "F")
  
  barplot( Counts2, 
     beside = TRUE,
     xlab = "Number of possums",
     #ylab = "Season",
     main = "Roadkill possums",
     space = c(0.15, 1),
     ylim = c(0, 16),
     xlim = c(0, 80),
     horiz = TRUE,
     axes = FALSE,
     names.arg = c("", "", ""),
     las = 1)
  grid(ny = NA, 
       nx = 8)
  
  out <- barplot( Counts2, # Plot again to cover the lines added by  grid()
     beside =TRUE,
     xlab = "Number of possums",
     #ylab = "Season",
     main = "Bar chart of possum road deaths\nby season and sex of possum",
     space = c(0.15, 1),
     ylim = c(0, 16),
     xlim = c(0, 80),
     las = 1,
     horiz = TRUE,
     names.arg = c("", "", ""),
     add = TRUE)
  
  mtext( "Au", 
         at = out[1, ], 
         side = 2, 
         line = 1, 
         las = 1)
  mtext( "Wn",
          at = out[2, ], 
  	side = 2, 
  	line = 1, 
  	las = 1)
  mtext( "Sp", 
         at = out[3, ], 
         side = 2, 
         line = 1, 
         las = 1)
  mtext( "Su", 
          at = out[4, ], 
  	side = 2, 
  	line = 1, 
  	las = 1)
  
  mtext( "Females", 
          at = 14, 
  	side = 2, 
  	line = 3, 
  	las = 3)
  mtext( "Males",   
         at = 8, 
         side = 2, 
         line = 3, 
         las = 3)
  mtext( "Unknown", 
         at = 3, 
         side = 2, 
         line = 3, 
         las = 3)
  
  
  Percentages <- t( apply(Counts2, 1,
                          function(x){ round(x / sum(x) * 100, 1)}))
  colnames(Percentages) <- c("Unknown", "Males", "Females")
  
  barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     main = "Roadkill possums",
     #xlab = "Season",
     space = c(0.15,1),
     names.arg = rep("", 3),
     ylim = c(0,80),
     las = 1)
  grid(nx = NA, 
       ny = 8)
  
  out <- barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     #xlab =" Season",
     main = "Bar chart of possum road deaths\nby season and sex of possum",
     space = c(0.15, 1),
     names.arg = rep("", 3),
     ylim=c(0, 80),
     las = 1,
     add = TRUE)
  
  mtext( "Au", 
         at = out[1, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Wn", 
         at = out[2, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Sp", 
         at = out[3, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Su", at=out[4,], side=1, line=1, cex=0.8, padj=0.5)
  
  mtext( "Unknown", 
         at = mean(out[, 1]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Males",   
         at = mean(out[, 2]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Females", 
         at = mean(out[, 3]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  
  barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     main = "Roadkill possums",
     #xlab = "Season",
     space = c(0.15, 1),
     names.arg = rep("", 3),
     ylim = c(0, 80),
     las = 1)
  grid(nx = NA, 
       ny = 8)
  
  out <- barplot( Percentages, 
     beside = TRUE,
     ylab = "Percentages of possums",
     #xlab = "Season",
     main = "Bar chart of possum road deaths\nby season and sex of possum",
     space = c(0.15, 1),
     names.arg = rep("", 3),
     #horiz = TRUE,
     ylim = c(0, 80),
     las = 1,
     add = TRUE)
  
  mtext( "Au", 
         at = out[1, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Wn", 
         at = out[2, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Sp", 
         at = out[3, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  mtext( "Su", 
         at = out[4, ], 
         side = 1, 
         line = 1, 
         cex = 0.8, 
         padj = 0.5)
  
  mtext( "Unknown", 
         at = mean(out[, 1]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Males",   
         at = mean(out[, 2]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  mtext( "Females", 
         at = mean(out[, 3]), 
         side = 1, 
         line = 3, 
         cex = 1, 
         padj = 0.5)
  
  
  
  colnames(Counts) <- c("Unknown", 
                      "Males", 
                      "Females")

dotchart( t(Counts),
         main = "Dot chart of possum road deaths\nby season and sex of possum",
         xlim = c(0, 80),
 	       xlab = "Number of possums",
	       pch = c(1, 2, 3, 4))
}
```


<!-- tables are side-by-side in LaTeX, so oly need separate one for HTML -->
```{r PossumRoadkillPercentagesHTML}
colnames(Percentages) <- c("Unknown", "M", "F")

if( knitr::is_html_output() ) {
  knitr::kable(Percentages,
               format = "html",
               longtable = FALSE,
               align = c("r", "r", "r"),
               booktabs = TRUE,
               caption = "The percentages of possums found as road kill by sex, within each season (i.e., rows sum to 100%)") %>%
    kable_styling(full_width = FALSE)
}
```

::: {.thinkBox .think data-latex="{iconmonstr-light-bulb-2-240.png}"}
Of all these displays, which one do you think best communicates the message in the data?
:::


<iframe src="https://learningapps.org/watch?v=pgn3q7ptv22" style="border:0px;width:100%;height:500px" allowfullscreen="true" webkitallowfullscreen="true" mozallowfullscreen="true"></iframe>



## Other types of graphs {#GraphsOther}

Many types of graphs have been studied, but sometimes other plots are useful.
The same principle applies:

> The purpose of a graph is to display the information in the clearest, simplest possible way, to help the reader understand the message(s) in the data.

In this section, graphs for some other situations are discussed:

* *Geographic plots* (Sect.\ \@ref(GeographicPlots)): 
  Useful when the RQ concerns geographical regions.
* *Case-profile plots* (Sect.\ \@ref(CaseProfilePlot)): 
  Useful when the same units are measured over a small number of time points, or are otherwise connected.
* *Histogram of differences* (Sect.\ \@ref(HistoDiffPlot)): 
  Useful for within-subject comparisons for *two* comparisons.
* *Time plots* (Sect.\ \@ref(TimePlots)): 
  Useful when the units vary over a large number of time points.


### Geographic plots {#GeographicPlots}

Data with a geographic component can be plotted accordingly.


::: {.example #GeographicPlots name="Geographics plots"}
A study examining lower-limb amputation incidence in Australia (@data:Dillon2017:Geographic) produced the graphs in Fig.\ \@ref(fig:AmputationAustralia).
Clearly, the incidence of amputation is higher in the Northern Territory than other parts of Australia for both females and males.
Furthermore, males have higher incidence of amputation than females in every jurisdiction.
:::


```{r AmputationAustralia, fig.cap="Age-adjusted incidence of lower limb amputations in Australia, from August 2007 to December 2011: females (left) and males (right). Numbers are incidents per 100 000.",fig.width=11, fig.height=5, out.width='99%'}
library(ozmaps)

colList <- viridis::viridis(6)[3:6]

colF <- rep( colList[1], 9)
colF[7] <- colList[3]  # NT
colM <- rep( colList[2], 9)
colM[7] <- colList[4]  # NT

par( mfrow = c(1, 2),
     mar = c(4, 4, 4, 0.5))

# FEMALES
plot( x = c(110, 155),
      y = c(-45, -5),
      type = "n",
      asp = 1,
      las = 1,
      xlim = c(110, 155),
      ylim = c(-45, -5),
      main = "Limb amputations per 100,000:\nFemales",
      las = 1,
      xlab = "Longitude",
      ylab = "Latitude")

abline( v = seq(110, 160, by = 10),
         lwd = 1,
         col = "grey",
         lty = 2)
abline( h = seq(-10, -40, by = -10),
         lwd = 1,
         col = "grey",
         lty = 2)
ozmap(x = "states",
      add = TRUE,
      col = colF)

text(x = 133,
     y = -20,
     labels = "NT",
     col = grey(0.9))
text(x = 143, 
     y = -22, 
     labels = "Qld",
     col = grey(0.9))
text(x = 121, 
     y = -25, 
     labels = "WA",
     col = grey(0.9))
text(x = 147,
     y = -32,
     labels = "NSW",
     col = grey(0.9))
text(x = 143,
     y = -37,
     labels = "Vic",
     col = grey(0.9))
text(x = 151, 
    y = -44, 
    labels = "Tas", 
    col = grey(0))
text(x = 134, 
    y = -29, 
    labels = "SA", 
    col = grey(0.9))
legend("top",
       fill = colList,
       horiz = TRUE,
       box.col = "white",
       bg = "white",
       inset = 0.015,
       legend = c("15 to 35", "35 to 55", "55 to 75", "95 to 115"))

# MALES
plot( x = c(110, 155),
      y = c(-45, -5),
      type = "n",
      asp = 1,
      las = 1,
      xlim = c(110, 155),
      ylim = c(-45, -5),
      main = "Limb amputations per 100,000:\nMales",
      las = 1,
      xlab = "Longitude",
      ylab = "Latitude")

abline( v = seq(110, 160, by = 10),
         lwd = 1,
         col = "grey",
         lty = 2)
abline( h = seq(-10, -40, by = -10),
         lwd = 1,
         col = "grey",
         lty = 2)
ozmap(x = "states", 
      add = TRUE,
      col = colM)
text(x = 133,
     y = -20,
     labels = "NT",
     col = grey(0))
text(x = 143, 
     y = -22, 
     labels = "Qld",
     col = grey(0.9))
text(x = 121, 
     y = -25, 
     labels = "WA",
     col = grey(0.9))
text(x = 147,
     y = -32,
     labels = "NSW",
     col = grey(0.9))
text(x = 143,
     y = -37,
     labels = "Vic",
     col = grey(0.9))
text(x = 151, 
    y = -44, 
    labels = "Tas", 
    col = grey(0))
text(x = 134, 
    y = -29, 
    labels = "SA", 
    col = grey(0.9))
legend("top",
       fill = colList,
       horiz = TRUE,
       box.col = "white",
       bg = "white",
       inset = 0.015,
       legend = c("15 to 35", "35 to 55", "55 to 75", "95 to 115"))

```



### Case-profile plots {#CaseProfilePlot}

Sometimes the same variable is measured on each unit of analysis more than once (i.e., *within*-individual comparisons for each unit of analysis) but only a small number of times.
Then a *case-profile plot* can be used: a plot showing how the response variable changes for each unit of analysis.
Examples of this type of data include:

* Measurements of household water consumption before and after installing water-saving devices, for many households.
* Blood pressure measurements from left arms and right arms, for many people.

In both cases, the data are *paired* (*two* observations per unit of analysis) as each unit of analysis gets a pair of similar measurements.


::: {.example #CaseProfilePlots name="Case-profile plots"}
A study of children with atopic asthma [@data:Lothian2006:Whey] included a graph showing the immunogoblin E concentrations (IgE) before and after an intervention 
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:PairedGraphCaseProfileLATEX), left panel).'
} else {
   '(Fig.\\ \\@ref(fig:PairedGraphCaseProfileHTML), top panel).'
}`
Each line on the graph represents one person.
:::


```{r PairedGraphCaseProfileLATEX, fig.cap="The IgE data. Left: A case-profile plot. Each line represents one subject, joining that person's pre-intervention score to their post-intervention score; right: A histogram of the differences", fig.align="center", fig.width=9.25, fig.height=3.5, out.width='100%'}
Pre <-  c( 83, 292, 293, 623, 792, 1543, 1668, 1960, 2877, 2961, 5504)
Post <- c( 83, 292, 292, 542, 709, 1000, 1000, 1626, 2502, 2711, 4504)


WhichTime <- c(
  rep(1, length(Pre)),
  rep(2, length(Post))
)


if( knitr::is_latex_output() ) {
  par(mfrow  = c(1, 2) )
  
  plot( c(Pre, Post) ~ WhichTime,
        pch = 19,
        las = 1,
        axes = FALSE,
        ylim = c(0, 6000),
        xlim = c(0.8, 2.2),
        main = "Case-profile plot\nof the IgE data",
        xlab = "Time",
        ylab = "IgE (in micrograms/L)")
  axis(side = 1, 
       at = 1:2, 
       labels = c("Pre-intervention", "Post-intervention"))
  box()
  
  axis(side = 2, 
       las = 1)
  for (i in 1:length(Pre)){
    lines( c(1, 2),
           c(Pre[i], Post[i]),
           col = grey(0.6)
    )
  }
  
  
  Reduction <- Pre - Post
  
  IgE <- data.frame( IgE.Pre = Pre, 
                     IgE.Post = Post,
                     IgE.Reduction = Reduction)
  
  hist(Reduction,
       col = plot.colour,
       las = 1,
       breaks = seq(0, 1000, by = 100),
       ylab = "Number of subjects",
       xlab = "Reduction in IgE (in micrograms/L)",
       main = "Reduction in IgE\n pre and post supplementation")
  box()
}
```





```{r PairedGraphCaseProfileHTML, fig.cap="The IgE data. Top: A case-profile plot. Each line represents one subject, joining that person's pre-intervention score to their post-intervention score; bottom: A histogram of the differences", fig.align="center", fig.width=5, fig.height=7.5}
Pre <-  c( 83, 292, 293, 623, 792, 1543, 1668, 1960, 2877, 2961, 5504)
Post <- c( 83, 292, 292, 542, 709, 1000, 1000, 1626, 2502, 2711, 4504)


WhichTime <- c(
  rep(1, length(Pre)),
  rep(2, length(Post))
)

if( knitr::is_html_output() ) {
  
  par(mfrow  = c(2, 1) )
  
  plot( c(Pre, Post) ~ WhichTime,
        pch = 19,
        las = 1,
        axes = FALSE,
        ylim = c(0, 6000),
        xlim = c(0.8, 2.2),
        xlab = "Time",
        ylab = "IgE (in micrograms/L)")
  axis(side = 1, 
       at = 1:2, 
       labels = c("Pre-intervention", "Post-intervention"))
  box()
  
  axis(side = 2, 
       las = 1)
  for (i in 1:length(Pre)){
    lines( c(1, 2),
           c(Pre[i], Post[i]),
           col = grey(0.6)
    )
  }
  
  
  Reduction <- Pre - Post
  
  IgE <- data.frame( IgE.Pre = Pre, 
                     IgE.Post = Post,
                     IgE.Reduction = Reduction)
  
  hist(Reduction,
       col = plot.colour,
       las = 1,
       breaks = seq(0, 1000, by = 100),
       ylab = "Number of subjects",
       xlab = "Reduction in IgE (in micrograms/L)",
       main = "Reduction in IgE\n pre and post supplementation")
  box()
}
```



### Histogram of differences {#HistoDiffPlot}

An alternative way to present paired data (two observations made for each unit of analysis) is to produce a histogram of the *changes* for each individual.
Consider the person in the case-profile plot whose line is at the top in 
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:PairedGraphCaseProfileLATEX)).'
} else {
   '(Fig.\\ \\@ref(fig:PairedGraphCaseProfileHTML)).'
}`
Their 'pre' IgE level is about 5500 micrograms/L, and their 'post' IgE level is about 4500 micrograms/L, a *change* (or more descriptively, a *reduction*)  of about 1000 micrograms/L.
These reductions could be computed for each person (Table\ \@ref(tab:IgEChanges)), and a histogram constructed based on these IgE *reductions*s, with one reduction for each person
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:PairedGraphCaseProfileLATEX), right panel).'
} else {
   '(Fig.\\ \\@ref(fig:PairedGraphCaseProfileHTML), bottom panel).'
}`


```{r IgEChanges}
Reduction <- Pre - Post

IgE <- data.frame( IgE.Pre = Pre, 
                   IgE.Post = Post,
                   IgE.Reduction = Reduction)


if( knitr::is_latex_output() ) {
  kable(IgE,
        format = "latex",
        booktabs = TRUE,
        col.names = c("IgE (before) in micrograms/L", 
                      "IgE (after) in micrograms/L", 
                      "IgE reduction in micrograms/L"),
        linesep = c("", "", "", "\\addlinespace", "", "", "\\addlinespace", "", "", "", ""), # Otherwise adds a space after five lines... 
        caption = "The IgE before and after an intervention, and the change in IgE (in micrograms/L)",
        digits = 0)  %>%
    kable_styling(font_size = 10) %>%
    row_spec(0, bold = TRUE)
}
if( knitr::is_html_output() ) {
  kable(IgE,
        format = "html",
        longtable = TRUE,
        booktabs = TRUE,
        col.names = c("IgE (before) in micrograms/L", 
                      "IgE (after) in micrograms/L", 
                      "IgE reduction in micrograms/L"),
        linesep = c("", "", "", "\\addlinespace", "", "", "\\addlinespace", "", "", "", ""), # Otherwise addes a space after five lines... 
        caption = "The IgE before and after an intervention, and the change in IgE (in micrograms/L)",
        digits = 0)
}
```



::: {.example #InsulationGraph name="Graphing paired data"}
The Electricity Council in Bristol examined if a certain type of wall-cavity insulation was effective in reducing energy consumption in winter [@data:OpenUni:insulationBA].
Their RQ was: "Is there an *average reduction* in energy consumption after adding insulation?"

The data collected, shown below, can be graphed using a case-profile plot (Fig.\ \@ref(fig:DataInsulationGraphPlots), left panel).
A dashed line has been used to show an *increase* in energy usage, and a solid line for houses where energy was *saved* after installing insulation (*encoding* extra information.)
:::

```{r DataInsulationGraph}
data(InsulationBeforeAfter)

InsulationBeforeAfter$Diff <- InsulationBeforeAfter$Before - InsulationBeforeAfter$After
```

For these data, finding the difference in energy consumption for each house is possible, since the same unit of analysis is measured twice on the same variable: energy consumption *before* and *after* adding insulation.
The difference (or *saving*) in energy consumption for each house can be graphed using a histogram, stemplot, or dot chart (Fig.\ \@ref(fig:DataInsulationGraphPlots), right panel).


```{r DataInsulationGraphPlots, fig.cap="Two plots of the energy consumption data. The dotted line in the left panel shows the one home where energy consumption increased.", fig.align="center", fig.width=10, fig.height=3, out.width='95%'}
par( mfrow = c(1, 2))


WhichTime <- c(
  rep(1, length(InsulationBeforeAfter$Before)),
  rep(2, length(InsulationBeforeAfter$After))
)

plot( range(c(InsulationBeforeAfter$Before, InsulationBeforeAfter$After)) ~ c(1, 2),
      pch = 19,
      las = 1,
      ylim = c(8, 16),
      axes = FALSE,
      type = "n",
      main = "Energy consumption before and\nafter installing insulation",
      xlab = "Time",
      ylab = "Energy saving (in MWh)")
axis(side = 1, 
     at = 1:2,
     labels = c("Before", "After"))
axis(side = 2, 
     las = 1)
box()

for (i in 1:length(InsulationBeforeAfter$Before)){
  lines( c(1, 2),
         c(InsulationBeforeAfter$Before[i], InsulationBeforeAfter$After[i]),
         lty = ifelse( (InsulationBeforeAfter$Before[i] - InsulationBeforeAfter$After[i]) > 0, 1, 2)
#         col = grey(0.6)
  )
  points( c(1, 2),
         c(InsulationBeforeAfter$Before[i], InsulationBeforeAfter$After[i]),
         pch = 19
  )
}


###

stripchart( InsulationBeforeAfter$Diff,
	method = "jitter",
	xlab = "Energy saving (MWh)",
	las = 1,
	pch = 19,
	main = "",
	xlim = c(-2, 3),
	axes = FALSE,
	col = "blue")
axis(side = 1)
abline(v = 0, 
	col = "grey",
	lwd = 2)

text(-1, 1.30,
     labels = "Energy use\nincreased",
     col = "grey")
text(1.75, 1.30,
     labels = "Energy use\ndecreased",
     col = "grey")
```


<div style="float:right; width: 222x; border: 1px; padding:10px">
<img src="Illustrations/cake-916253_640.jpg" width="200px"/>
</div>



::: {.example #GraphPairedData name="Graphing paired data"}
A study [@data:Dawson2017:BirthdayCakes] examined the average number of bacteria on birthday cakes *before* and *after* blowing out the candles.
This RQ could be studied by taking two measurements from each cake: before and after blowing out candles.
The *change* in the number of bacteria could be computed for each cake, and a histogram of the differences plotted.
:::


### Time plots {#TimePlots}
   
Sometimes, a variable is measured over many time points.
A *time plot* can be used, which shows how the variable changes over time.


::: {.example #TimePlots name="Time plots"}
The baby-birth data (in Sect.\ \@ref(StemAndLeafPlots)) recorded the time of each birth.
A time plot shows how the weights varied over time 
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:BabyBoomTimeLATEX), left panel).'
} else {
   '(Fig.\\ \\@ref(fig:BabyBoomTimeHTML), top panel).'
}`
:::


```{r BabyBoomTimeLATEX, fig.cap="Time plots. Left: The weight of babies born on 21 December 1997 at a Brisbane hospital, plotted over time (encoded with gender); right: The number of lynx trapped in the Mackenzie River district of north-west Canada, from 1821 to 1934", fig.align="center", fig.height=4, fig.width=9, out.width='100%'}
data(Lynx)

if( knitr::is_latex_output() ) {
  
  par(mfrow = c(1, 2) )
    
  plot( Weight ~ Mins.Since.Midnight,
        data = BB,
        main = "Birthweight for babies born on\n21 December 1997",
        type = "b",
  #      pch = 19,
        ylim = c(1.5, 4.5),
        xlim = c(0, 1500),
        xlab = "Minutes after midnight",
        ylab = "Weight (in kg)",
        pch = ifelse(Gender == "Female", 19, 3),
        las = 1)
  legend("bottom", 
         pch = c(19, 3),
         legend = c("Female", "Male"))


  plot( Trappings ~ Year, 
        data = Lynx,
        type = "b",
        xlab = "Year",
        ylab = "Number of lynx trappings in Mackenzie River\nfrom 1821 to 1934",
        main = "Lynx trappings",
        pch = 19,
        xlim = c(1820, 1935),
        ylim = c(0, 7000),
        las = 1)
  
}
```




```{r BabyBoomTimeHTML, fig.cap="Time plots. Top: The weight of babies born on 21 December 1997 at a Brisbane hospital, plotted over time (encoded with gender); bottom: The number of lynx trapped in the Mackenzie River district of north-west Canada, from 1821 to 1934", fig.align="center", fig.height=4, fig.width=5}
if( knitr::is_html_output() ) {

  par(mfrow = c(2, 1) )
  
  plot( Weight ~ Mins.Since.Midnight,
      data = BB,
      main = "Birthweight for babies born on\n21 December 1997",
      type = "b",
#      pch = 19,
      ylim = c(1.5, 4.5),
      xlim = c(0, 1500),
      xlab = "Minutes after midnight",
      ylab = "Weight (in kg)",
      pch = ifelse(Gender == "Female", 19, 3),
      las = 1)
legend("bottom", 
       pch = c(19, 3),
       bty = "n",
       legend = c("Female", "Male"))


  plot( Trappings ~ Year, 
        data = Lynx,
        type = "b",
        xlab = "Year",
        ylab = "Number of lynx trappings",
        pch = 19,
        xlim = c(1820, 1935),
        ylim = c(0, 7000),
        las = 1)
  
}
```



::: {.example #TimePlots2 name="Time plots"}
The number of lynx trapped in the Mackenzie River area of Canada [@elton1942ten] each year from 1821 to 1934 is shown in 
`r if (knitr::is_latex_output()) {
   'Fig.\\ \\@ref(fig:BabyBoomTimeLATEX) (right panel).'
} else {
   'Fig.\\ \\@ref(fig:BabyBoomTimeHTML) (bottom panel).'
}`
A regular cycle may be present, when the number of trapping is very large.
:::



## Notes on constructing graphs {#GraphsConstructing}

Always remember:

::: {.importantBox .important data-latex="{iconmonstr-warning-8-240.png}"}
The purpose of a graph is to display the information in the simplest and clearest possible way, to help the reader understand the message(s) in the data.
:::


Helping readers to understand the data is the essence of producing a good graph.
You should be able to construct the graphs by hand, but usually software is used (e.g., jamovi or SPSS).
Graphs need to be clear and well-labelled with captions.

Ensure that you:

* *Do* add units of measurement and value labels where appropriate.
* *Do* add titles and axis labels.
* *Do* ensure the axis scales are appropriate.
* *Do* add any necessary explanations.
* *Do* make it easy for the reader to be see how the data answer the RQ (for example, to easily make the comparison of greatest interest).
* *Do not* add artificial third dimensions, or other 'chart junk' [@su2008s].
* *Do not* add optical illusions, such as an artificial third dimension.
* *Do not* make errors.


::: {.example #VerticalTruncation name="Truncating bar charts"}
One common optical illusion occurs when the frequency (or percentage) axis on a bar chart (or dot chart) does not start at zero.
For example, consider data recording the number of lung cancer cases in various Danish cities
`r if (knitr::is_latex_output()) {
   '(Fig.\\ \\@ref(fig:LCgraphsPDF)).'
} else {
   '(Fig.\\ \\@ref(fig:LCgraphsHTML)).'
}`

`r if (knitr::is_latex_output()) {
   'Figure\\ \\@ref(fig:BarchartTruncated1) (left panel) shows the original bar chart with the count (vertical) axis starting at zero; the counts in each age group look very similar. In contrast, if the vertical axis starts at 9.75, the counts look very different (Fig.\\ \\@ref(fig:BarchartTruncated1), right panel) for two age categories, suggesting large difference between the number of lung cancer cases. (The online version has an animation.)'
} else {
   'The animation below shows a bar chart with the count (vertical axis) starting at zero (when the counts in each age group look similar), and then gradually changing where the vertical axis starts... so that the final bar chart make the number of cases in each age groups look very different.'
}`
The graph is visually misleading when the graph does not start at a count of zero.
:::



```{r animation.hook='gifski', interval=0.5, dev=if (is_latex_output()){"pdf"}else{"png"}}
data(DanishLC)

ylim.lower.vec <- c( rep(0, 10), 
                     rep(5, 10),
		                 rep(9.75, 10) )
#                     seq(0, 9.75, by = 0.25) )
ylim.lower.vec <- c( ylim.lower.vec, 
                     rep( max(ylim.lower.vec), 10))

draw.bar <- function( ylim.lower) {
  egdata <- subset( DanishLC, 
                    City == "Fredericia")
  
  upperLimitYAxis <- 11 + (11 - ylim.lower) * 0.05 # Try to keep the gap between top bar and graph edge constant
  out <- barplot( egdata$Cases, 
                  las = 1,
                  col = plot.colour,
                  xlab = "Age groups",
                  ylab = "Number of cases",
                  main = paste("Number of lung cancer cases in Fredericia\n",
		                           ifelse(ylim.lower == 0,
		                                  "(Vertical axis NOT truncated)", 
				                              paste0("(Vertical axis truncated at ", ylim.lower,")") 
				      )
				),
                  xpd = FALSE, # BARS **DON'T** GO outside plotting region
                  ylim = c( ylim.lower, upperLimitYAxis),
                  axes = FALSE)
		  
  # ADD WHITE DOT to ensure many plots are created...
  for (i in (1:10)){
  points(2 + (i / 10), 
         mean( upperLimitYAxis, 11),
         pch = 19,
	 cex = 0.2,
	 col = "white")
  }
  axis(side = 1, 
       at = out[, 1], 
       labels = egdata$Age, 
       las = 2)
  axis(side = 2, 
       at = seq( floor(ylim.lower), 12, by = 1), 
       las = 1)
  box()
}

if (knitr::is_html_output()){
    for (i in 1:length(ylim.lower.vec) ){
      ylim.lower <- ylim.lower.vec[i]
      draw.bar( ylim.lower)
    }
  }
```



```{r BarchartTruncated1, fig.align="center", fig.width=10, fig.height=4, out.width='95%', fig.cap="The same bar chart, without truncating the vertical axis (left) and truncating the vertical axis (right)" }
if (knitr::is_latex_output()){
  len.ylv <- length(ylim.lower.vec)
  par( mfrow = c(1, 2))
  
  draw.bar( ylim.lower.vec[1] )
  draw.bar( ylim.lower.vec[ len.ylv ])
}
```


Graphs should focus on clear communication, not trying to appear fancy. 
One way people try to be fancy is to use an unnecessary third dimension.
This is poor: the resulting graphs can be misleading and hard to read [@siegrist1996use].


::: {.example #Bar2D3D name="Two- and three-dimensional plots"}
In the NHANES study [@data:NHANES3:Data], the age and sex of each participant were recorded (among other variables).
Using Fig.\ \@ref(fig:Bar3D2D) (left panel), can you easily determine if more females or more males are in each age group?

The artificial third dimension makes determining the heights of the bars hard.
In contrast, a 2-D graph (a side-by-side bar graph; Fig.\ \@ref(fig:Bar3D2D), right panel) makes it clear whether each age group has more females or more males.
:::


```{r Bar3D2D, out.width=c('50%', '46%'), fig.width=c(5, 4), fig.height=c(3.75, 4.5), fig.align="center", fig.cap="Two plots of the NHANES participants, divided by age group and sex. Left: A three-dimensional bar chart; right: a side-by-side bar chart", fig.show="hold"}
AgeD.Gender <- xtabs( ~ AgeD + NHANES$Gender)

TwoDData <- data.frame(
  Counts = c(
    AgeD.Gender[, 1],  
    AgeD.Gender[, 2]),
    Gender = c( rep("F", 7), 
                rep("M", 7)),
    Age = levels(AgeD)
  
)
TwoDData <- data.frame(
  Counts = c(
    AgeD.Gender[, 1],  
    AgeD.Gender[, 2]),
    Gender = c( rep("F", 7), 
                rep("M", 7)),
    Age = levels(AgeD)
  
)
TwoDData.mat <- array( dim = c(2, 7))
TwoDData.mat[1, ] <- TwoDData$Counts[1:7]
TwoDData.mat[2, ] <- TwoDData$Counts[8:14]
rownames(TwoDData.mat) <- c("Female", "Male")
colnames(TwoDData.mat) <- unique(TwoDData$Age)

col.mat <- TwoDData.mat
col.mat[1, ] <- rep(1, 7)
col.mat[2, ] <- rep(2, 7)

cols <- c('#FFFAF0', '#5844ff')

pmat <- plot3D::hist3D(z = TwoDData.mat, 
       x = c(0, 1),
       y = 1:7,
       col = cols,
       colvar = col.mat,
       border = "black",
       colkey = FALSE,
       facets = TRUE,
       xlab = "Sex",
       ylab = "Age group",
       zlab = "Number",
       #main = "Participants by group",
       space = 0.4,
       axes = FALSE,
       ticktype = "detailed",
       zlim = c(0, 900),
       ylim = c(0.5, 7.5),
       xlim = c(-0.5, 1.5),
       expand = 1,
       d = 25,
       mar = c(2, 2, 0, 0.5))


# Set vars
x.axis <- 0:1
min.x <- -0.5
max.x <- 1.5
y.axis <- 1:7
min.y <- 0.5
max.y <- 7.5
z.axis <- seq(0, 900, 
              by = 100)
min.z <- 0
max.z <- 900


# LINES
lines(trans3d(x = x.axis, 
              y = min.y, 
              z = min.z, 
              pmat) , 
      col = "black")
lines(trans3d(x = max.x, 
              y = y.axis, 
              z = min.z, 
              pmat) , 
      col = "black")
lines(trans3d(x = min.x, 
              y = min.y, 
              z = z.axis, 
              pmat) , 
      col = "black")



# TICKS
### See: http://entrenchant.blogspot.com/2014_03_01_archive.html
tick.start <- trans3d(x = x.axis, 
                      y = min.y, 
                      z = min.z, 
                      pmat)
tick.end <- trans3d(x = x.axis, 
                    y = (min.y - 0.25), 
                    z = min.z, 
                    pmat)
segments(x0 = tick.start$x, 
         y0 = tick.start$y, 
         x1 = tick.end$x, 
         y1 = tick.end$y)

tick.start <- trans3d(x = max.x,
                      y = y.axis, 
                      z = min.z, 
                      pmat)
tick.end <- trans3d(x = max.x + 0.20, 
                    y = y.axis, 
                    z = min.z, 
                    pmat)
segments(x0 = tick.start$x, 
         y0 = tick.start$y, 
         x1 = tick.end$x, 
         y1 = tick.end$y)

tick.start <- trans3d(x = min.x, 
                      y = min.y, 
                      z = z.axis, 
                      pmat)
tick.end <- trans3d(x = min.x, 
                    y = (min.y - 0.20), 
                    z = z.axis, 
                    pmat)
segments(x0 = tick.start$x, 
         y0 = tick.start$y, 
         x1 = tick.end$x, 
         y1 = tick.end$y)



# LABELS
### See: http://entrenchant.blogspot.com/2014_03_01_archive.html
labels <- c("F", 
            "M")
label.pos <- trans3d(x = x.axis, 
                     y = (min.y - 0.95),
                     z = min.z, 
                     pmat)
text(x = label.pos$x, 
     y = label.pos$y, 
     labels = labels, 
     adj = c(0, NA), 
     #srt = 270, 
     cex = 0.65)


labels <- levels(factor(TwoDData$Age))
label.pos <- trans3d(x = (max.x + 0.25), 
                     y = y.axis,
                     z = min.z, 
                     pmat)
text(x = label.pos$x, 
     y = label.pos$y, 
     labels = labels, 
     adj = c(0, NA), 
     srt = 0,
     cex = 0.65)

labels <- as.character(z.axis)
label.pos <- trans3d(x = min.x, 
                     y = (min.y - 0.5), 
                     z = z.axis, 
                     pmat)
text(x = label.pos$x, 
     y = label.pos$y, 
     labels = labels, 
     adj = c(1, NA), 
     cex = 0.65)


### AXIS NAMES
label.pos <- trans3d(x = mean(x.axis), 
                     y = (min.y - 1.25),
                     z = min.z, 
                     pmat)
text(x = label.pos$x, 
     y = label.pos$y, 
     labels = "Sex",
     cex = 0.75)

label.pos <- trans3d(x = (max.x + 0.80), 
                     y = 4.5,
                     z = min.z, 
                     pmat)
text(x = label.pos$x, 
     y = label.pos$y, 
     labels = "Age group", 
     cex = 0.75,
     srt = 42)

label.pos <- trans3d(x = min.x + 0.03, 
                     y = (min.y + 0.5), 
                     z = 650, 
                     pmat)
text(x = label.pos$x, 
     y = label.pos$y, 
     labels = "Number", 
     srt = 90,
     cex = 0.75)

title("Participants by group",
      cex.main = 0.9)

######################################################################################################

# Restore defaults
par( mar = c(5, 4, 4, 2) + 0.1 )

locatn <- barplot(t( AgeD.Gender ), 
        beside = TRUE,
        legend = FALSE, 
        las = 2, 
        cex.names = 0.8,
        ylim = c(0, 900),
        args.legend = list(x = "bottom", 
	                         bg = "white", 
			  cex = 0.75), 
        ylab = "Number",
        xlab = "Age group",
        main = "Participants by group",
        col = c("#FFFAF0", "#5844ff"))
#text(locatn[1, 4], 700, "F")
#text(locatn[2, 4], 700, "M")
legend("top", 
       legend = c("Females                 ", 
                  "Males"), 
       #levels(NHANES$Gender),
       bty = "n", 
       xpd = TRUE, 
       horiz = TRUE, 
       cex = 0.9, # For the *text* in the legend 
       fill = c( '#FFFAF0', '#5844ff') )

```





## Case study: the NHANES data {#NHANESGraphs}


<div style="float:right; width: 222x; border: 1px; padding:10px">
<img src="Illustrations/hamburger-2253349_640.jpg" width="200px"/>
</div>


To demonstrate how graphs help answer RQs, consider the following RQ:

> Among Americans, is the average direct HDL cholesterol concentration different for current smokers and non-smokers?

The Population is 'Americans', the Outcome is the 'average direct HDL cholesterol concentration', and the Comparison is 'Between those who currently do and do not smoke'.
There is no intervention.


::: {.thinkBox .think data-latex="{iconmonstr-light-bulb-2-240.png}"}
As with any study, possible confounding variables should be considered.\label{thinkBox:ConfoundingVars}
Can you think of any possible confounding variables?
`r if (!knitr::is_html_output()) '<!--'`
`r webexercises::hide()`
Maybe weight... Weight may be related to both the HDL cholesterol concentration *and* the smoking status.
`r webexercises::unhide()`
`r if (!knitr::is_html_output()) '-->'`
:::


Data to answer this RQs are obtained from the American *National Health and Nutrition Examination Survey* (NHANES) [@data:NHANES3:Data; @data:NHANES3].
From the NHANES webpage, this survey:

> ... examines a nationally *representative* sample of about 5,000 persons each year... located in counties across the country, 15 of which are visited each year.
>
> --- NHANES webpage (emphasis added): http://www.cdc.gov/nchs/nhanes/about_nhanes.htm

The data are equivalent to a "a simple random sample from the American population" [@data:NHANES:Rpackage].
In total, 10\,000 observations on scores of variables are available (from the 2009/2010 and the 2011/2012 surveys).
`r if (knitr::is_html_output()) {
   'Fig.\\ \\@ref(fig:NHANESFullDataTable) shows the first 1000 observations on the first 40 variables only.'
}`


```{r NHANESFullDataTable, fig.cap="The full NHANES dataset", fig.align="center"}
if( knitr::is_html_output() ) {
  DT::datatable(NHANES[1:1000, 1:40 ],
               list(scrollX = TRUE, 
	             scrollY = TRUE, 
	             ordering = FALSE),
               caption = "Part of the large NHANES dataset")
}
```


To answer any RQ, exploring and understanding the data are important. 
Graphs are an excellent way to do so (especially with large datasets).
Begin by graphing both the response and explanatory variables involved in the RQ.

For the NHANES data, the *response variable* is direct HDL cholesterol concentration (quantitative continuous).
The histogram (Fig.\ \@ref(fig:NHANESHDLHistBar), left panel) shows that the direct HDL cholesterol concentration is usually around 1.5mmol/L, and varies between about 0.5 and 3 mmol/L.
The data are slightly skewed right.


```{r NHANESHDLHistBar, fig.cap="Plots from the NHANES data. Left: histogram of the direct HDL cholesterol concetration (large outliers exist but are hard to see); centre: bar chart of current smoking status for the NHANES data; right: bar chart of current smoking status, including missing data", fig.align="center", fig.height = 3.5, fig.width=3.5, out.width=c('33%', '32%', '32%'), fig.show="hold"}

hist(NHANES$DirectChol,
     xlab = "Direct HDL cholesterol conc. (mmol/L)",
     ylab = "Frequency",
     main = "Distribution of direct\nHDL cholesterol concentration",
     las = 1,
     axes = FALSE,
     ylim = c(0, 2500),
     xlim = c(0, 4.5),
     breaks = seq(0, 4.5, 
                  by = 0.25),
     col = plot.colour)
axis(side = 1, 
     at = seq(0, 4.5, 
              by = 0.5))
axis(side = 2, 
     las = 1)
box()


barplot( table(NHANES$SmokeNow),
         xlab = "Does respondent currently smoke?",
         ylab = "Count",
         main = "Smoking status",
         ylim = c(0, 2000),
         las = 1,
         col = plot.colour)
box()

tabNA <- table(NHANES$SmokeNow, 
               exclude = NULL)
names(tabNA) <- c("No",
                  "Yes",
                  "Missing")

barplot( tabNA,
         xlab = "Does respondent currently smoke?",
         ylab = "Count",
         main = "Smoking status",
         ylim = c(0, 7000),
         las = 1,
         col = plot.colour)
box()
```


The *explanatory variable* is 'current smoking status', which could be graphed using a bar chart (Fig.\ \@ref(fig:NHANESHDLHistBar), centre panel) or pie chart, for example.
Since only *one* piece of information is present, a graph probably isn't necessary: `r round( table(NHANES$SmokeNow)["Yes"] / sum( table(NHANES$SmokeNow) ) * 100, 1)`% (just under half) of the respondents currently smoke.

Even though the NHANES data contains 10\,000 respondents, the bar chart clearly does not have 10\,000 responses: Many respondents did not answer this question.
A bar could be added to the bar chart to show the number of non-responses (though it probably isn't necessary): Fig.\ \@ref(fig:NHANESHDLHistBar), right panel.
  
The main RQ involves the *relationship* between average direct HDL cholesterol concentration and current smoking status, so a graph displaying this relationship is needed (Fig.\ \@ref(fig:NHANESboxplot)).
From the plot of the *sample* data, is there a difference in the average HDL cholesterol concentrations between the two smoking groups in the *population*?

    
```{r NHANESboxplot, fig.cap="Boxplot of direct HDL cholesterol concentration against current smoking status, for the NHANES data", fig.align="center", fig.height=4, fig.width=5}
boxplot(NHANES$DirectChol ~ NHANES$SmokeNow, 
        col = plot.colour, 
        las = 1,
        ylim = c(0, 4),
        ylab = "Direct HDL cholesterol conc. (mmol/L)",
        xlab = "Does respondent smoke now?",
        main = "Direct HDL cholesterol concentration by smoking\nstatus, for the NHANES data"
)
```


```{r animation.hook="gifski", dev=if (is_latex_output()){"pdf"}else{"png"}}
data(NHANES)

# Grab some of the necessary stats
Qs <- quantile(NHANES$Age, 
               seq(0, 1, by = 0.25)) # See below for better way when there are outliers

left <- 0.75
right <- 1.25
axis.edge <- axis.edge.left <- 0.5
axis.edge.right <- 1.5
angle <- 15
mid <- mean( c(left, 
               axis.edge.left) )
mid.r <- mean( c(right, 
                 axis.edge.right) )

# The box is at x=1, and goes from x=0.5 to 1.5 on the x-axis.
# See this by following the boxplot()) call with axis(side-1)

if (knitr::is_html_output()){
  for (i in 1:6){
    boxplot( NHANES$Age, 
             col = plot.colour,
             las = 1,
             ylim = c(0, 85),
             main = "Boxplot of Age (NHANES data)",
             ylab = "Age (in years)")
    
    if ( i == 2){
      arrows(left, Qs[5],
             axis.edge, 
	     Qs[5],
             angle = angle,
             lwd = 2)
      text(mid, Qs[5], 
           label = "Maximum", 
	   pos = 1)
    }
    if ( i == 3 ){
      arrows(left, Qs[4],
             axis.edge, 
	     Qs[4],
             angle = angle,
             lwd = 2)
      text(mid, Qs[4], 
           label = expression( italic(Q)[3]), 
	   pos = 1)
    }
    if (i == 4){
      arrows(left, 
             Qs[3],
             axis.edge, 
	     Qs[3],
             angle = angle,
             lwd = 2)
      text(mid, 
           Qs[3], 
           label = expression( italic(Q)[2]), 
	   pos = 1)
    }
    if (i == 5){  
      arrows(left, 
             Qs[2],
             axis.edge, 
	     Qs[2],
             angle = angle,
             lwd = 2)
      text(mid, 
           Qs[2], 
	   label = expression( italic(Q)[1]), 
	   pos = 1)
    }
    if ( i == 6){    
      arrows(left, 
             Qs[1],
             axis.edge, 
	     Qs[1],
             angle = angle,
             lwd = 2)
      text(mid, 
           Qs[1], 
	   label = "Minimum", 
	   pos = 3)
    }
  }
}
```


The NHANES study is an observational study, so confounding is a potential problem (Sect.\ \@ref(ExtraneousVariables)), so exploring the relationships involving extraneous variables is useful (Fig.\ \@ref(fig:NHANESplots)).


```{r NHANESplots, fig.cap="Some plots from the NHANES data. Top panels: the response variables and some extraneous variables; bottom panels: the explanatory variable and some extraneous variables", fig.align="center", fig.height=5, fig.width=6.5, out.width='100%'}
par( mfrow = c(2, 2),
     mar = c(5, 5.5, 4, 2) + 0.1 ) # DEFAULT: c(5, 4, 4, 2) + 0.1.

plot( DirectChol ~ Weight, 
      data = NHANES, 
	    xlab = "Weight (in kg)", 
	    ylab = "Direct HDL cholesterol\nconc. (mmol/L)",
      main = "NHANES: HDL cholesterol\nconc. vs Weight",
	    ylim = c(0, 4),
      col = blueTransparent,
	    las = 1)


plot( DirectChol ~ Gender, 
      data = NHANES, 
	    xlab = "Gender", 
	    ylab = "Direct HDL cholesterol\nconc. (mmol/L)",
      main = "NHANES: HDL cholesterol\nconc. vs Gender",
	    names = c("Female", "Male"),
      col = plot.colour,
	    ylim = c(0, 4),
	    las = 1)

plot( Weight ~ SmokeNow, 
      data = NHANES, 
	    xlab = "Current smoker?", 
	    ylab = "Weight (in kg)",
      main = "NHANES: Weight\nvs smoking status",
	    horiz = TRUE,
      col = plot.colour,
	    ylim = c(0, 250),
	    las = 1)

par(mar = c(5, 4, 4, 2) + 0.1)  # DEFAULT: c(5, 4, 4, 2) + 0.1.

out <- barplot( prop.table( xtabs(~ SmokeNow + Gender, 
                                  data = NHANES),  
				  margin = 2) * 100, 
	ylim = c(0, 100),
	xlab = "Gender", 
	ylab = "Percentage of people",
  main = "NHANES: Smokers\nby Gender",
	names = c("Female", 
	          "Male"),
	col = c("white", 
	        plot.colour),
	las = 1)
text( out[2], 25, 
      "Smokers"  	 )
text( out[2], 75, 
      "Non-\nsmokers")
```


::: {.thinkBox .think data-latex="{iconmonstr-light-bulb-2-240.png}"}
What do these graphs suggest about possible confounding relationships?
:::



## Summary {#Chap12-Summary}

For *one quantitative variable*, possible graphs include a histogram (for medium to large amounts of data), a dot chart (for small to moderate amounts of data), or a stemplot (for small amounts of data, and only useful in some circumstances).

For *one qualitative variable*, possible graphs include a bar chart (usually a good choice), a dot chart (usually a good choice), or a pie chart (only useful in special circumstances, and can be harder to interpret).

For two variables, the graphs shown in Table\ \@ref(tab:GraphsToUse) can be used.
In general, the 'response' variable is shown on the vertical axis.


```{r GraphsToUse}
GraphTab <- array( dim = c(2, 2) )
colnames(GraphTab) <- c( "Explanatory: Qualitative", 
                         "Explanatory: Quantitative")
rownames(GraphTab) <- c("Response: Qualitative",
                        "Response: Quantitative")

GraphTab[1, ] <- c("Stacked bar; side-by-side bar", 
                   "Boxplot; back-to-back stem; 2-D dot")
GraphTab[2, ] <- c("Boxplot; back-to-back stem; 2-D dot", 
                   "Scatterplot")

if( knitr::is_latex_output() ) {
    kable(GraphTab,
        format = "latex",
        longtable = FALSE,
        booktabs = TRUE,
        escape = FALSE, # For latex to work in \rightarrow
        caption = "Graphs to use for most situations",
        linesep = "\\addlinespace",
        align = c("l", "l"))   %>%
   kable_styling(full_width = FALSE, font_size = 10) %>%
  row_spec(0, bold = TRUE) %>% # Columns headings (row 0) in bold
  column_spec(1, bold = TRUE) %>%     # Row headings in bold
  column_spec(1, width = "20mm") %>%
  column_spec(2, width = "40mm") %>%
  column_spec(3, width = "40mm")
}

if( knitr::is_html_output() ) {
    kable(GraphTab,
        format = "html",
        align = c("l", "l"),
        longtable = FALSE,
        caption = "Graphs to use for most situations",
        booktabs = TRUE)
}
```



::: {.thinkBox .think data-latex="{iconmonstr-light-bulb-2-240.png}"}
In the following scenarios, select the most appropriate graph to use.
\vspace{-2ex}

1. To compare the average number of hours of sunlight exposure per day for female and male teachers. \tightlist
`r if( knitr::is_html_output()) {longmcq( c(
   "A histogram of the differences", 
   answer = "A boxplot",
   "A stacked bar chart",
   "A histogram") )}`
2. To compare the proportion of trees with koalas, for trees more than 10 metres tall with trees 10 metres or shorter.
`r if( knitr::is_html_output()) {longmcq( c(
   "A pie chart", 
   answer = "A stacked bar chart",
   "A boxplot",
   "A histogram of the differences") )}`
3. To compare the average number of hours spent on social media for people aged over 30, people aged 30 to under 50, and people over 50.
`r if( knitr::is_html_output()) {longmcq( c(
   "A stacked bar chart",
   "A pie chart", 
   answer = "A boxplot",
   "A histogram",
   "A histogram of the differences") )}`
4. To compare the average number of hours of exercise per week for office workers, both in summer and in winter for each person.
`r if( knitr::is_html_output()) {longmcq( c(
   "A side-by-side bar chart",
   "A boxplot",
   "A barchart",
   answer = "A histogram of the differences for each person") )}`
:::



`r if (knitr::is_html_output()){
  'The following short video may help explain some of these concepts:'
}`

<div style="text-align:center;">
<iframe width="560" height="315" src="https://www.youtube.com/embed/FMprCwua-ew" frameborder="0" allow="accelerometer; encrypted-media; gyroscope; picture-in-picture"></iframe>
<iframe width="560" height="315" src="https://www.youtube.com/embed/FMprCwua-ew" frameborder="0" allow="accelerometer; encrypted-media; gyroscope; picture-in-picture"></iframe>
</div>


<!-- <iframe src="https://www.oppia.org/embed/exploration/hGescxPJnaHe" width="700" height="1000"> -->


## Quick review questions {#Chap12-QuickReview}

::: {.webex-check .webex-box}
A study on the bruising of apples [@doosti2016development] purposefully hit apples with three different forces (200, 700 and 1200 mJ), in three different regions of the apple (lower; middle; upper). 
The researchers then recorded the depth of the bruising (the *response* variable), to determine if it could be estimated from the surface temperature (the *explanatory* variable) at the bruise location (an *extraneous* variable).

1. What graph would be appropriate to display the *response* variable?\tightlist  
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	answer = "A histogram",
	"A boxplot",
	"A scatterplot",
	"A side-by-side bar chart"))}`

2. What graph would be appropriate to display the *explanatory* variable?  
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	answer = "A histogram",
	"A boxplot",
	"A scatterplot",
	"A side-by-side bar chart"))}`

3. What graph could be used to display the *relationship* between the response and explanatory variables?  
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	"A histogram",
	"A boxplot",
	answer = "A scatterplot",
	"A side-by-side bar chart"))}`

4. What graph could be used to display the *relationship* between the 'surface temperature' and 'bruise location'?  
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	"A histogram",
	answer = "A boxplot",
	"A scatterplot",
	"A side-by-side bar chart"))}`
:::


`r if (!knitr::is_html_output()) '<!--'`
`r webexercises::hide()`
1. **Histogram**: One quantitative variable.
2. **Histogram**: One quantitative variable.
3. **Scatterplot**: Two quantitative variables.
4. **Boxplot**: One quantitative and one qualitative variable.
`r webexercises::unhide()`
`r if (!knitr::is_html_output()) '-->'`



## Exercises {#GraphsExercises}

Selected answers are available in Sect.\ \@ref(GraphsAnswer).


::: {.exercise #Graphs123}
What would be the best graph for displaying the data for these situations?

1. To explore the relationship between the pH of water and the temperature of the water, in various creeks.
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	"A histogram",
	"A boxplot",
	"A histogram of the differences",
	answer = "A scatterplot",
	"A side-by-side bar chart"))}`
2. Suppose, in a research study, the researchers measures the difference between a swimmer's fastest 100m time and their fastest 200m time.
The researchers were interested in the average time *difference*.
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	"A histogram",
	"A boxplot",
	answer = "A histogram of the differences",
	"A scatterplot",
	"A side-by-side bar chart"))}`
3. A research study examined the way in which students usually came to university (bus; private car; car pooling; etc.) and their postcode.
`r if( knitr::is_html_output() ) {longmcq( c(
	"A bar chart",
	"A histogram",
	"A boxplot",
	"A histogram of the differences",
	"A scatterplot",
	answer = "A side-by-side bar chart"))}`
:::


::: {.exercise #GraphsCars}
A study [@henderson1981building] recorded the number of cylinders in many models of cars: eleven cars had four cylinders, seven had six cylinders, and fourteen had eight cylinders.
The *number* of cylinders is quantitative discrete, but with so few different values, this data could be plotted with some of the graphs used for qualitative data.
For these data:

1. Produce a dot chart.
1. Produce a histogram.
1. Produce a bar chart.
1. Produce a pie chart.

What graph do you think is best?
Why?
:::


::: {.exercise #GraphsLimeTrees}
A study of lime trees (*Tilia cordata*) recorded these variables for 385 lime trees in Russia [@data:ForestBiomass2017]:

* the foliage biomass, in kg (the response variable);
* the tree diameter, in cm (the explanatory variable);
* the age of the tree, in years; and
* the origin of the tree, one of Coppice, Natural, or Planted.

The purpose of the study is to estimate the foliage biomass from the tree diameter, and perhaps the other extraneous variables.
What graphs would be useful?
:::


::: {.exercise #GraphsOrthoses}
In a study of the influence of using ankle-foot orthoses in children with cerebral palsy [@data:Swinnen2017:orthoses], the data in Table\ \@ref(tab:DescribeAnkleFoot) describe the 15 subjects.
(GMFCS is an *ordinal* variable used to describe the impact of cerebral palsy on their motor function: the
`r if (knitr::is_latex_output()) {
   'Gross Motor Function Classification System.)'
} else {
   '[Gross Motor Function Classification System](https://en.wikipedia.org/wiki/Gross_Motor_Function_Classification_System).)'
}`
Sketch some graphs to explore the *relationships* between these variables.
:::


::: {.exercise #GraphNitrogenInSoil}
A study of fertilizer use [@data:Lane2002:GLMsoilscience] recorded the soil nitrogen after applying different fertilizer doses.
The researchers recorded:

* the fertilizer dose, in kilograms of nitrogen per hectare; 
* the soil nitrogen, in kilograms of nitrogen per hectare; and
* the fertilizer source; one of 'inorganic' or 'organic'.

What graphs would be useful for understanding the data?
:::


::: {.exercise #GraphSurveyData}
A 
`r if (knitr::is_latex_output()) {
   'survey of voice assistants'
} else {
   '[survey of voice assistants](https://www.nielsen.com/us/en/insights/news/2018/smart-speaking-my-language-despite-their-vast-capabilities-smart-speakers-all-about-the-music.print.html)'
}`
(such as Amazon Echo; Google Home; etc.) conducted by
`r if (knitr::is_latex_output()) {
   'Nielsen'
} else {
   '[Nielsen](https://www.nielsen.com/au/en.html)'
}`
asked respondents to indicate how they used their voice assistant.
The options were:

* Listening to music;
* Search for real-time info (e.g., traffic; weather);
* Search for factual info (e.g., trivia; history);
* Listen to news;
* Chat with voice assistant for fun;
* Use alarms, timer.

What would be the best graph for displaying respondents answers?
Would a pie chart be suitable? Explain your answer.
:::


::: {.exercise #GraphsAIS}
A study of athletes at the *Australian Institute of Sport* (AIS) measured numerous physical and blood measurements from high performance athletes [@data:Telford1991:sexsportsize].
The graph in Fig.\ \@ref(fig:AISfemalesport) compares the heights of females in two similar sports: basketball and netball.
(Netball was derived from basketball.)
How would you compare the heights of the athletes in the two sports?
:::

```{r}
data(AISsub)

subAIS <- subset(AISsub, 
                 Sex=="F" & (Sport=="BBall" | Sport=="Netball"))
subAIS$Sport <- factor(subAIS$Sport)
```

```{r AISfemalesport, fig.cap="The heights of female basketball and netball players attending the AIS", fig.align="center", fig.width=5, fig.height=3.5}
boxplot(Ht ~ Sport, 
        data = subAIS,
        las = 1,
        ylim = c(165, 200),
        col = plot.colour,
        xlab = "Sport",
        ylab = "Height (in cm)",
        main = "Boxplot of heights:\nbasketballers and netballers",
        names = c("Basketball", 
                  "Netball")
)
```


::: {.exercise #GraphNoisyMiners}
A study of noisy miners (a small Australian bird) counted the number of noisy miners and the number of eucalyptus trees in random quadrats [@data:Maron:eucthreshold; @mypapers:dunnsmyth:glms].
Critique the graph of the data (Fig.\ \@ref(fig:MinerPlot)).
:::


```{r, MinerPlot, fig.cap="The number of noisy miners and the number of eucalyptus trees", fig.align="center", fig.width=5, fig.height=3.5}
data(NMiner)

plot(Minerab ~ Eucs,
     data = NMiner,
     ylim = c(0, 35),
     xlab = "E_trees",
     ylab = "num_miners",
     pch = 1:9)
```


::: {.exercise #GraphHorseshoeCrabs}
A study of 173 female horseshoe crabs [@data:brockmann:crabs; @mypapers:dunnsmyth:glms] recorded, among other things, the colour of the carapace ('Light medium', 'Medium', 'Dark medium' or 'Dark') and the condition of the carapace  ('Both OK', 'One OK', 'None OK').
Critique the scatterplot (Fig.\ \@ref(fig:CrabPlot)) used to explore the data.
:::



```{r CrabPlot, fig.cap="A scatterplot of the colour of female horseshoe crabs and the condition of their spines. There are no missing values.", fig.align="center", fig.width=5, fig.height=3.5}
data(HCrabs)

plot( as.numeric( factor(Col) ) ~ as.numeric( factor(Spine) ), 
      data = HCrabs,
      xlab = "Condition of the spine",
      ylab = "Carapace colour",
      main = "Dot chart of carapace colouyr\nagainst spine condition",
      las = 1)
```



::: {.exercise #GraphsMADRS}
A study [@data:Danielsson2014:Depression] examined the change in MADRS (a *quantitative* scale measuring level of depression) and treatment group (whether each person was treated using: exercise; body awareness; or advice).

1. What is the response variable?
1. What is the explanatory variable?
1. What graphs would be useful for exploring the data and the relationships of interest?
:::


::: {.exercise #GreenBuilding}
In a study of the temperature in offices, @data:Paul2008:Comfort compared the temperature in three offices (during working hours) at Charles Sturt University (Australia); the data are summarised in Table\ \@ref(tab:OfficeTemps).
Using this information, draw the boxplot comparing the three offices.
What do we learn from this graph? 
:::


```{r OfficeTemps}
OfficeTemps <- array( dim = c(6, 3) )
colnames(OfficeTemps) <- c( "Office A", 
                            "Office B", 
                            "Office C")
rownames(OfficeTemps) <- c("Mean", 
                           "Minimum", 
                           "$Q_1$", 
                           "Median", 
                           "$Q_3$", 
                           "Maximum")

OfficeTemps[1, ] <- c(24.1, 25.3, 25.7)
OfficeTemps[2, ] <- c(16.4, 15.9, 20.1)
OfficeTemps[3, ] <- c(22.8, 23.8, 24.6)
OfficeTemps[4, ] <- c(24.4, 25.5, 26.1)
OfficeTemps[5, ] <- c(25.5, 26.9, 27.2)
OfficeTemps[6, ] <- c(27.4, 31.0, 30.3)

if( knitr::is_latex_output() ) {
    kable(OfficeTemps,
        format = "latex",
        longtable = FALSE,
        booktabs = TRUE,
        escape = FALSE, # For latex to work
        linesep = c( "", "", "", "", "", ""), # Add a bit of space between all rows. 
        caption = "A summary of the temperature (in degrees C) in three offices at CSU during working hours
according to current smoking status",
        align = c("c", "c"))   %>%
   kable_styling(full_width = FALSE, 
                 font_size = 10) %>%
  row_spec(0, bold = TRUE) %>% # Columns headings (row 0) in bold
  column_spec(0, bold = TRUE)     # Row headings in bold
}

if( knitr::is_html_output() ) {
  kable(OfficeTemps,
        format = "html",
        align = c("c", "c"),
        longtable = FALSE,
        caption = "A summary of the temperature (in degrees C) in three offices at CSU during working hours according to current smoking status",
        booktabs = TRUE) 
}
```


::: {.exercise #GraphsSkewBar}
A study of high-performance athletes at the *Australian Institute of Sport* (AIS) [@data:Telford1991:sexsportsize] recorded numerous variables about athletes.
A plot for the sports played by the athletes is shown in Fig.\ \@ref(fig:AISSportBarchart).

How would you describe the data: Left skewed, right skewed, approximately symmetrical?
Or something else?
:::


```{r AISSportBarchart, fig.cap="Sports played by athletes in the AIS study", fig.align="center", fig.width=5, fig.height=4}
data(AISsub) 

SPT <- sort( table(AISsub$Sport) ) # Sorted from smallest to largest category

par(mar = c(6, 4, 4, 2) + 0.1 ) # DEFAULT: c(5, 4, 4, 2) + 0.1.
barplot(SPT, 
        las = 2,
        ylab = "Number of AIS athletes",
        xlab = "Sport",
        main = "Sports played by AIS athletes",
        col = plot.colour)
```


<!-- QUICK REVIEW ANSWERS -->
`r if (knitr::is_html_output()) '<!--'`
::: {.EOCanswerBox .EOCanswer data-latex="{iconmonstr-check-mark-14-240.png}"}
**Answers to in-chapter questions:**

- Sect. \ref{thinkBox:DescribeBrainFreeze}:
*Average*: Hard to be sure... maybe between 10 or 15. (More observations appear at the smaller values (as the bars are higher).) 
*Variation*: From about 0 to about 40. 
*Shape*: Slightly skewed right. 
*Outliers*: No outliers or unusual observations. The observation between 35 and 40 *may* be an outlier.
I suspect it is *not* an outlier, as a larger sample may very well have observations between 30 and 35.
Every sample is a bit different, of course.
Of course, I could be wrong.
- Sect. \ref{thinkBox:WhenPiechart}:
**1.** A pie chart is **not suitable**. 
The three browsers are not *mutually exclusive* (people can use more than one) nor *exhaustive* (some people may use other browsers). 
**2.** A pie chart is **not suitable**: the percentages are not parts of a whole.
Each individual has information recorded on *two* qualitative variables: (a) state of residence; and (b) whether or not they own an iPhone.

- Sect. \ref{thinkBox:ConfoundingVars}: Maybe weight... Weight may be related to both the HDL cholesterol concentration *and* the smoking status.

- \textbf{\textit{Quick Revision} questions:}
**1.** Histogram.
**2.** Histogram.
**3.** Scatterplot.
**4.** Boxplot.
:::
`r if (knitr::is_html_output()) '-->'`