Chapter 5.2.1 Exercises

Untitled
In [2]:
# Some meaningless calculations. Not important

a <- (100 + 3) - 2
mean(c(a / 100, 642564624.34))

# t.test comparing revenue of sequels v non-sequels

t.test(formula = revenue.all ~ sequel,
       data = movies)

# A scatterplot of budget and dvd revenue. 
#  Hard to see a relationship

plot(x = movies$budget,
     y = movies$dvd.usa,
     main = "myplot")
321282312.675
	Welch Two Sample t-test

data:  revenue.all by sequel
t = -10.876, df = 614.22, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -149.0457 -103.4524
sample estimates:
mean in group 0 mean in group 1 
       83.63355       209.88262 
In [5]:
library(ggplot2)
library(dplyr)

# Missing values: 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.
# Drop the entire row with the strange values:

    diamonds2 <- diamonds %>% 
      filter(between(y, 3, 20))

# I don’t recommend this option because just because one measurement is invalid, doesn’t mean all the measurements are. Additionally, if you 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!
In [6]:
# Instead, I recommend replacing the unusual values with missing values. The easiest way to do this is to use mutate() to replace the variable with a modified copy. You can use the ifelse() function to replace unusual values with NA:

    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.
In [7]:
# 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:

ggplot(data = diamonds2, mapping = aes(x = x, y = y)) + 
  geom_point()
Warning message:
“Removed 9 rows containing missing values (geom_point).”
In [8]:
# To suppress that warning, set na.rm = TRUE:

ggplot(data = diamonds2, mapping = aes(x = x, y = y)) + 
  geom_point(na.rm = TRUE)
In [10]:
# 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().

nycflights13::flights %>% 
  mutate(
    cancelled = is.na(dep_time),
    sched_hour = sched_dep_time %/% 100,
    sched_min = sched_dep_time %% 100,
    sched_dep_time = sched_hour + sched_min / 60
  ) %>% 
  ggplot(mapping = aes(sched_dep_time)) + 
    geom_freqpoly(mapping = aes(colour = cancelled), binwidth = 1/4)
In [2]:
library(ggplot2)

# 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:

ggplot(data = diamonds) +
  geom_bar(mapping = aes(x = cut))
In [3]:
# The height of the bars displays how many observations occurred with each x value. You can compute these values manually with dplyr::count():

library(dplyr)
diamonds %>% 
  count(cut)
cutn
Fair 1610
Good 4906
Very Good12082
Premium 13791
Ideal 21551
In [3]:
# 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:

ggplot(data = diamonds) +
  geom_histogram(mapping = aes(x = carat), binwidth = 0.5)
In [4]:
# You can compute this by hand by combining dplyr::count() and ggplot2::cut_width():

diamonds %>% 
  count(cut_width(carat, 0.5))
cut_width(carat, 0.5)n
[-0.25,0.25] 785
(0.25,0.75] 29498
(0.75,1.25] 15977
(1.25,1.75] 5313
(1.75,2.25] 2002
(2.25,2.75] 322
(2.75,3.25] 32
(3.25,3.75] 5
(3.75,4.25] 4
(4.25,4.75] 1
(4.75,5.25] 1
In [4]:
# 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.

smaller <- diamonds %>% 
  filter(carat < 3)
  
ggplot(data = smaller, mapping = aes(x = carat)) +
  geom_histogram(binwidth = 0.1)
In [6]:
# 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.

ggplot(data = smaller, mapping = aes(x = carat, colour = cut)) +
  geom_freqpoly(binwidth = 0.1)

# There are a few challenges with this type of plot, which we will come back to in visualising a categorical and a continuous variable.
In [7]:
# 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?

ggplot(data = smaller, mapping = aes(x = carat)) +
  geom_histogram(binwidth = 0.01)
In [8]:
# 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.

ggplot(data = faithful, mapping = aes(x = eruptions)) + 
  geom_histogram(binwidth = 0.25)

# 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.
In [9]:
ggplot(diamonds) + 
  geom_histogram(mapping = aes(x = y), binwidth = 0.5)
In [12]:
# 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():

ggplot(diamonds) + 
  geom_histogram(mapping = aes(x = y), binwidth = 0.5) +
  coord_cartesian(ylim = c(0, 50))
In [15]:
# (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:
diamonds
unusual <- diamonds %>% 
  filter(y < 3 | y > 20) %>% 
  select(price, x, y, z) %>%
  arrange(y)
