---
title: "Summary statistics"
author: "Jonathan Harrington, Johanna Cronenberg"
output: 
  bookdown::html_document2:
    number_sections: TRUE
    toc: true
    theme: flatly
    highlight: pygments
---

<style>
div.col {background-color: #ffd5c0; border-radius: 5px; padding: 20px;}
body {font-size: 16pt;}
h1 {font-size: 24pt;}
h2 {font-size: 22pt;}
p.caption {font-size: 12pt; text-align: justify;}
code.sourceCode {font-size: 16pt;}
</style>

If not already done, carry out parts 1-5 of the setup,
[as described here](https://www.phonetik.uni-muenchen.de/studium_lehre/lehrmaterialien/R_speech_processing/setup.html).

Load this library:

```{r, message=FALSE}
library(tidyverse)
```

Load these objects:

```{r}
urla = "https://www.phonetik.uni-muenchen.de/studium_lehre/"
urlb = "lehrmaterialien/R_speech_processing/Rdf"
url = paste0(urla, urlb)
int = read.table(file.path(url, "intdauer.txt"), 
                 stringsAsFactors = T)
vdata = read.table(file.path(url, "vdata.txt"), 
                   stringsAsFactors = T)
```

The module is about using summary statistics that are almost always applied before doing any statistical test. Summary statistics are used to provide a general overview of the data.  Different types of summary statistics are applied to *numeric* and *categorical* data. The important differences between these two types of data 
[should be reviewed here](https://www.phonetik.uni-muenchen.de/studium_lehre/lehrmaterialien/R_speech_processing/01_Rintro.html#numeric-and-categorical-variables) before continuing with this module.

# Summary statistics for numeric variables

Most of the functions for summary statistics are evident in the output of the `summary()` function, e.g. as applied to the data-frame `vdata`:

```{r}
summary(vdata)
```

Some further details are as follows:

```{r}
# create a vector of values
vec = c(10, 12, -8, 4, 2, 0, 19, -3, 18)
# mean 
mean(vec)      
# variance
var(vec) 
# standard deviation
sd(vec)       
# minimum
min(vec)       
# maximum
max(vec)     
# range
range(vec)
# first quartile (25% quantile)
quantile(vec, 0.25)
# third quartile
quantile(vec, 0.75) 
# median or 2nd quartile
median(vec)         
# as above
quantile(vec, 0.5)
# interquartile range
IQR(vec)
# as above
quantile(vec, 0.75) - quantile(vec, 0.25) 
```

<div class="col">

1. The mean is the sum of values divided by the number of values:

```{r}
sum(vec)/length(vec)
# the same
mean(vec)
```

2. The median or 0.5 quantile is the middle ordered value *after the values have been sorted*. The middle ordered value is (`N`-1)/2 + 1 where `N` is the length of the vector. The length of the vector in this example is 9:

```{r}
length(vec)
```

so the middle ordered value is ((9-1) * 0.5) + 1 = 5. Now sort the vector in ascending order and find the 5th value:

```{r}
vec.sort = sort(vec)
vec.sort[5]
# the same
median(vec)
```

The median is less sensitive to outliers than the median. For example:

```{r}
# make another vector, the same as `vec`
vec2 = vec
# make the first value of `vec2` an outlier
vec2[1] = 20000
vec
vec2
# mean of vec and vec2: big difference
mean(vec)
mean(vec2)
# median of vec and vec2: the same
median(vec)
median(vec2)
```

3. The 25% or 0.25 quantile is obtained in the same way as the median, but this time multiply by 0.25 to find which sorted value is needed: ((9-1) * 0.25) + 1 = 3. So the 0.25 quantile is:

```{r}
vec.sort
vec.sort[3]
# the same
quantile(vec, 0.25)
```

4. The variance is a measure of how dispersed the values are.

```{r}
vec
# this has a smaller variance than thos in `vec` 
# because they cluster more tightly around the mean
vec2 = c(10, 12, -1, 4, 2, 0, 10, -1, 11)
var(vec)
var(vec2)

# the variance of this must be zero
var(rep(1, 10))
```

5. The standard deviation is the square root of the variance

```{r}
sqrt(var(vec))
# or
(var(vec))^0.5
# the same:
sd(vec)
```
</div>

There are numerous ways of applying any of the above functions to numeric variables of data-frames. Here are 4 of them used to calculate the mean of the numeric variable `F1` in the data-frame `vdata`:

```{r}
mean(vdata$F1)
with(vdata, mean(F1))
# pull(F1) extracts the F1 values from the data-frame
vdata %>% pull(F1) %>% mean()
# this makes another data-frame with the variable F1mean
vdata %>% summarise(F1mean = mean(F1))
```

By far the most useful of these for plotting and applying statistical functions will be the last involving `summarise()` because it can be applied to several variables at once:

```{r}
# mean of F1 and mean of F2. These summary variables
# are stored in mF1 and mF2; but you can choose other names
vdata %>%
  summarise(mF1 = mean(F1), mF2 = mean(F2))
```

```{r}
# as above but include the standard deviations and interquartile ranges
vdata %>%
  summarise(mF1 = mean(F1), 
            mF2 = mean(F2),
            sF1 = sd(F1),
            sF2 = sd(F2),
            iF1 = IQR(F1), 
            iF2 = IQR(F2))
```

<div class="col">
The function `mutate()` 
[that was discussed in an earlier module](https://www.phonetik.uni-muenchen.de/studium_lehre/lehrmaterialien/R_speech_processing/02_Rintro.html#adding-extra-variables-to-data-frames-mutate)
produces the same numerical output as `summarise()`. The main difference is that whereas the output from `summarise()` is in a *separate* data-frame consisting of just the newly created variable, `mutate()` attaches the new variable to an *existing* data-frame. For example, the earlier command: 

```{r}
vdata %>%
  summarise(mF1 = mean(F1), mF2 = mean(F2))
```

could also be accomplished with:

```{r}
mdata = vdata %>%
  mutate(mF1 = mean(F1), mF2 = mean(F2))
head(mdata)
```

from which it becomes clear that the mean F1 and mean F2 values have been attached to every observation of the data-frame `vdata`. In general, **`summarise()` is the preferred function when (as the name suggests) providing summary statistics across several observations** i.e. when applying summary functions like `mean()`, `median()`, `var()` and so on. `mutate()` **should instead be used when a separately calculated value is needed for each observation**. 

Consider in this regard a common transformation in speech science of converting numeric (especially formant) data into so-called *z*-scores which show how many standard deviations a value is from the mean of the distribution. Applied to the earlier example of the vector of numbers, this is given by the following:

```{r}
(vec - mean(vec))/sd(vec)
```

The same can be applied to a data-frame using `mutate()`. In this case, *z*-scores are calculated for `F1` and the output is additionally stored as a new data frame:

```{r}
# store the output in the new data-frame `zdata`
zdata = 
  vdata %>%
  mutate(F1z  = (F1 - mean(F1)/sd(F1)))
```

The above *could* also be done with `summarise()` but  `mutate()` is for the reason stated above the preferred function in this case. Finally, the above can also be accomplished by writing a function then passing that function to `mutate()`. Here is a function called `zscore()` for computing *z*-scores from any vector:

```{r}
zscore = function(x) {
  (x - mean(x))/sd(x)
}
```

Here is the *z*-score for `vec`:

```{r}
vec
zscore(vec)
```

Here are the *z*-scores for F1 and for F2 of `vdata`:

```{r}
zdata = 
  vdata %>%
  mutate(zF1 = zscore(F1), zF2 = zscore(F2))
```
</div>

# Summary statistics for categorical variables (factors)

The simplest way to obtain this kind of information for a categorical variable in a data-frame is with the functions `group_by()` and `n()` inside `summarise()`. For example, to make a table to show how often the separate levels of the factor `V` in the data-frame `vdata` occur:

```{r}
vdata %>%
  group_by(V) %>%
  summarise(count = n()) %>%
  ungroup()

# do the same but store it as a new data-frame:
Vcount = 
  vdata %>%
  group_by(V) %>%
  summarise(count = n()) %>%
  ungroup()
Vcount
```

from which it becomes immediately apparent that:

- the factor `V` is made up of 7 levels
- the levels are "%" (=/ø/ as in German 'S**öh**ne'), "A", "E", "I", "O", "U", "Y"
- they are each over 400 times in the data-frame.

<div class="col">
The output of the above code with `group_by()` returns an object called a *tibble* as shown here:

```{r}
class(Vcount)
```

A tibble is a simplified version of a data-frame and the difference between the two will be largely unimportant for the subsequent statistical analyses to be considered.

Conversion back to a data-frame is possible as follows:

```{r}
Vcount.df = Vcount %>% as.data.frame()
class(Vcount.df)
```

The function `as_tibble()` can be used to convert a data-frame to a tibble. More information about tibbles is in
[chapter 10, *R for Data Science*](https://r4ds.had.co.nz/tibbles.html).

It's always a good idea to use the function `ungroup()` after using the `group_by()` function because otherwise R makes assumptions about the grouping of variables in a data-frame that can influence subsequent calculations in unexpected ways. Using `ungroup()` in the above example was in fact unnecessary because there was only one variable `V` but when there is more than one variable (as in the example below), `ungroup()` should always be used.

</div>

A table of levels and their frequency of occurrence for the *combination* of any factors can be made by passing the names of the factors to the `group_by()` function together, thus:

```{r}
vdata %>%
  group_by(V, Cons) %>%
  summarise(count = n()) %>%
  ungroup()
```

For example, the fourth line means that there are 147 observations which have  level `A` (of factor `V`) and level `K` (of factor `Cons`). The first few observations above show that there are three consonant levels per vowel level. The function `n_distinct()` can be used to count how many levels of one factor occur relative to any other. For example:

```{r}
vdata %>% 
  group_by(V) %>% 
  summarise(count = n_distinct(Cons)) %>%
  ungroup()
```

means that for each of the 7 vowel levels in factor `V` there are three consonant levels (of factor `Cons`). The following:

```{r}
vdata %>% 
  group_by(V, Cons) %>% 
  summarise(count = n_distinct(Subj)) %>% 
  ungroup()
```

means that for each combination of consonant and vowel levels, there are 7 subject levels.

# Summary statistics for numeric variables grouped by factor.

The `group_by()` function can be also used to calculate summary statistics for numeric variables separately for the levels of a factor. For example,  instead of calculating mean F1 for the entire data-frame, separate F1 means might instead be calculated for each of the 7 vowel levels in `V`. This is done as follows:

```{r}
vdata %>%
  # for each level in `V`
  group_by(V) %>%
  # calculate the mean of F1
  summarise(mF1 = mean(F1)) %>%
  ungroup()
```

The above can be extended to include any number of factors. Thus to calculate mean F1 for each level when combining consonant and vowel factors:

```{r}
cvmeans = 
  vdata %>%
  group_by(Cons, V) %>%
  summarise(mF1 = mean(F1)) %>%
  ungroup()
cvmeans
```

The meaning of the 4th line in the above is:  the mean of F1 in /KA/ (i.e. the mean of observations for which `Cons` is `K` and for which `V` is `A`) is 307 Hz.

