---
title: "Data Manipulation with dplyr"
author: "Evan L. Ray"
date: "September 15, 2017"
output: ioslides_presentation
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE)
```

## The HELPrct data

```{r preprocess_HELPrct, echo=FALSE, warning=FALSE, message=FALSE}
library(mosaicData)
library(dplyr)
library(tidyr)
library(ggplot2)
library(readr)

HELPrct <- HELPrct %>%
  select(
    age,
    homeless,
    substance,
    cesd,
    mcs
  ) %>%
  rename(
    dep_score = cesd,
    mental_score = mcs)
```

```{r head_HELPrct, echo = TRUE}
head(HELPrct)
```

The depression score is a number between 0 and 60, with high scores indicating greater depressive symptoms.

The mental component score is a number between 0 and 100, with 0 indicating disability and 100 indicating no disability

## mutate: add a new variable

Let's create a combined measure of mental health, equal to (60 - depression score) + mental component score:

```{r, echo = TRUE}
HELPrct <- mutate(HELPrct,
  mental_health = (60 - dep_score) + mental_score)
head(HELPrct)
```

## arrange: sort according to a variable

```{r, echo = TRUE}
HELPrct_sorted <- arrange(HELPrct, mental_health)
head(HELPrct_sorted)
```

```{r, echo = TRUE}
HELPrct_sorted <- arrange(HELPrct, desc(mental_health))
head(HELPrct_sorted)
```

## Statistics for numeric variables

```{r, echo = TRUE}
mean(HELPrct$mental_health)
sd(HELPrct$mental_health)
min(HELPrct$mental_health)
max(HELPrct$mental_health)
```

## Statistics for numeric variables, cont.

```{r, echo = TRUE}
median(HELPrct$mental_health)
quantile(HELPrct$mental_health, p = 0.5)
quantile(HELPrct$mental_health, p = c(0.25, 0.75))
```

## Statistics for numeric variables, cont.

```{r, echo = TRUE}
quantile(HELPrct$mental_health, p = 0.75) -
  quantile(HELPrct$mental_health, p = 0.25)

IQR(HELPrct$mental_health)
```

## summarize: statistics in a data frame

```{r, echo = TRUE}
quantile(HELPrct$mental_health, p = 0.25)
median(HELPrct$mental_health)
HELPrct_summary <- summarize(HELPrct,
  Q1_mh = quantile(mental_health, p = 0.25),
  median_mh = median(mental_health),
  Q3_mh = quantile(mental_health, p = 0.75))
HELPrct_summary
```

## Tangent: Putting the ideas together

```{r, echo = TRUE, fig.height=2.7}
ggplot() +
  geom_histogram(mapping = aes(x = mental_health, y = ..density..),
    data = HELPrct) +
  geom_vline(mapping = aes(xintercept = Q1_mh), data = HELPrct_summary) +
  geom_vline(mapping = aes(xintercept = median_mh), data = HELPrct_summary) +
  geom_vline(mapping = aes(xintercept = Q3_mh), data = HELPrct_summary)
```

## group_by: set up a grouping

Nothing visible happens -- it's all behind the scenes!

```{r, echo = TRUE}
HELPrct_grouped <- group_by(HELPrct, substance)
head(HELPrct_grouped)
```

## group_by continued...

...but now summaries are calculated separately for each group:

```{r, echo = TRUE}
HELPrct_grouped_summary <- summarize(HELPrct_grouped,
  mean_mh = mean(mental_health),
  sd_mh = sd(mental_health),
  iqr_mh = IQR(mental_health))
HELPrct_grouped_summary
```

## %>% -- the pipe operator

Combine multiple steps that happen one after the other.

```{r, echo = TRUE}
HELPrct_grouped_summary2 <- HELPrct %>%
  group_by(substance) %>%
  summarize(
    mean_mh = mean(mental_health),
    sd_mh = sd(mental_health))

HELPrct_grouped_summary2
```

Read as "... and then ..."

## %>% -- how does it work?

The result of the calculation on the left hand side of `%>%` is used as the first argument to the function on the right hand side.

```{r, echo = TRUE, eval = FALSE}
HELPrct %>%
  group_by(substance)
```

is the same as

```{r, echo = TRUE, eval = FALSE}
group_by(HELPrct, substance)
```

## Back to statistics...

 * If the distribution is **unimodal and symmetric**
    * Measure *center* with the **mean**: $\bar{y} = \frac{1}{n} \sum_{i = 1}^n y_i$
    * Measure *spread* with the **standard deviation** -- think "average distance from the mean" (though this is not formally quite right): $s = \sqrt{\frac{\sum_{i = 1}^n (y_i - \bar{y})^2}{n - 1}}$
 * If the distribution has **multiple modes**, or is **skewed**, or has **outliers**:
    * Measure *center* with the **median** -- half of the data are less than the median and half of the data are greater than the median.
    * Measure *spread* with the **Interquartile Range (IQR)** -- the width of the interval that contains the middle half of the data.

## Why?  Recall Bechdel data

```{r, echo = FALSE, fig.height=4.5, message=FALSE}
movies <- read_csv("https://mhc-stat140-2017.github.io/data/bechdel/bechdel.csv") %>%
  mutate(
    bechdel_test = factor(
      bechdel_test,
      ordered = TRUE,
      levels = c("nowomen", "notalk", "men", "dubious", "ok")),
    bechdel_test_binary = factor(
      bechdel_test_binary,
      ordered = TRUE,
      levels = c("FAIL", "PASS")),
    mpaa_rating = factor(
      mpaa_rating,
      ordered = TRUE,
      levels = c("UNRATED", "NOT RATED", "G", "PG", "TV-PG", "PG-13", "TV-14", "R", "NC-17"))
  ) %>%
  filter(!is.na(intgross_2013))
```



```{r, fig.height=3.5, echo = TRUE, message=FALSE}
movie_stats <- summarize(movies, mean_intgross = mean(intgross_2013), median_intgross = median(intgross_2013))
ggplot() +
  geom_histogram(aes(x = intgross_2013), data = movies) + 
  geom_vline(aes(xintercept = mean_intgross), color = "orange", data = movie_stats) +
  geom_vline(aes(xintercept = median_intgross), color = "blue", data = movie_stats)
```

## Outlier effects on summary statistics

```{r, echo = TRUE}
summary_all_obs <- movies %>%
  summarize(mean_intgross = mean(intgross_2013), median_intgross = median(intgross_2013), sd_mpg = sd(intgross_2013), iqr_intgross = IQR(intgross_2013))
summary_all_obs
```

```{r, echo = TRUE}
summary_drop_outliers <- movies %>%
  filter(intgross_2013 < 3 * 10^9) %>%
  summarize(mean_intgross = mean(intgross_2013), median_intgross = median(intgross_2013), sd_intgross = sd(intgross_2013), iqr_intgross = IQR(intgross_2013))
summary_drop_outliers
```
