Lecture 4

class: center, middle, inverse, title-slide

# Lecture 4
## Getting basic statistics / Probability distributions
### Jung-Jin Lee
### Jan 28, 2020

---

## group_by() / summarize()

**Example**: find the mean `Sepal.Length` for each `Species`.

```r
setosa_df <- iris %>% 
  filter(Species == "setosa")
mean(setosa_df$Sepal.Length)
```

```
## [1] 5.006
```

```r
versi_df <- iris %>% 
  filter(Species == "versicolor")
mean(versi_df$Sepal.Length)
```

```
## [1] 5.936
```

```r
virgi_df <- iris %>% 
  filter(Species == "virginica")
mean(virgi_df$Sepal.Length)
```

```
## [1] 6.588
```

---

## group_by() / summarize(), continued

`group_by()` and `summarize()` combination: summary of data (e.g. mean, max, min, median) by each category of a categorical variable.

**Example**: find the mean `Sepal.Length` for each `Species`.

```r
# For each 'Species' category in the 'iris' data set
# compute the mean of 'Sepal.Length'
*group_by(iris, Species) %>%
* summarize(mean.Sepal.Length = mean(Sepal.Length))
```

```
## # A tibble: 3 x 2
##   Species    mean.Sepal.Length
##   <fct>                  <dbl>
## 1 setosa                  5.01
## 2 versicolor              5.94
## 3 virginica               6.59
```

Note that the result is a data frame with two columns: `Species` and `mean.Sepal.Length`.

---

## group_by() / summarize(), continued

**Exercise**: What is the maximum (and minimum) of `Sepal.Length` in each `Species`?

```r
iris %>% group_by(Species) %>%
* summarize(max.Sepal.Length = max(Sepal.Length), # max(): maximum
*           min.Sepal.Length = min(Sepal.Length)) # min(): minimum
```

```
## # A tibble: 3 x 3
##   Species    max.Sepal.Length min.Sepal.Length
##   <fct>                 <dbl>            <dbl>
## 1 setosa                  5.8              4.3
## 2 versicolor              7                4.9
## 3 virginica               7.9              4.9
```

---

## group_by() / summarize(), continued

One can get summary data using _multiple_ "group by" variables.

```r
diamonds %>% 
* group_by(color, cut) %>%
  summarize(avg.price = mean(price), top.price = max(price)) 
```

```
## # A tibble: 35 x 4
## # Groups:   color [7]
##    color cut       avg.price top.price
##    <ord> <ord>         <dbl>     <int>
##  1 D     Fair          4291.     16386
##  2 D     Good          3405.     18468
##  3 D     Very Good     3470.     18542
##  4 D     Premium       3631.     18575
##  5 D     Ideal         2629.     18693
##  6 E     Fair          3682.     15584
##  7 E     Good          3424.     18236
##  8 E     Very Good     3215.     18731
##  9 E     Premium       3539.     18477
## 10 E     Ideal         2598.     18729
## # … with 25 more rows
```

---

## Exercise

**Exercise**: `starwars`, which comes with the package `dplyr`, is a data set that contains data related with characters in the movie Star Wars.

- Explorer the data set. Subset the data so that the resulting data set contains only _Human_ species with `height` data available (_Hint_: look up `is.na()`). How many observations?

- Using the previous subset, compute the mean height for each gender.

- Find the homeworld of the tallest Human male character.

---

## Normal distribution

Consider the data set `heights.txt` which consists of mother/daughter height pairs.

```r
ht <- read.table(file = "heights.txt", header = T, sep = " ")
dim(ht)
```

```
## [1] 1375    2
```

```r
head(ht)
```

```
##   Mheight Dheight
## 1    59.7    55.1
## 2    58.2    56.5
## 3    60.6    56.0
## 4    60.7    56.8
## 5    61.8    56.0
## 6    55.5    57.9
```

---

## Normal distribution, continued

.pull-left[

```r
ggplot(ht, aes(Mheight)) + 
  geom_histogram()
```

<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-8-1.png" style="display: block; margin: auto;" />
]

.pull-right[

```r
ggplot(ht, aes(Dheight)) + 
  geom_histogram()
```

<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-9-1.png" style="display: block; margin: auto;" />
]

Note that the histograms look like a bell.

---

## Normal distribution, continued

Bell shapes are observed in other data sets, too.

.pull-left[

```r
ggplot(diamonds, aes(depth)) + 
  geom_histogram()
```

<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-10-1.png" style="display: block; margin: auto;" />
]

.pull-right[

```r
ggplot(iris, aes(Sepal.Width)) + 
  geom_histogram() 
```

<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-11-1.png" style="display: block; margin: auto;" />
]

---

## Rigorous definition of normal distribution

- A continuous random variable `$X$` is said to follow a _normal distribution_ `$N(\mu, \sigma^2)$` if the probability `$P(a\leq X \leq b)$` that `$X$` is between `$a$` and `$b$` is given by 
`$$P(a\leq X \leq b) = \int_a^b \phi_{\mu,\sigma}(x)\,dx,$$`
where `$$\phi_{\mu,\sigma}(x) = \frac{1}{\sqrt{2\pi}\sigma} e^{-\frac{(x-\mu)^2}{2\sigma^2}}.$$`

- In particular, the area under the curve equals 1.

- `$\phi_{\mu,\sigma}(x)$` is a bell-shaped function called the probability distribution function (pdf) of `$X$`. `$\mu$` is called the _mean_ and `$\sigma$` is called the _standard deviation_.