unusual
caratcutcolorclaritydepthtablepricexyz
0.23 Ideal E SI2 61.5 55 326 3.95 3.98 2.43
0.21 Premium E SI1 59.8 61 326 3.89 3.84 2.31
0.23 Good E VS1 56.9 65 327 4.05 4.07 2.31
0.29 Premium I VS2 62.4 58 334 4.20 4.23 2.63
0.31 Good J SI2 63.3 58 335 4.34 4.35 2.75
0.24 Very GoodJ VVS2 62.8 57 336 3.94 3.96 2.48
0.24 Very GoodI VVS1 62.3 57 336 3.95 3.98 2.47
0.26 Very GoodH SI1 61.9 55 337 4.07 4.11 2.53
0.22 Fair E VS2 65.1 61 337 3.87 3.78 2.49
0.23 Very GoodH VS1 59.4 61 338 4.00 4.05 2.39
0.30 Good J SI1 64.0 55 339 4.25 4.28 2.73
0.23 Ideal J VS1 62.8 56 340 3.93 3.90 2.46
0.22 Premium F SI1 60.4 61 342 3.88 3.84 2.33
0.31 Ideal J SI2 62.2 54 344 4.35 4.37 2.71
0.20 Premium E SI2 60.2 62 345 3.79 3.75 2.27
0.32 Premium E I1 60.9 58 345 4.38 4.42 2.68
0.30 Ideal I SI2 62.0 54 348 4.31 4.34 2.68
0.30 Good J SI1 63.4 54 351 4.23 4.29 2.70
0.30 Good J SI1 63.8 56 351 4.23 4.26 2.71
0.30 Very GoodJ SI1 62.7 59 351 4.21 4.27 2.66
0.30 Good I SI2 63.3 56 351 4.26 4.30 2.71
0.23 Very GoodE VS2 63.8 55 352 3.85 3.92 2.48
0.23 Very GoodH VS1 61.0 57 353 3.94 3.96 2.41
0.31 Very GoodJ SI1 59.4 62 353 4.39 4.43 2.62
0.31 Very GoodJ SI1 58.1 62 353 4.44 4.47 2.59
0.23 Very GoodG VVS2 60.4 58 354 3.97 4.01 2.41
0.24 Premium I VS1 62.5 57 355 3.97 3.94 2.47
0.30 Very GoodJ VS2 62.2 57 357 4.28 4.30 2.67
0.23 Very GoodD VS2 60.5 61 357 3.96 3.97 2.40
0.23 Very GoodF VS1 60.9 57 357 3.96 3.99 2.42
0.70 Premium E SI1 60.5 58 2753 5.74 5.77 3.48
0.57 Premium E IF 59.8 60 2753 5.43 5.38 3.23
0.61 Premium F VVS1 61.8 59 2753 5.48 5.40 3.36
0.80 Good G VS2 64.2 58 2753 5.84 5.81 3.74
0.84 Good I VS1 63.7 59 2753 5.94 5.90 3.77
0.77 Ideal E SI2 62.1 56 2753 5.84 5.86 3.63
0.74 Good D SI1 63.1 59 2753 5.71 5.74 3.61
0.90 Very GoodJ SI1 63.2 60 2753 6.12 6.09 3.86
0.76 Premium I VS1 59.3 62 2753 5.93 5.85 3.49
0.76 Ideal I VVS1 62.2 55 2753 5.89 5.87 3.66
0.70 Very GoodE VS2 62.4 60 2755 5.57 5.61 3.49
0.70 Very GoodE VS2 62.8 60 2755 5.59 5.65 3.53
0.70 Very GoodD VS1 63.1 59 2755 5.67 5.58 3.55
0.73 Ideal I VS2 61.3 56 2756 5.80 5.84 3.57
0.73 Ideal I VS2 61.6 55 2756 5.82 5.84 3.59
0.79 Ideal I SI1 61.6 56 2756 5.95 5.97 3.67
0.71 Ideal E SI1 61.9 56 2756 5.71 5.73 3.54
0.79 Good F SI1 58.1 59 2756 6.06 6.13 3.54
0.79 Premium E SI2 61.4 58 2756 6.03 5.96 3.68
0.71 Ideal G VS1 61.4 56 2756 5.76 5.73 3.53
0.71 Premium E SI1 60.5 55 2756 5.79 5.74 3.49
0.71 Premium F SI1 59.8 62 2756 5.74 5.73 3.43
0.70 Very GoodE VS2 60.5 59 2757 5.71 5.76 3.47
0.70 Very GoodE VS2 61.2 59 2757 5.69 5.72 3.49
0.72 Premium D SI1 62.7 59 2757 5.69 5.73 3.58
0.72 Ideal D SI1 60.8 57 2757 5.75 5.76 3.50
0.72 Good D SI1 63.1 55 2757 5.69 5.75 3.61
0.70 Very GoodD SI1 62.8 60 2757 5.66 5.68 3.56
0.86 Premium H SI2 61.0 58 2757 6.15 6.12 3.74
0.75 Ideal D SI2 62.2 55 2757 5.83 5.87 3.64
pricexyz
51390.00 0.0 0.00
63810.00 0.0 0.00
128000.00 0.0 0.00
156860.00 0.0 0.00
180340.00 0.0 0.00
21300.00 0.0 0.00
21300.00 0.0 0.00
20755.15 31.8 5.12
122108.09 58.9 8.06
In [16]:
#  Covariation - a categorical and continuous variable; let’s explore how the price of a diamond varies with its quality.

