The empirical QQ plot

Introduction

The empirical quantile-quantile plot (QQ plot) is probably one of the most underused and least appreciated plots in univariate analysis. It is used to compare two distributions across their full range of values. It is a generalization of the boxplot in that it does not limit the comparison to just the median and upper and lower quartiles. In fact, it compares all values by matching each value in one batch to its corresponding quantile in the other batch. The sizes of each batch need not be the same. If they differ, the larger batch is interpolated down to the smaller batch’s set of quantiles.

A QQ plot does not only help visualize the differences in distributions, but it can also model the relationship between both batches. Note that this is not to be confused with modeling the relationship between a bivariate dataset where the latter pairs up the points by observational units whereas a QQ plot pairs up the values by their matching quantiles.

Anatomy of the `eda_qq` plot

#>  [1] "Suggested offsets:y = x * 1.4573 + (0.9914)"

Each point represents matching quantiles from each batch.
The shaded boxes represent each batch’s interquartile range (mid 50% of values).
The Solid dashed lines inside the shaded boxes represent each batch’s medians.
The lightly shaded dashed dots represent each batch’s 12.5^th and 87.5^th quantiles (i.e. they show the ends of the mid 80% of values).
The upper right-hand text indicates the power transformation applied to the both batches (default is a power of 1 which is original measurement scale). If a formula is applied to one or both batches, it too will appear in the upper right-hand text.
The eda_qq will also output the suggested relationship between the y variable and the x variable in the console. It bases this on each batch’s interquartile values.
If the output is assigned to a new object, that object will store a list of the following values: the x value (interpolated if needed), the y value (interpolated if needed), the power parameter, the formula applied to x variable, and the formula applied to y variable.

An overview of some of the function arguments

Data type

The function will accept a dataframe with a values column and a group column, or it will accept two separate vector objects. For example, to pass two separate vector object, x and y, type:

library(tukeyedar)
set.seed(207)
x <- rnorm(30)
y <- rnorm(30) + 0.5
eda_qq(x, y)

If the data are in a dataframe, type:

dat <- data.frame(val = c(x, y), cat = rep(c("x", "y"), each = 30))
eda_qq(dat, val, cat)

Suppressing the plot

You can suppress the plot and have the x and y values outputted to a list. If the batches did not match in size, the output will show their interpolated values such that the output batches match in size.

The output will also include the power parameter applied to both batches as well as any formula applied to one or both batches (fx is the formula applied to the x variable and fy is the formula applied to the y variable).

out <- eda_qq(x,y, plot = FALSE)
#>  [1] "Suggested offsets:y = x * 1.12 + (0.5619)"
out
#>  $x
#>   [1] -2.0207122 -1.6048333 -1.5620907 -1.5128732 -1.3126378 -1.1770882
#>   [7] -1.0871906 -0.9258832 -0.8896555 -0.6152073 -0.3140113 -0.2996734
#>  [13] -0.2954234 -0.2199849 -0.2108781  0.1202060  0.2608893  0.2680445
#>  [19]  0.2910663  0.4239690  0.4262605  0.4301416  0.5176361  0.6085180
#>  [25]  0.6880919  0.6929772  0.7640838  0.9037644  1.0124869  1.0503544
#>  
#>  $y
#>   [1] -2.10710669 -1.30465821 -1.17618932 -1.17253191 -0.86423268 -0.40162859
#>   [7] -0.37002087 -0.19629536 -0.07210822 -0.01829722  0.05826287  0.09105884
#>  [13]  0.09371398  0.13698992  0.22318119  0.43006689  0.52597363  0.72665767
#>  [19]  0.81351407  1.00612388  1.01831440  1.06713353  1.23708449  1.24360530
#>  [25]  1.33232007  1.43973056  1.59125312  1.64852115  1.77464625  1.93823390
#>  
#>  $p
#>  [1] 1
#>  
#>  $fx
#>  NULL
#>  
#>  $fy
#>  NULL

Setting the grey box and dashed line parameters

The grey box highlights the interquartile ranges for both batches. Its boundary can be modified via the b.val argument. Likewise, the lightly shaded dashed dots that highlight the mid 80% of values can be modified via the l.val argument.

For example, to highlight the mid 68% of values using the grey boxes and the mid 95% of values using the lightly shaded dashed dots, type:

eda_qq(x, y, b.val = c(0.16, 0.84), l.val = c(0.025, 0.975))