- In particular, if `$\mu = 0$` and `$\sigma = 1$`, `$X$` is said to follow the _standard_ normal distribution.

---

## Normal pdf

.pull-left[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-13-1.png" style="display: block; margin: auto;" />
]

.pull-right[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-14-1.png" style="display: block; margin: auto;" />
]

_mean_ determines the center of the bell shape, while _standard deviation_ measures the spread.

---

## Why is normal distribution important?

- It has mathematically important properties (e.g. it plays a crucial role in the Central Limit Theorem).

- Related with other important distributions (e.g. `$t$`, `$F$`, `$\chi^2$`).

- Provides a good model for many real-world data (e.g. A normal distribution is a good model for data such as heights or weights).

---

## Normal pdf, continued

**Example**: the life of a fully-charged cell phone battery is normally distributed with a mean of 13 hours with a standard deviation of 1 hour. `$X$`, the life of a randomly selected fully-charged phone battery, can be assumed to follow `$N(13, 1^2)$`.

---

## Normal pdf, continued

- What is the probability that `$X \leq 13$`?

- As the curve is symmetric, the desired probability is `$0.5$`.

---

## Normal pdf, continued

- What is the probability that `$X \leq 12$`?

.pull-left[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-17-1.png" style="display: block; margin: auto;" />
]

.pull-right[
`pnorm()` computes the desired probability:

```r
pnorm(12, mean = 13, sd = 1)
```

```
## [1] 0.1586553
```

```r
pnorm(13, mean = 13, sd = 1)
```

```
## [1] 0.5
```
]

**Exercise**: compute the probability that X `$\geq$` 14.

---

## Normal pdf, continued

- What is the probability that `$11 \leq X \leq 14$`?

.pull-left[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-19-1.png" style="display: block; margin: auto;" />
]

.pull-right[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-20-1.png" style="display: block; margin: auto;" />
]

---

## Normal pdf, continued

**Exercise**:

- Compute the probability that `$X\geq 15$`, where `$X$` follows `$N(13,1^2)$`.
- Compute the probability that `$Z \leq -1$`, where `$Z$` follows `$N(0,1^2)$`.
- Determine `$z_{0.1}$` such that the probability that `$Z\leq z_{0.1}$` equals 0.9 (_hint_: look up `qnorm()`).

.pull-left[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-21-1.png" width="90%" style="display: block; margin: auto;" />
]

.pull-right[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-22-1.png" width="90%" style="display: block; margin: auto;" />
]

---

## _t_-distribution

- Class of bell-shaped distribution determined by a parameter called _degrees of freedom_.

---

## _t_-distribution, continued

**Example**: suppose `$X$` follows `$t$`-distribution with `$5$` degrees of freedom. What is the probability that `$X\leq 0$`? `$X\leq -1.5$`? `$-1 \leq X \leq 1.2$`? `$X \geq 1.6$`?

.pull-left[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-24-1.png" style="display: block; margin: auto;" />
]

.pull-right[

```r
pt(0, df = 5)
```

```
## [1] 0.5
```

```r
pt(-1.5, df = 5)
```

```
## [1] 0.09695184
```

```r
pt(1.2, df = 5) - pt(-1, df = 5)
```

```
## [1] 0.6764457
```

```r
1-pt(1.6, df = 5)
```

```
## [1] 0.08524762
```
]

---

## _t_-distribution, continued

**Exercise**: suppose `$X$` follows `$t$`-distribution with `$5$` degrees of freedom. Determine `$t_{5, 0.2}$` so that the probability that `$X\geq t_{5, 0.2} = 0.2$` (_Hint_: look up `qt()`).

---

## _t_-distribution, continued

**Exercise**: suppose `$X$` follows `$t$`-distribution with `$4$` degrees of freedom. Determine `$t_{4, 0.025}$` so that the probability that `$-t_{4, 0.025} \leq X \leq t_{4, 0.025} = 0.95$`.

---

## Using data to make inference

- Point estimate

- Confidence interval

- Statistical tests

---

## Point estimate

- Assume that the heights of all UK women follow a normal distribution `$N(\mu, \sigma^2)$`.
- Interested in determining `$\mu$` and `$\sigma$`.
- Sample mean ans sample standard deviation are known to be _good_ estimates for `$\mu$` and `$\sigma$`, respectively.

.pull-left[
<img src="Drexel_2020_Lecture_4_files/figure-html/unnamed-chunk-28-1.png" width="95%" style="display: block; margin: auto;" />
]

.pull-right[

```r
mean(ht$Mheight)
```

```
## [1] 62.4528
```

```r
sd(ht$Mheight)
```

```
## [1] 2.355103
```

**Exercise**: compute the sample mean and sample standard deviation of `Sepal.Width` in `iris` data set.

]

---

## Confidence interval

Point estimate is too naive -- how much can we trust 62.4528 as a reliable estimator for `$\mu$`, the true population mean?

**Example**: suppose that we are measuring the battery life which is assumed to follow a normal distribution. Assume that we collected two sets of samples of size 5:
`$$\text{Sample 1: }\{12, 13, 13, 14, 13\}$$`
`$$\text{Sample 2: }\{10, 11, 16, 13, 15\}$$`

```r
sample1 <- c(12, 13, 13, 14, 13); sample2 <- c(10, 11, 16, 13, 15)
c(mean(sample1), mean(sample2))
```

```
## [1] 13 13
```

```r
c(sd(sample1), sd(sample2))
```

```
## [1] 0.7071068 2.5495098
```

---