ggplot(data = diamonds, mapping = aes(x = price)) + 
  geom_freqpoly(mapping = aes(colour = cut), binwidth = 500)
In [17]:
# It’s hard to see the difference in distribution because the overall counts differ so much:

ggplot(diamonds) + 
  geom_bar(mapping = aes(x = cut))
In [18]:
# 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.

ggplot(data = diamonds, mapping = aes(x = price, y = ..density..)) + 
  geom_freqpoly(mapping = aes(colour = cut), binwidth = 500)
In [5]:
# Let’s take a look at the distribution of price by cut using geom_boxplot():

ggplot(data = diamonds, mapping = aes(x = cut, y = price)) +
  geom_boxplot()
In [20]:
# For example, take the class variable in the mpg dataset. You might be interested to know how highway mileage varies across classes:

ggplot(data = mpg, mapping = aes(x = class, y = hwy)) +
  geom_boxplot()
In [21]:
# To make the trend easier to see, we can reorder class based on the median value of hwy:

ggplot(data = mpg) +
  geom_boxplot(mapping = aes(x = reorder(class, hwy, FUN = median), y = hwy))
In [22]:
# If you have long variable names, geom_boxplot() will work better if you flip it 90°. You can do that with coord_flip().

ggplot(data = mpg) +
  geom_boxplot(mapping = aes(x = reorder(class, hwy, FUN = median), y = hwy)) +
  coord_flip()
In [23]:
# Two categorical variables: 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():

ggplot(data = diamonds) +
  geom_count(mapping = aes(x = cut, y = color))
In [26]:
# Then visualise with geom_tile() and the fill aesthetic:

diamonds %>% 
  count(color, cut) %>%  
  ggplot(mapping = aes(x = color, y = cut)) +
    geom_tile(mapping = aes(fill = n))

# 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.
In [27]:
# Two continuous variables: 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.

ggplot(data = diamonds) +
  geom_point(mapping = aes(x = carat, y = price))
In [28]:
# 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.

ggplot(data = diamonds) + 
  geom_point(mapping = aes(x = carat, y = price), alpha = 1 / 100)

# 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.
In [29]:
# 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().

ggplot(data = smaller) +
  geom_bin2d(mapping = aes(x = carat, y = price))

# install.packages("hexbin")
ggplot(data = smaller) +
  geom_hex(mapping = aes(x = carat, y = price))
In [30]:
# 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:

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.
In [32]:
# Another approach is to display approximately the same number of points in each bin. That’s the job of cut_number():

ggplot(data = smaller, mapping = aes(x = carat, y = price)) + 
  geom_boxplot(mapping = aes(group = cut_number(carat, 20)))
In [33]:
ggplot(data = faithful) + 
  geom_point(mapping = aes(x = eruptions, y = waiting))
In [34]:
# 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.

library(modelr)

mod <- lm(log(price) ~ log(carat), data = diamonds)

diamonds2 <- diamonds %>% 
  add_residuals(mod) %>% 
  mutate(resid = exp(resid))

ggplot(data = diamonds2) + 
  geom_point(mapping = aes(x = carat, y = resid))
In [35]:
# 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.

ggplot(data = diamonds2) + 
  geom_boxplot(mapping = aes(x = cut, y = resid))

# You’ll learn how models, and the modelr package, work in the final part of the book, model. 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.
In [36]:
# 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:

ggplot(data = faithful, mapping = aes(x = eruptions)) + geom_freqpoly(binwidth = 0.25)
In [37]:
# Rewriting the previous plot more concisely yields:

ggplot(faithful, aes(eruptions)) + 
  geom_freqpoly(binwidth = 0.25)