You can suppress the plotting of the grey box and the lightly shaded dashed dots by setting q = FALSE. This does not affect the median dashed lines.

Applying a formula to one of the batches

You can apply a formula to a batch via the fx argument for the x-variable and the fy argument for the y-variable. The formula is passed as a text string. For example, to add 0.5 to the x values, type:

eda_qq(x, y, fx = "x + 0.5")

Quantile type

There are many different quantile algorithms available in R. To see the full list of quantile types, refer to the quantile help page: ?quantile. By default, eda_qq() adopts q.type = 5. In general, the choice of quantiles will not really matter, especially for large datasets. If you want to adopt R’s default type, set q.type = 7.

Point symbols

The point symbol type, color and size can be modified via the pch, p.col (and/or p-fill) and size arguments. The color can be either a built-in color name (you can see the full list by typing colors()) or the rgb() function. If you define the color using one of the built-in color names, you can adjust its transparency via the alpha argument. An alpha value of 0 renders the point completely transparent and a value of 1 renders the point completely opaque. The point symbol can take on two color parameters depending on point type. If pch is any number between 21 and 25, p.fill will define its fill color and p.col will define its border color. For any other point symbol type, the p.fill argument is ignored.

Here are a few examples:

eda_qq(x, y, p.fill = "bisque", p.col = "red", size = 1.2)

eda_qq(x, y, pch = 16, p.col = "tomato2", size = 1.5, alpha = 0.5)

eda_qq(x, y, pch = 3, p.col = "tomato2", size = 1.5)

Interpreting a QQ plot

To help interpret the following QQ plots, we’ll compare each plot to matching kernel density plots.

Identical distributions

In this first example, we will generate a QQ plot of two identical distributions.

library(tukeyedar)
set.seed(543)
x <- rnorm(100)
y <- x
eda_qq(x, y)
eda_dens(x, y)

When two distributions are identical, the points will line up along the x=y line as shown above. This will also generate overlapping density plots as seen on the right plot.

Additive offset

We will work off of the same batches, but this time we will offset the second batch, y, by 2. This case is referred to as an additive offset.

library(tukeyedar)
set.seed(543)
x <- rnorm(100)
y <- x + 2
eda_qq(x, y)
eda_dens(x, y)

You’ll note that the points are parallel to the x=y line. This indicates that the distributions have the exact same shape. But, they do not fall on the x=y line–they are offset by +2 units when measured along the y-axis, as expected. We can confirm this by adding 2 to the x batch:

eda_qq(x, y, fx = "x + 2")
eda_dens(x, y, fx = "x + 2")

The points overlap the x=y line perfectly. The density distribution overlap exactly as well.

Multiplicative offset

Next, we explore two batches that share the same central value but where the second batch is 0.5 times that of the first. This case is referred to as a multiplicative offset.

y <- x * 0.5
eda_qq(x, y)
eda_dens(x, y)

Here, the series of points are at an angle of the x=y line yet, you’ll note that the points follow a perfectly straight line. This suggests a multiplicative offset with no change in location. This indicates that the “shape” between the batches are similar, but that one is “wider” than the other. Here, y is half as wide as x. We can also state that x is twice as wide as y

We know the multiplicative offset since we synthetically generated the values x and y. But in practice eyeballing the multiplier from the plot is not straightforward. We can use the suggested offset of 0.5 displayed in the console to help guide us. We can also use the angle between the points and the x=y line to judge the direction to take when choosing a multiplier. If the points make up an angle that is less than the x=y line, we want to choose an x multiplier that is less than 1. If the angle is greater than the x=y line, we want to choose a multiplier that is greater than 1.

Here, we know that the multiplier is 0.5. Let’s confirm this with the following code chunk:

eda_qq(x, y, fx = "x * 0.5")
eda_dens(x, y, fx = "x * 0.5")

Both additive and multiplicative offset

In this next example, we add both a multiplicative and an additive offset to the data.

y <- x * 0.5 + 2
eda_qq(x, y)
eda_dens(x, y)

We now see both a multiplicative offset (the points form an angle to the x=y line) and an additive offset (the points do not intersect the x=y line). This suggests that the width of the batches differ and that their values are offset by a constant value across the full range of values. It’s usually best to first identify the multiplicative offset such that the points are rendered parallel to the x=y line. Once the multiplier is identified, we can then identify the additive offset.

