r4ds/EDA.Rmd

596 lines
32 KiB
Plaintext
Raw Normal View History

2016-07-19 22:39:00 +08:00
# Exploratory Data Analysis
2016-07-15 22:42:51 +08:00
## Introduction
This chapter will show you how to use visualisation and transformation to explore your data in a systematic way, a task that statisticians call exploratory data analysis, or EDA for short. EDA is an iterative cycle. You:
2016-07-15 22:42:51 +08:00
2016-07-23 05:00:11 +08:00
1. Generate questions about your data.
2016-07-15 22:42:51 +08:00
1. Search for answers by visualising, transforming, and modelling your data.
2016-07-15 22:42:51 +08:00
2016-10-03 21:36:51 +08:00
1. Use what you learn to refine your questions and/or generate new questions.
2017-01-13 10:56:09 +08:00
EDA is not a formal process with a strict set of rules. More than anything, EDA is a state of mind. During the initial phases of EDA you should feel free to investigate every idea that occurs to you. Some of these ideas will pan out, and some will be dead ends. As your exploration continues, you will home in on a few particularly productive areas that you'll eventually write up and communicate to others.
2016-07-27 03:36:53 +08:00
EDA is an important part of any data analysis, even if the questions are handed to you on a platter, because you always need to investigate the quality of your data. Data cleaning is just one application of EDA: you ask questions about whether your data meets your expectations or not. To do data cleaning, you'll need to deploy all the tools of EDA: visualisation, transformation, and modelling.
2016-07-15 22:42:51 +08:00
### Prerequisites
2016-07-27 03:36:53 +08:00
In this chapter we'll combine what you've learned about dplyr and ggplot2 to interactively ask questions, answer them with data, and then ask new questions.
2016-07-15 22:42:51 +08:00
2016-07-23 05:00:11 +08:00
```{r setup, message = FALSE}
2016-10-04 01:30:24 +08:00
library(tidyverse)
2016-07-15 22:42:51 +08:00
```
2016-05-03 05:11:00 +08:00
## Questions
2016-05-03 05:11:00 +08:00
2016-07-12 06:32:36 +08:00
> "There are no routine statistical questions, only questionable statistical
> routines." --- Sir David Cox
2016-05-03 05:11:00 +08:00
2016-07-15 22:42:51 +08:00
> "Far better an approximate answer to the right question, which is often
> vague, than an exact answer to the wrong question, which can always be made
2016-07-23 05:00:11 +08:00
> precise." --- John Tukey
2016-05-03 05:11:00 +08:00
Your goal during EDA is to develop an understanding of your data. The easiest way to do this is to use questions as tools to guide your investigation. When you ask a question, the question focuses your attention on a specific part of your dataset and helps you decide which graphs, models, or transformations to make.
2016-05-03 05:11:00 +08:00
2016-10-03 21:36:51 +08:00
EDA is fundamentally a creative process. And like most creative processes, the key to asking _quality_ questions is to generate a large _quantity_ of questions. It is difficult to ask revealing questions at the start of your analysis because you do not know what insights are contained in your dataset. On the other hand, each new question that you ask will expose you to a new aspect of your data and increase your chance of making a discovery. You can quickly drill down into the most interesting parts of your data---and develop a set of thought-provoking questions---if you follow up each question with a new question based on what you find.
2016-07-15 22:42:51 +08:00
There is no rule about which questions you should ask to guide your research. However, two types of questions will always be useful for making discoveries within your data. You can loosely word these questions as:
2016-05-03 05:11:00 +08:00
2016-10-05 00:08:12 +08:00
1. What type of variation occurs within my variables?
2016-05-05 02:41:46 +08:00
2016-10-05 00:08:12 +08:00
1. What type of covariation occurs between my variables?
2016-05-22 10:14:30 +08:00
The rest of this chapter will look at these two questions. I'll explain what variation and covariation are, and I'll show you several ways to answer each question. To make the discussion easier, let's define some terms:
2016-07-15 22:42:51 +08:00
* A __variable__ is a quantity, quality, or property that you can measure.
2016-07-15 22:42:51 +08:00
* A __value__ is the state of a variable when you measure it. The value of a
variable may change from measurement to measurement.
2016-05-03 05:11:00 +08:00
2016-07-15 22:42:51 +08:00
* An __observation__ is a set of measurements made under similar conditions
(you usually make all of the measurements in an observation at the same
time and on the same object). An observation will contain several values,
each associated with a different variable. I'll sometimes refer to
an observation as a data point.
2016-07-23 05:00:11 +08:00
* __Tabular data__ is a set of values, each associated with a variable and an
2016-07-15 22:42:51 +08:00
observation. Tabular data is _tidy_ if each value is placed in its own
"cell", each variable in its own column, and each observation in its own
row.
2016-11-04 02:18:13 +08:00
So far, all of the data that you've seen has been tidy. In real-life, most data isn't tidy, so we'll come back to these ideas again in [tidy data].
2016-05-22 10:14:30 +08:00
## Variation
2016-10-04 21:28:46 +08:00
**Variation** is the tendency of the values of a variable to change from measurement to measurement. You can see variation easily in real life; if you measure any continuous variable twice, you will get two different results. This is true even if you measure quantities that are constant, like the speed of light. Each of your measurements will include a small amount of error that varies from measurement to measurement. Categorical variables can also vary if you measure across different subjects (e.g. the eye colors of different people), or different times (e.g. the energy levels of an electron at different moments).
2016-11-04 02:18:13 +08:00
Every variable has its own pattern of variation, which can reveal interesting information. The best way to understand that pattern is to visualise the distribution of the variable's values.
### Visualising distributions
2016-08-16 20:22:41 +08:00
How you visualise the distribution of a variable will depend on whether the variable is categorical or continuous. A variable is **categorical** if it can only take one of a small set of values. In R, categorical variables are usually saved as factors or character vectors. To examine the distribution of a categorical variable, use a bar chart:
```{r}
ggplot(data = diamonds) +
geom_bar(mapping = aes(x = cut))
```
2016-07-23 05:00:11 +08:00
The height of the bars displays how many observations occurred with each x value. You can compute these values manually with `dplyr::count()`:
```{r}
2016-08-01 00:32:16 +08:00
diamonds %>%
count(cut)
```
2016-10-03 21:36:51 +08:00
A variable is **continuous** if it can take any of an infinite set of ordered values. Numbers and date-times are two examples of continuous variables. To examine the distribution of a continuous variable, use a histogram:
2016-07-15 22:42:51 +08:00
```{r}
ggplot(data = diamonds) +
geom_histogram(mapping = aes(x = carat), binwidth = 0.5)
```
2016-07-23 05:00:11 +08:00
You can compute this by hand by combining `dplyr::count()` and `ggplot2::cut_width()`:
```{r}
2016-08-01 00:32:16 +08:00
diamonds %>%
count(cut_width(carat, 0.5))
2016-07-23 05:00:11 +08:00
```
2016-11-04 02:18:13 +08:00
A histogram divides the x-axis into equally spaced bins and then uses the height of a bar to display the number of observations that fall in each bin. In the graph above, the tallest bar shows that almost 30,000 observations have a `carat` value between 0.25 and 0.75, which are the left and right edges of the bar.
2016-07-27 05:01:19 +08:00
You can set the width of the intervals in a histogram with the `binwidth` argument, which is measured in the units of the `x` variable. You should always explore a variety of binwidths when working with histograms, as different binwidths can reveal different patterns. For example, here is how the graph above looks when we zoom into just the diamonds with a size of less than three carats and choose a smaller binwidth.
2016-07-15 22:42:51 +08:00
```{r}
2016-08-01 00:32:16 +08:00
smaller <- diamonds %>%
filter(carat < 3)
2016-07-15 22:42:51 +08:00
ggplot(data = smaller, mapping = aes(x = carat)) +
geom_histogram(binwidth = 0.1)
```
2016-07-15 22:42:51 +08:00
If you wish to overlay multiple histograms in the same plot, I recommend using `geom_freqpoly()` instead of `geom_histogram()`. `geom_freqpoly()` performs the same calculation as `geom_histogram()`, but instead of displaying the counts with bars, uses lines instead. It's much easier to understand overlapping lines than bars.
2016-07-15 22:42:51 +08:00
```{r}
2016-08-01 00:32:16 +08:00
ggplot(data = smaller, mapping = aes(x = carat, colour = cut)) +
2016-07-15 22:42:51 +08:00
geom_freqpoly(binwidth = 0.1)
```
2016-08-16 20:22:41 +08:00
There are a few challenges with this type of plot, which we will come back to in [visualising a categorical and a continuous variable](#cat-cont).
2016-08-01 00:32:16 +08:00
2016-10-05 00:08:12 +08:00
Now that you can visualise variation, what should you look for in your plots? And what type of follow-up questions should you ask? I've put together a list below of the most useful types of information that you will find in your graphs, along with some follow-up questions for each type of information. The key to asking good follow-up questions will be to rely on your curiosity (What do you want to learn more about?) as well as your skepticism (How could this be misleading?).
2016-07-15 22:42:51 +08:00
### Typical values
2016-07-23 05:00:11 +08:00
In both bar charts and histograms, tall bars show the common values of a variable, and shorter bars show less-common values. Places that do not have bars reveal values that were not seen in your data. To turn this information into useful questions, look for anything unexpected:
2016-07-15 22:42:51 +08:00
* Which values are the most common? Why?
* Which values are rare? Why? Does that match your expectations?
2016-07-15 22:42:51 +08:00
* Can you see any unusual patterns? What might explain them?
As an example, the histogram below suggests several interesting questions:
* Why are there more diamonds at whole carats and common fractions of carats?
* Why are there more diamonds slightly to the right of each peak than there
are slightly to the left of each peak?
* Why are there no diamonds bigger than 3 carats?
```{r}
ggplot(data = smaller, mapping = aes(x = carat)) +
geom_histogram(binwidth = 0.01)
```
2016-07-15 22:42:51 +08:00
Clusters of similar values suggest that subgroups exist in your data. To understand the subgroups, ask:
2016-07-15 22:42:51 +08:00
* How are the observations within each cluster similar to each other?
2016-07-15 22:42:51 +08:00
* How are the observations in separate clusters different from each other?
* How can you explain or describe the clusters?
* Why might the appearance of clusters be misleading?
2016-11-04 02:18:13 +08:00
The histogram below shows the length (in minutes) of 272 eruptions of the Old Faithful Geyser in Yellowstone National Park. Eruption times appear to be clustered into two groups: there are short eruptions (of around 2 minutes) and long eruptions (4-5 minutes), but little in between.
2016-07-15 22:42:51 +08:00
```{r}
ggplot(data = faithful, mapping = aes(x = eruptions)) +
geom_histogram(binwidth = 0.25)
```
2016-07-23 05:00:11 +08:00
Many of the questions above will prompt you to explore a relationship *between* variables, for example, to see if the values of one variable can explain the behavior of another variable. We'll get to that shortly.
2016-07-15 22:42:51 +08:00
### Unusual values
Outliers are observations that are unusual; data points that don't seem to fit the pattern. Sometimes outliers are data entry errors; other times outliers suggest important new science. When you have a lot of data, outliers are sometimes difficult to see in a histogram. For example, take the distribution of the `y` variable from the diamonds dataset. The only evidence of outliers is the unusually wide limits on the x-axis.
2016-07-15 22:42:51 +08:00
```{r}
ggplot(diamonds) +
geom_histogram(mapping = aes(x = y), binwidth = 0.5)
2016-07-15 22:42:51 +08:00
```
2017-05-03 20:42:40 +08:00
There are so many observations in the common bins that the rare bins are so short that you can't see them (although maybe if you stare intently at 0 you'll spot something). To make it easy to see the unusual values, we need to zoom to small values of the y-axis with `coord_cartesian()`:
2016-07-15 22:42:51 +08:00
```{r}
ggplot(diamonds) +
geom_histogram(mapping = aes(x = y), binwidth = 0.5) +
2016-07-15 22:42:51 +08:00
coord_cartesian(ylim = c(0, 50))
```
(`coord_cartesian()` also has an `xlim()` argument for when you need to zoom into the x-axis. ggplot2 also has `xlim()` and `ylim()` functions that work slightly differently: they throw away the data outside the limits.)
This allows us to see that there are three unusual values: 0, ~30, and ~60. We pluck them out with dplyr:
2016-10-04 05:08:44 +08:00
```{r, include = FALSE}
old <- options(tibble.print_max = 10, tibble.print_min = 10)
```
2016-07-15 22:42:51 +08:00
```{r}
unusual <- diamonds %>%
filter(y < 3 | y > 20) %>%
select(price, x, y, z) %>%
2016-07-15 22:42:51 +08:00
arrange(y)
unusual
```
2016-10-04 05:08:44 +08:00
```{r, include = FALSE}
options(old)
```
2016-07-27 03:36:53 +08:00
The `y` variable measures one of the three dimensions of these diamonds, in mm. We know that diamonds can't have a width of 0mm, so these values must be incorrect. We might also suspect that measurements of 32mm and 59mm are implausible: those diamonds are over an inch long, but don't cost hundreds of thousands of dollars!
2016-07-15 22:42:51 +08:00
2016-08-01 00:32:16 +08:00
It's good practice to repeat your analysis with and without the outliers. If they have minimal effect on the results, and you can't figure out why they're there, it's reasonable to replace them with missing values, and move on. However, if they have a substantial effect on your results, you shouldn't drop them without justification. You'll need to figure out what caused them (e.g. a data entry error) and disclose that you removed them in your write-up.
2016-07-15 22:42:51 +08:00
### Exercises
1. Explore the distribution of each of the `x`, `y`, and `z` variables
in `diamonds`. What do you learn? Think about a diamond and how you
might decide which dimension is the length, width, and depth.
1. Explore the distribution of `price`. Do you discover anything unusual
or surprising? (Hint: Carefully think about the `binwidth` and make sure
2016-08-13 01:43:09 +08:00
you try a wide range of values.)
2016-07-18 21:47:52 +08:00
2016-08-26 23:15:59 +08:00
1. How many diamonds are 0.99 carat? How many are 1 carat? What
2016-07-27 05:01:19 +08:00
do you think is the cause of the difference?
1. Compare and contrast `coord_cartesian()` vs `xlim()` or `ylim()` when
2016-07-15 22:42:51 +08:00
zooming in on a histogram. What happens if you leave `binwidth` unset?
What happens if you try and zoom so only half a bar shows?
2016-07-15 22:42:51 +08:00
## Missing values
2016-05-05 02:41:46 +08:00
2016-07-15 22:42:51 +08:00
If you've encountered unusual values in your dataset, and simply want to move on to the rest of your analysis, you have two options.
1. Drop the entire row with the strange values:
2016-08-01 00:32:16 +08:00
```{r, eval = FALSE}
diamonds2 <- diamonds %>%
filter(between(y, 3, 20))
```
2016-07-15 22:42:51 +08:00
I don't recommend this option because just because one measurement
is invalid, doesn't mean all the measurements are. Additionally, if you
2016-07-23 05:00:11 +08:00
have low quality data, by time that you've applied this approach to every
variable you might find that you don't have any data left!
2016-07-15 22:42:51 +08:00
1. Instead, I recommend replacing the unusual values with missing values.
2016-10-03 21:36:51 +08:00
The easiest way to do this is to use `mutate()` to replace the variable
2016-07-15 22:42:51 +08:00
with a modified copy. You can use the `ifelse()` function to replace
unusual values with `NA`:
```{r}
diamonds2 <- diamonds %>%
mutate(y = ifelse(y < 3 | y > 20, NA, y))
```
`ifelse()` has three arguments. The first argument `test` should be a logical vector. The result will contain the value of the second argument, `yes`, when `test` is `TRUE`, and the value of the third argument, `no`, when it is false. Alternatively to ifelse, use `dplyr::case_when()`. `case_when()` is particularly useful inside mutate when you want to create a new variable that relies on a complex combination of existing variables.
2016-07-23 05:00:11 +08:00
Like R, ggplot2 subscribes to the philosophy that missing values should never silently go missing. It's not obvious where you should plot missing values, so ggplot2 doesn't include them in the plot, but it does warn that they've been removed:
2016-07-18 22:52:55 +08:00
```{r, dev = "png"}
2016-07-15 22:42:51 +08:00
ggplot(data = diamonds2, mapping = aes(x = x, y = y)) +
geom_point()
```
2016-05-03 05:11:00 +08:00
To suppress that warning, set `na.rm = TRUE`:
2016-07-18 21:47:52 +08:00
```{r, eval = FALSE}
2016-07-15 22:42:51 +08:00
ggplot(data = diamonds2, mapping = aes(x = x, y = y)) +
geom_point(na.rm = TRUE)
```
2017-05-03 20:42:40 +08:00
Other times you want to understand what makes observations with missing values different to observations with recorded values. For example, in `nycflights13::flights`, missing values in the `dep_time` variable indicate that the flight was cancelled. So you might want to compare the scheduled departure times for cancelled and non-cancelled times. You can do this by making a new variable with `is.na()`.
2016-07-15 22:42:51 +08:00
```{r}
nycflights13::flights %>%
mutate(
cancelled = is.na(dep_time),
sched_hour = sched_dep_time %/% 100,
2016-07-23 05:00:11 +08:00
sched_min = sched_dep_time %% 100,
2016-07-15 22:42:51 +08:00
sched_dep_time = sched_hour + sched_min / 60
) %>%
ggplot(mapping = aes(sched_dep_time)) +
geom_freqpoly(mapping = aes(colour = cancelled), binwidth = 1/4)
```
2016-07-18 21:47:52 +08:00
However this plot isn't great because there are many more non-cancelled flights than cancelled flights. In the next section we'll explore some techniques for improving this comparison.
2016-07-15 22:42:51 +08:00
### Exercises
2016-07-18 21:47:52 +08:00
1. What happens to missing values in a histogram? What happens to missing
2016-10-03 21:36:51 +08:00
values in a bar chart? Why is there a difference?
2016-07-18 21:47:52 +08:00
1. What does `na.rm = TRUE` do in `mean()` and `sum()`?
## Covariation
If variation describes the behavior _within_ a variable, covariation describes the behavior _between_ variables. **Covariation** is the tendency for the values of two or more variables to vary together in a related way. The best way to spot covariation is to visualise the relationship between two or more variables. How you do that should again depend on the type of variables involved.
2016-08-01 00:32:16 +08:00
### A categorical and continuous variable {#cat-cont}
2016-10-03 21:36:51 +08:00
It's common to want to explore the distribution of a continuous variable broken down by a categorical variable, as in the previous frequency polygon. The default appearance of `geom_freqpoly()` is not that useful for that sort of comparison because the height is given by the count. That means if one of the groups is much smaller than the others, it's hard to see the differences in shape. For example, let's explore how the price of a diamond varies with its quality:
```{r}
2016-07-15 22:42:51 +08:00
ggplot(data = diamonds, mapping = aes(x = price)) +
geom_freqpoly(mapping = aes(colour = cut), binwidth = 500)
```
2016-07-15 22:42:51 +08:00
It's hard to see the difference in distribution because the overall counts differ so much:
2016-07-23 05:07:07 +08:00
```{r, fig.width = "50%", fig.width = 4}
ggplot(diamonds) +
geom_bar(mapping = aes(x = cut))
```
To make the comparison easier we need to swap what is displayed on the y-axis. Instead of displaying count, we'll display __density__, which is the count standardised so that the area under each frequency polygon is one.
2016-07-15 22:42:51 +08:00
```{r}
ggplot(data = diamonds, mapping = aes(x = price, y = ..density..)) +
geom_freqpoly(mapping = aes(colour = cut), binwidth = 500)
2016-07-15 22:42:51 +08:00
```
2016-07-27 03:36:53 +08:00
There's something rather surprising about this plot - it appears that fair diamonds (the lowest quality) have the highest average price! But maybe that's because frequency polygons are a little hard to interpret - there's a lot going on in this plot.
2016-07-15 22:42:51 +08:00
Another alternative to display the distribution of a continuous variable broken down by a categorical variable is the boxplot. A **boxplot** is a type of visual shorthand for a distribution of values that is popular among statisticians. Each boxplot consists of:
2016-07-15 22:42:51 +08:00
* A box that stretches from the 25th percentile of the distribution to the
2016-07-27 05:01:19 +08:00
75th percentile, a distance known as the interquartile range (IQR). In the
2016-07-15 22:42:51 +08:00
middle of the box is a line that displays the median, i.e. 50th percentile,
of the distribution. These three lines give you a sense of the spread of the
distribution and whether or not the distribution is symmetric about the
median or skewed to one side.
2016-07-15 22:42:51 +08:00
* Visual points that display observations that fall more than 1.5 times the
IQR from either edge of the box. These outlying points are unusual
2016-10-03 21:36:51 +08:00
so are plotted individually.
2016-07-15 22:42:51 +08:00
* A line (or whisker) that extends from each end of the box and goes to the
farthest non-outlier point in the distribution.
2016-07-18 22:52:55 +08:00
```{r, echo = FALSE, out.width = "100%"}
2016-07-15 22:42:51 +08:00
knitr::include_graphics("images/EDA-boxplot.png")
```
2016-07-15 22:42:51 +08:00
Let's take a look at the distribution of price by cut using `geom_boxplot()`:
```{r fig.height = 3}
2016-07-15 22:42:51 +08:00
ggplot(data = diamonds, mapping = aes(x = cut, y = price)) +
geom_boxplot()
```
2016-07-27 03:36:53 +08:00
We see much less information about the distribution, but the boxplots are much more compact so we can more easily compare them (and fit more on one plot). It supports the counterintuitive finding that better quality diamonds are cheaper on average! In the exercises, you'll be challenged to figure out why.
2016-07-15 22:42:51 +08:00
2016-08-26 23:15:59 +08:00
`cut` is an ordered factor: fair is worse than good, which is worse than very good and so on. Many categorical variables don't have such an intrinsic order, so you might want to reorder them to make a more informative display. One way to do that is with the `reorder()` function.
2016-08-01 00:32:16 +08:00
For example, take the `class` variable in the `mpg` dataset. You might be interested to know how highway mileage varies across classes:
2016-07-15 22:42:51 +08:00
```{r}
ggplot(data = mpg, mapping = aes(x = class, y = hwy)) +
geom_boxplot()
```
2016-08-01 00:32:16 +08:00
To make the trend easier to see, we can reorder `class` based on the median value of `hwy`:
```{r fig.height = 3}
ggplot(data = mpg) +
geom_boxplot(mapping = aes(x = reorder(class, hwy, FUN = median), y = hwy))
```
2016-07-15 22:42:51 +08:00
If you have long variable names, `geom_boxplot()` will work better if you flip it 90°. You can do that with `coord_flip()`.
```{r}
ggplot(data = mpg) +
geom_boxplot(mapping = aes(x = reorder(class, hwy, FUN = median), y = hwy)) +
coord_flip()
```
2016-07-15 22:42:51 +08:00
#### Exercises
2016-07-18 21:47:52 +08:00
1. Use what you've learned to improve the visualisation of the departure times
of cancelled vs. non-cancelled flights.
1. What variable in the diamonds dataset is most important for predicting
the price of a diamond? How is that variable correlated with cut?
2016-07-23 05:00:11 +08:00
Why does the combination of those two relationships lead to lower quality
diamonds being more expensive?
2016-07-27 03:36:53 +08:00
1. Install the ggstance package, and create a horizontal boxplot.
2016-07-18 21:47:52 +08:00
How does this compare to using `coord_flip()`?
2016-07-15 22:42:51 +08:00
1. One problem with boxplots is that they were developed in an era of
much smaller datasets and tend to display a prohibitively large
2016-07-15 22:42:51 +08:00
number of "outlying values". One approach to remedy this problem is
the letter value plot. Install the lvplot package, and try using
2016-07-27 03:36:53 +08:00
`geom_lv()` to display the distribution of price vs cut. What
2016-07-15 22:42:51 +08:00
do you learn? How do you interpret the plots?
1. Compare and contrast `geom_violin()` with a facetted `geom_histogram()`,
2016-07-23 05:00:11 +08:00
or a coloured `geom_freqpoly()`. What are the pros and cons of each
2016-07-15 22:42:51 +08:00
method?
1. If you have a small dataset, it's sometimes useful to use `geom_jitter()`
to see the relationship between a continuous and categorical variable.
The ggbeeswarm package provides a number of methods similar to
`geom_jitter()`. List them and briefly describe what each one does.
2016-07-15 22:42:51 +08:00
2016-07-23 05:00:11 +08:00
### Two categorical variables
2016-07-15 22:42:51 +08:00
2016-07-23 05:00:11 +08:00
To visualise the covariation between categorical variables, you'll need to count the number of observations for each combination. One way to do that is to rely on the built-in `geom_count()`:
2016-07-15 22:42:51 +08:00
```{r}
ggplot(data = diamonds) +
geom_count(mapping = aes(x = cut, y = color))
```
2016-07-23 05:00:11 +08:00
The size of each circle in the plot displays how many observations occurred at each combination of values. Covariation will appear as a strong correlation between specific x values and specific y values.
Another approach is to compute the count with dplyr:
2016-07-15 22:42:51 +08:00
```{r}
2016-08-01 00:32:16 +08:00
diamonds %>%
count(color, cut)
2016-07-15 22:42:51 +08:00
```
2016-07-23 05:00:11 +08:00
Then visualise with `geom_tile()` and the fill aesthetic:
2016-07-15 22:42:51 +08:00
```{r}
diamonds %>%
count(color, cut) %>%
ggplot(mapping = aes(x = color, y = cut)) +
2016-08-01 00:47:51 +08:00
geom_tile(mapping = aes(fill = n))
2016-07-15 22:42:51 +08:00
```
2016-07-27 03:36:53 +08:00
If the categorical variables are unordered, you might want to use the seriation package to simultaneously reorder the rows and columns in order to more clearly reveal interesting patterns. For larger plots, you might want to try the d3heatmap or heatmaply packages, which create interactive plots.
2016-07-15 22:42:51 +08:00
#### Exercises
2016-07-15 22:42:51 +08:00
2016-07-23 05:00:11 +08:00
1. How could you rescale the count dataset above to more clearly show
the distribution of cut within colour, or colour within cut?
2016-07-15 22:42:51 +08:00
2016-07-23 05:00:11 +08:00
1. Use `geom_tile()` together with dplyr to explore how average flight
delays vary by destination and month of year. What makes the
plot difficult to read? How could you improve it?
2016-07-15 22:42:51 +08:00
1. Why is it slightly better to use `aes(x = color, y = cut)` rather
2016-07-18 21:47:52 +08:00
than `aes(x = cut, y = color)` in the example above?
2016-07-23 05:00:11 +08:00
### Two continuous variables
2016-07-18 21:47:52 +08:00
You've already seen one great way to visualise the covariation between two continuous variables: draw a scatterplot with `geom_point()`. You can see covariation as a pattern in the points. For example, you can see an exponential relationship between the carat size and price of a diamond.
2016-07-18 22:52:55 +08:00
```{r, dev = "png"}
ggplot(data = diamonds) +
geom_point(mapping = aes(x = carat, y = price))
```
2016-10-04 05:08:44 +08:00
Scatterplots become less useful as the size of your dataset grows, because points begin to overplot, and pile up into areas of uniform black (as above).
You've already seen one way to fix the problem: using the `alpha` aesthetic to add transparency.
2016-07-18 21:47:52 +08:00
2016-07-18 22:52:55 +08:00
```{r, dev = "png"}
ggplot(data = diamonds) +
2016-10-04 05:08:44 +08:00
geom_point(mapping = aes(x = carat, y = price), alpha = 1 / 100)
2016-07-18 21:47:52 +08:00
```
2016-10-04 22:01:45 +08:00
But using transparency can be challenging for very large datasets. Another solution is to use bin. Previously you used `geom_histogram()` and `geom_freqpoly()` to bin in one dimension. Now you'll learn how to use `geom_bin2d()` and `geom_hex()` to bin in two dimensions.
2016-10-04 05:08:44 +08:00
`geom_bin2d()` and `geom_hex()` divide the coordinate plane into 2d bins and then use a fill color to display how many points fall into each bin. `geom_bin2d()` creates rectangular bins. `geom_hex()` creates hexagonal bins. You will need to install the hexbin package to use `geom_hex()`.
```{r, fig.asp = 1, out.width = "50%", fig.align = "default", message = FALSE}
ggplot(data = smaller) +
geom_bin2d(mapping = aes(x = carat, y = price))
2016-05-22 10:14:30 +08:00
# install.packages("hexbin")
ggplot(data = smaller) +
geom_hex(mapping = aes(x = carat, y = price))
```
2016-07-31 23:41:35 +08:00
Another option is to bin one continuous variable so it acts like a categorical variable. Then you can use one of the techniques for visualising the combination of a categorical and a continuous variable that you learned about. For example, you could bin `carat` and then for each group, display a boxplot:
```{r}
ggplot(data = smaller, mapping = aes(x = carat, y = price)) +
geom_boxplot(mapping = aes(group = cut_width(carat, 0.1)))
```
`cut_width(x, width)`, as used above, divides `x` into bins of width `width`. By default, boxplots look roughly the same (apart from number of outliers) regardless of how many observations there are, so it's difficult to tell that each boxplot summarises a different number of points. One way to show that is to make the width of the boxplot proportional to the number of points with `varwidth = TRUE`.
2016-05-23 01:35:52 +08:00
2016-07-18 21:47:52 +08:00
Another approach is to display approximately the same number of points in each bin. That's the job of `cut_number()`:
2016-05-23 01:35:52 +08:00
```{r}
2016-07-18 21:47:52 +08:00
ggplot(data = smaller, mapping = aes(x = carat, y = price)) +
geom_boxplot(mapping = aes(group = cut_number(carat, 20)))
2016-05-23 01:35:52 +08:00
```
2016-07-18 21:47:52 +08:00
#### Exercises
2016-05-23 01:35:52 +08:00
2016-07-18 21:47:52 +08:00
1. Instead of summarising the conditional distribution with a boxplot, you
could use a frequency polygon. What do you need to consider when using
2016-07-27 03:36:53 +08:00
`cut_width()` vs `cut_number()`? How does that impact a visualisation of
2016-07-18 21:47:52 +08:00
the 2d distribution of `carat` and `price`?
2016-05-23 01:35:52 +08:00
2016-07-18 21:47:52 +08:00
1. Visualise the distribution of carat, partitioned by price.
1. How does the price distribution of very large diamonds compare to small
diamonds? Is it as you expect, or does it surprise you?
2016-07-18 21:47:52 +08:00
1. Combine two of the techniques you've learned to visualise the
combined distribution of cut, carat, and price.
1. Two dimensional plots reveal outliers that are not visible in one
dimensional plots. For example, some points in the plot below have an
2016-07-27 05:01:19 +08:00
unusual combination of `x` and `y` values, which makes the points outliers
even though their `x` and `y` values appear normal when examined separately.
2016-07-18 21:47:52 +08:00
2016-07-18 22:52:55 +08:00
```{r, dev = "png"}
2016-07-18 21:47:52 +08:00
ggplot(data = diamonds) +
geom_point(mapping = aes(x = x, y = y)) +
2016-07-18 21:47:52 +08:00
coord_cartesian(xlim = c(4, 11), ylim = c(4, 11))
```
Why is a scatterplot a better display than a binned plot for this case?
2016-05-23 01:35:52 +08:00
2016-07-18 21:47:52 +08:00
## Patterns and models
Patterns in your data provide clues about relationships. If a systematic relationship exists between two variables it will appear as a pattern in the data. If you spot a pattern, ask yourself:
2016-05-23 01:35:52 +08:00
+ Could this pattern be due to coincidence (i.e. random chance)?
2016-05-23 01:35:52 +08:00
+ How can you describe the relationship implied by the pattern?
2016-05-23 01:35:52 +08:00
+ How strong is the relationship implied by the pattern?
2016-05-23 01:35:52 +08:00
+ What other variables might affect the relationship?
2016-05-23 01:35:52 +08:00
+ Does the relationship change if you look at individual subgroups of the data?
2016-05-23 01:35:52 +08:00
A scatterplot of Old Faithful eruption lengths versus the wait time between eruptions shows a pattern: longer wait times are associated with longer eruptions. The scatterplot also displays the two clusters that we noticed above.
2016-07-18 21:47:52 +08:00
```{r fig.height = 2}
ggplot(data = faithful) +
geom_point(mapping = aes(x = eruptions, y = waiting))
```
2016-05-23 01:35:52 +08:00
2016-07-18 21:47:52 +08:00
Patterns provide one of the most useful tools for data scientists because they reveal covariation. If you think of variation as a phenomenon that creates uncertainty, covariation is a phenomenon that reduces it. If two variables covary, you can use the values of one variable to make better predictions about the values of the second. If the covariation is due to a causal relationship (a special case), then you can use the value of one variable to control the value of the second.
2016-05-23 01:35:52 +08:00
2016-10-07 00:19:55 +08:00
Models are a tool for extracting patterns out of data. For example, consider the diamonds data. It's hard to understand the relationship between cut and price, because cut and carat, and carat and price are tightly related. It's possible to use a model to remove the very strong relationship between price and carat so we can explore the subtleties that remain. The following code fits a model that predicts `price` from `carat` and then computes the residuals (the difference between the predicted value and the actual value). The residuals give us a view of the price of the diamond, once the effect of carat has been removed.
2016-05-23 01:35:52 +08:00
2016-07-18 22:52:55 +08:00
```{r, dev = "png"}
library(modelr)
2016-05-23 01:35:52 +08:00
mod <- lm(log(price) ~ log(carat), data = diamonds)
2016-05-23 01:35:52 +08:00
diamonds2 <- diamonds %>%
add_residuals(mod) %>%
mutate(resid = exp(resid))
2016-05-23 01:35:52 +08:00
ggplot(data = diamonds2) +
geom_point(mapping = aes(x = carat, y = resid))
```
2016-07-23 05:00:11 +08:00
Once you've removed the strong relationship between carat and price, you can see what you expect in the relationship between cut and price: relative to their size, better quality diamonds are more expensive.
```{r}
ggplot(data = diamonds2) +
geom_boxplot(mapping = aes(x = cut, y = resid))
```
2016-07-19 22:39:00 +08:00
2016-10-07 00:19:55 +08:00
You'll learn how models, and the modelr package, work in the final part of the book, [model](#model-intro). We're saving modelling for later because understanding what models are and how they work is easiest once you have tools of data wrangling and programming in hand.
2016-07-20 21:47:17 +08:00
## ggplot2 calls
As we move on from these introductory chapters, we'll transition to a more concise expression of ggplot2 code. So far we've been very explicit, which is helpful when you are learning:
```{r, eval = FALSE}
ggplot(data = faithful, mapping = aes(x = eruptions)) +
2016-07-22 22:15:55 +08:00
geom_freqpoly(binwidth = 0.25)
2016-07-20 21:47:17 +08:00
```
2016-07-23 05:00:11 +08:00
Typically, the first one or two arguments to a function are so important that you should know them by heart. The first two arguments to `ggplot()` are `data` and `mapping`, and the first two arguments to `aes()` are `x` and `y`. In the remainder of the book, we won't supply those names. That saves typing, and, by reducing the amount of boilerplate, makes it easier to see what's different between plots. That's a really important programming concern that we'll come back in [functions].
Rewriting the previous plot more concisely yields:
2016-07-20 21:47:17 +08:00
```{r, eval = FALSE}
ggplot(faithful, aes(eruptions)) +
2016-07-22 22:15:55 +08:00
geom_freqpoly(binwidth = 0.25)
2016-07-20 21:47:17 +08:00
```
2016-08-16 20:22:41 +08:00
Sometimes we'll turn the end of a pipeline of data transformation into a plot. Watch for the transition from `%>%` to `+`. I wish this transition wasn't necessary but unfortunately ggplot2 was created before the pipe was discovered.
2016-07-20 21:47:17 +08:00
```{r, eval = FALSE}
diamonds %>%
count(cut, clarity) %>%
ggplot(aes(clarity, cut, fill = n)) +
2016-08-01 00:32:16 +08:00
geom_tile()
2016-07-20 21:47:17 +08:00
```
2016-07-23 05:00:11 +08:00
## Learning more
2017-05-03 20:50:03 +08:00
If you want to learn more about the mechanics of ggplot2, I'd highly recommend grabbing a copy of the ggplot2 book: <https://amzn.com/331924275X>. It's been recently updated, so it includes dplyr and tidyr code, and has much more space to explore all the facets of visualisation. Unfortunately the book isn't generally available for free, but if you have a connection to a university you can probably get an electronic version for free through SpringerLink.
2016-08-16 21:31:09 +08:00
Another useful resource is the [_R Graphics Cookbook_](https://amzn.com/1449316956) by Winston Chang. Much of the contents are available online at <http://www.cookbook-r.com/Graphs/>.
2016-10-03 21:36:51 +08:00
I also recommend [_Graphical Data Analysis with R_](https://amzn.com/1498715230), by Antony Unwin. This is a book-length treatment similar to the material covered in this chapter, but has the space to go into much greater depth.