Managing graphics with R

A very basic plot can be thought of as a succession of layers. In the example below, we first draw a scatterplot using variable lwt (on the x-axis) and bwt (on the y-axis) from the MASS::birthwt dataset, and then add a scatterplot smoother (a lowess curve, with default smoothing parameters) on top of the scatterplot:

p <- ggplot() +
    layer(data = MASS::birthwt,
          stat = "identity",
          geom = "point",
          mapping = aes(x = lwt, y = bwt),
          position = "identity") +
    layer(data = MASS::birthwt,
          stat = "smooth",
          geom = "line",
          mapping = aes(x = lwt, y = bwt),
          position = "identity",
          params = list(method = "auto"))

Here is the simplified version using the approach described above:

library(MASS)
p <- ggplot(data = birthwt, aes(x = lwt, y = bwt))
p + geom_point() + geom_smooth(method = "auto")

This follows the same principled approach: You add each graphical piece together using the + operateur – which is specific to ggplot2 – and you print the final graphical object to the graphical device. This is the reason why we don't use any assignment operateur (<-) in the last expression. It would be possible to write the following expression:

ggplot(data = birthwt, aes(x = lwt, y = bwt)) + geom_point() + geom_smooth(method = "auto")

However, it is better practice to save the graphical commands in a variable, and even to accumulate (+) the instructions as they go along since this allows to build the final plot in an incremental way.

Histogram and density curve are two common estimators for the distribution of a continuous random variable. They do indeed have their counterparts in ggplot2, namely geom_histogram and geom_density, although in the latter case it is also possible to use geom_line using the a density estimator that ggplot2 will compute for us. Here are two examples of use:

p <- ggplot(data = d, aes(x = age)) +
  geom_histogram(binwidth = 5) +
  labs(x = "Age", y = "Counts")
p

p <- ggplot(data = d, aes(x = age)) +
  geom_line(stat = "density") +
  labs(x = "Age", y = "Density")
p

The syntax is almost identical in both case, pending the options specific to each geometrical object: In the case of geom_histogram, it is possible to indicate the size (binwidth) or the number (bins) of intervals to use when discretizing the continuous variable x; with the density estimator, we can specify the size of the smoothing window (bw or adjust) and/or the kernel (kernel).

Boxplots can be used to visually depict the output of the summary function, that is a five-number summary of the distribution of a numerical variable (location, shape and range). It is most useful when there are several distributions to plot side by side, hence the use of x and y in the mapping below:

p <- ggplot(data = d, aes(x = mar, y = logincome)) +
  geom_boxplot() +
  labs(x = "Marital Status", y = "Annual Income (log)") +
  coord_flip()
p

Let us summarize the distribution of average (log) income across socio-economic classes. First, we need to compute the mean and standard deviation of logincome for each level of egp. This is easily performed using aggregate (or tapply):

egp_stats <- aggregate(logincome ~ egp, data = d, mean)
egp_stats$sd <- aggregate(logincome ~ egp, data = d, sd)$logincome
names(egp_stats)[2] <- "mean"
egp_stats

   egp      mean        sd
1 High 10.210198 0.9861072
2  Mid  9.831354 0.9518870
3  Low  9.601578 0.9659040
4 None  8.621649 1.3534442

Unfortunately, aggregate returns a data frame with only one column for the computed results (named after the variable that appears in the LHS), no matter if the function returns a single valued vector or multivariate vectors. Hence, we are forced to extract the result computed with the second call to aggregate and append it as a distinct column. It is then possible to use a combination of geom_bar and geom_errorbar to produce the desired result:

p <- ggplot(data = egp_stats, aes(x = egp, y = mean)) +
  geom_bar(stat = "identity") +
  geom_errorbar(aes(ymin = mean - sd, ymax = mean + sd), width = .2, col = "black") +
  labs(x = "Socio-Economic Status", y = "Average income (log)")
p

Here is another example using two variables in the RHS, egp and sex:

egp_stats <- aggregate(logincome ~ egp, data = d, mean)
egp_stats <- aggregate(logincome ~ egp + sex, data = d, mean)
egp_stats$sd <- aggregate(logincome ~ egp + sex, data = d, sd)$logincome
names(egp_stats)[3] <- "mean"

p <- ggplot(data = egp_stats, aes(x = egp, y = mean, fill = sex)) +
  geom_bar(stat = "identity", position = position_dodge()) +
  geom_errorbar(aes(ymin = mean - sd, ymax = mean + sd), width = .2, col = "black", position = position_dodge(.9)) +
  scale_fill_manual("", values = c("cornflowerblue", "darkorange")) +
  labs(x = "Socio-Economic Status", y = "Average income (log)")
p

Scatterplot are a common way to display the joint variation of two series of observations, specified as x and y coordinates in a 2D plan. A scatterplot smoother, e.g. a lowess (or loess) curve, can be superimposed on the scatterplot to help in gauging local deviation to the linear relationship assumption commonly found in linear regression. In addition to the x and y variables, it is also possible to ask ggplot2 to highlight point using a color that depends on the level of a third categorical variable. This is controlled in the aesthetic directly. When we are not happy with the default color scheme, we can provide our own values (RGB, HEX or named colors) using a "manual" scale, which targets an existing aesthetic. Finally, facet_grid or facet_wrap can be used to build a treillis graphic, meaning that the overall plot is subsetted according to the facet defined in this specific instruction. Here is a full working example in action:

p <- ggplot(data = subset(d, egp != "None"), aes(x = age, y = logincome, color = sex)) +
  geom_point(alpha = 0.5) +
  geom_smooth(method = "loess", se = FALSE, span = 1.5) +
  scale_color_manual("", values = c("cornflowerblue", "darkorange")) +
  facet_grid(mar ~ egp) +
  labs(x = "Age", y = "Average income (log)")
p