eda_qq(x, y, fx = "x * 0.5")
eda_dens(x, y, fx = "x * 0.5")

No surprise, the multiplier of 0.5 renders the series of points parallel to the x=y line. We can now eyeball the offset of +2 when measuring the distance between any of the points and the x=y line when measured along the y-axis.

eda_qq(x, y, fx = "x * 0.5 + 2")
eda_dens(x, y, fx = "x * 0.5 + 2")

We can model the relationship between y and x as $y = x * 0.5 + 2$.

Batches need not be symmetrical

So far, we have worked with a normally distributed dataset. But note that any distribution can be used in a QQ plot. For example, two equally skewed distributions that differ by their central values will generate points that are perfectly lined up.

set.seed(540)
x2 <- rbeta(100, 1, 8)
y2 <- x2 + 0.2 
eda_qq(x2, y2)
eda_dens(x2, y2)

Since both distributions have an identical shape, the points will follow a straight line regardless of the nature of the shape (skewed, unimodal, bimodal, etc…)

Perfectly alligned points are rare

Note that with observational data, two batches pulled from the same underlying population will seldom follow a perfectly straight line. For example, note the meandering pattern generated by the following QQ plot when both batches are pulled from a same Normal distribution.

set.seed(702)
x <- rnorm(100)
y <- rnorm(100)
eda_qq(x, y)

Note how the points not only meander about the x=y line, but some of the points tail off near the end of the distribution.

Power transformation

The eda_qq function allows you to apply a power transformation to both batches. Note that transforming just one batch makes little sense since you end up comparing two batches measured on a different scale. In this example, we make use of R’s Indometh dataset where we compare indometacin plasma concentrations between two test subjects.

s1 <- subset(Indometh, Subject == 1, select = conc, drop = TRUE) # Test subject 1
s2 <- subset(Indometh, Subject == 2, select = conc, drop = TRUE) # Test subject 2

eda_qq(s1, s2)

The QQ plot’s grey boxes are shifted towards lower values. The grey boxes show the mid 50% of values for each batch (aka, the IQR range). When the IQR is shifted towards lower values or higher values, this suggests skewed data. Another telltale sign of a skewed dataset is the gradual dispersion of the points in any of the two diagonal directions. Here, we go from a relatively high density of points near the lower values and a lower density of points near higher values. This does not preclude us from identifying multiplicative/additive offsets, but many statistical procedures benefit from a symmetrical distribution. Given that the values are measures of concentration, we might want to adopt a log transformation. The power transformation is defined with the p argument (a power parameter value of 0 defines a log transformation).

eda_qq(s1, s2, p = 0)

#>  [1] "Suggested offsets:y = x * 0.8508 + (0.3914)"

This seems to do a decent job of symmetrizing both distributions.

Note that the suggested offset displayed in the console applies to the transformed dataset. We can verify this by applying the offset to the x batch.

eda_qq(s1, s2, p = 0, fx = "x * 0.8501 + 0.3902")

We can characterize the differences in indometacin plasma concentrations between subject 1 and subject 2 as $log(conc)_{s2} = log(conc)_{s1} * 0.8501 + 0.3902$

The Tukey mean-difference plot

The Tukey mean-difference plot is simply an extension of the QQ plot whereby the plot is rotated such that the x=y line becomes horizontal. This can be useful in helping identify subtle differences between the point pattern and the line. The plot is rotated 45° by mapping the difference between both batches to the y-axis, and mapping the mean between both batches to the x-axis. For example, the following figure on the left (the QQ plot) shows the additive offset between both batches but it fails to clearly identify the multiplicative offset. The latter can be clearly seen in the Tukey mean-difference plot (on the right) which is invoked by setting the argument md = TRUE.

y <- x * 0.97 + 0.3
eda_qq(x, y, title = "QQ plot")
eda_qq(x, y, md = TRUE, title = "M-D plot")

A working example

Many datasets will have distributions that differ not only additively and/or multiplicatively, but also by their general shape. This may create complex point patterns in your QQ plot. Such cases can be indicative of different processes at play for different ranges of values. We’ll explore such a case using wat95 and wat05 dataframes available in this package. The data represent derived normal temperatures for the 1981-2010 period (wat95) and the 1991-2020 period (wat05) for the city of Waterville, Maine (USA). We will subset the data to the daily average normals, avg:

old <- wat95$avg # legacy temperature normals
new <- wat05$avg # current temperature normals

We will now compare their distributions.

out <- eda_qq(old, new)

At first glance, the batches do not seem to differ. But how close are the points to the x=y line? We will rotate the plot to zoom in on the x=y line by setting md = TRUE.

out <- eda_qq(old, new, md = TRUE)

This view is proving to be far more insightful. Note that for most of the data, we have an overall offset of 0.5 suggesting that the new normals are about 0.5°F. warmer. We could stop there, but the pattern observed in the Tukey mean-difference plot is far from random. In fact, we can break the pattern into three distinct ones at around 35°F and 50°F. We will categorize these groups as low, mid and high values.

This next chunk of code will generate three separate QQ plots for each range of values.

labs <- c("low", "mid", "high")
out$avg <- (out$x + out$y) / 2
out <- as.data.frame(out[c(1:2,6)])
out2 <- split(out, cut(out$avg, c(min(out$avg), 35, 50, max(out$avg)), 
                       labels = labs, include.lowest = TRUE))
sapply(labs, FUN = \(x) {eda_qq(out2[[x]]$x, out2[[x]]$y ,
                                xlab = "old", ylab = "new", md = T)
                         title(x, line = 3, col.main = "orange")} )

#>  [1] "Suggested offsets:y = x * 0.9507 + (1.6599)"
#>  [1] "Suggested offsets:y = x * 1.047 + (-1.6502)"
#>  [1] "Suggested offsets:y = x * 0.9938 + (0.9268)"

The suggested offsets are displayed in the console for each group.

It is important to note that we are splitting the paired values generated from the earlier eda_qq function (out$avg) and not the values from the original datasets (old and new). Had we split the original data prior to combining them with eda_qq, we could have generated different patterns in the QQ and Tukey plots. This is because this approach could generate different quantile pairs.

Next, we’ll adopt the suggested offsets generated in the console.

xform <- c("x * 0.9506 + 1.661",
           "x * 1.0469 - 1.6469",
           "x * 0.9938 + 0.9268")
names(xform) <- labs
sapply(labs, FUN = \(x) {eda_qq(out2[[x]]$x, out2[[x]]$y, 
                                fx = xform[x],
                                xlab = "old", ylab = "new", md = T)
                         title(x, line = 3, col.main = "coral3")}
)

The proposed offsets seem to do a good job in characterizing the differences in temperatures.

The characterization of the differences in normal temperatures between the old and new set of normals can be formalized as follows:

\[ new = \begin{cases} old * 0.9506 + 1.661, & T_{avg} < 35 \\ old * 1.0469 - 1.6469, & 35 \le T_{avg} < 50 \\ old * 0.9938 + 0.9268, & T_{avg} \ge 50 \end{cases} \]

The key takeaways from this analysis can be summarized as follows:

Overall, the new normals are about 0.5°F warmer.
The offset is not uniform across the full range of temperature values.
For lower temperature (those less than 35°F), the new normals have a slightly narrower distribution and are about 1.7°F warmer.
For the mid temperature values, (between 35°F and 50°F the new normals have a slightly wider distribution and are overall about 1.6°F cooler.
For the higher range of temperature values (those greater than 50°F), the new normals have a slightly narrower distribution and are about 0.9°F warmer.

Manuel Gimond

2024-03-18

Introduction

Anatomy of the `eda_qq` plot

An overview of some of the function arguments

Data type

Suppressing the plot

Setting the grey box and dashed line parameters

Applying a formula to one of the batches

Quantile type

Point symbols

Interpreting a QQ plot

Identical distributions

Additive offset

Multiplicative offset

Both additive and multiplicative offset

Batches need not be symmetrical

Perfectly alligned points are rare

Power transformation

The Tukey mean-difference plot

A working example

The empirical QQ plot

Manuel Gimond

2024-03-18

Introduction

Anatomy of the eda_qq plot

An overview of some of the function arguments

Data type

Suppressing the plot

Setting the grey box and dashed line parameters

Applying a formula to one of the batches

Quantile type

Point symbols

Interpreting a QQ plot

Identical distributions

Additive offset

Multiplicative offset

Both additive and multiplicative offset

Batches need not be symmetrical

Perfectly alligned points are rare

Power transformation

The Tukey mean-difference plot

A working example

Anatomy of the `eda_qq` plot