-
Notifications
You must be signed in to change notification settings - Fork 25
/
Copy pathbivariate_rf_plots.qmd
152 lines (107 loc) · 7.28 KB
/
bivariate_rf_plots.qmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
---
output: html_document
editor_options:
chunk_output_type: console
---
# Bivariate residual-fit spread plot
```{r echo=FALSE}
source("libs/Common.R")
```
```{r echo = FALSE}
pkg_ver(c("dplyr", "tidyr", "ggplot2", "tukeyedar"))
```
------------------------------------------------------------------------
William Cleveland's *residual-fit spread* (RFS) plot was introduced earlier in this course in the univariate context. The plot provides a visual means of comparing the variability explained by a fitted model—such as group means—to the variability of the residuals. The RFS plot can also be applied to bivariate models offering insight into the relationship between a model's effects and its residual variability.
In bivariate analysis, an important objective is to identify a model that maximizes our ability to explain the variability in $Y$. The greater the variability in $Y$ explained by the model, the smaller the uncertainty in estimating $Y$ for each value of $X$. This objective translates into minimizing the spread of residuals relative to the variability captured by the model.
For example, in the figure below, the plot on the left illustrates a model that effectively explains the variability in $Y$. Notice the relatively small residuals, represented as vertical dashed lines connecting the fitted (predicted) values (red points) to the observed values (grey points), compared to the variability explained by the model. The model accounts for a range in $Y$ of approximately 11 units, while the residual magnitudes average about 1 unit, indicating minimal unexplained variability.
In contrast, the plot on the right shows a model with limited explanatory power. Here, the model accounts for only about 2 units of variability in $Y$, while the residual magnitudes average approximately 4 units. This suggests that a significant portion of the variability in $Y$ remains unexplained by $X$.
```{r echo = FALSE, results='hide', fig.width=8, fig.height=3}
library(tukeyedar)
n <- 20
set.seed(123)
x <- sort(runif(n, 1, 10))
y <- 1 + 1.5 * x + rnorm(n)
#y[20] <- 28
df <- data.frame(x,y)
OP <- par(mfrow = c(1,2))
M1 <- eda_lm(df, x, y, show.par = FALSE, mean.l = FALSE, sd = FALSE)
error_segs <- lapply(1:nrow(df),
FUN = \(i) data.frame(x = df$x[i], y = c(M1$fitted.values[i], df$y[i]),
row.names = NULL))
OP2 <- par(mar = c(0,0,0,0))
lapply(error_segs, \(x) lines(x, col = "grey", lty = 2))
points(x = df$x, y = M1$fitted.values, pch = 16, col = "darkred")
par(OP2)
set.seed(17)
y<- 10 + runif(n, -6, 6)
df <- data.frame(x,y)
M2 <-eda_lm(df,x,y, robust = TRUE, show.par = FALSE, mean.l = FALSE, sd = FALSE)
error_segs <- lapply(1:nrow(df),
FUN = \(i) data.frame(x = df$x[i], y = c(M2$fitted.values[i], df$y[i]),
row.names = NULL))
OP2 <- par(mar = c(0,0,0,0))
lapply(error_segs, \(x) lines(x, col = "grey", lty = 2))
points(x = df$x, y = M2$fitted.values, pch = 16, col = "darkred")
par(OP2)
par(OP)
```
While a scatter plot can provide insights into the model's ability to explain variability in extreme cases, it is not well-suited for assessing performance across a broader range of scenarios. In contrast, an RFS plot is specifically designed to evaluate this aspect making it a more effective tool for such assessments.
## Constructing a residual-fit spread plot
Consider the 1^st^ order polynomial fitted to the miles-per-gallon (`mpg`) vs. horsepower (`hp`) scatter plot from the `mtcars` dataset:
```{r echo=FALSE, fig.height=3, fig.width=3, results='hide'}
M <- eda_lm(mtcars, hp, mpg, show.par = FALSE, mean.l = FALSE, sd = FALSE)
error_segs <- lapply(1:nrow(mtcars),
FUN = \(i) data.frame(x = mtcars$hp[i], y = c(M$fitted.values[i], mtcars$mpg[i]),
row.names = NULL))
OP2 <- par(mar = c(0,0,0,0))
lapply(error_segs, \(x) lines(x, col = "grey", lty = 2))
points(x = mtcars$hp, y = M$fitted.values, pch = 16, col = "darkred")
par(OP2)
```
The estimatefd values ranges from about 7 miles-per-gallon to about 26 miles-per-gallon for an absolute range of around 19 mpg. The residuals are highlighted by the vertical dashed lines.
An RFS plot can help us assess how much of the variability in `mpg` can be explained by the fitted model and how it compares to the residuals.
```{r echo=FALSE, fig.height=3, fig.width=4, results='hide'}
eda_rfs(M)
```
The RFS plot consists of two side-by-side quantile plots. The left plot, known as the fit-minus-mean plot, displays the quantiles of the estimated `mpg` values after subtracting the overall mean (i.e., recentering the estimated `mpg` values at 0). The right plot shows the quantiles of the residuals.
Each point in the RFS plot can be tied back to the scatter plot. For example, the estimated `mpg` value for the largest `hp` value of 335 is 7.24 mpg and its residual is +7.75 mpg. These values are split between the fit-minus-mean quantile plot and the residuals quantile plot.
{width="586"}
By plotting the quantiles of both the estimated values and the residuals, the RFS plot allows for the comparison their respective magnitudes.
The fit-minus-mean plot shows estimates covering a range of approximately 19 mpg, which aligns with the observations from the scatter plot. The residuals span a range of about 15 mpg; however, 75% of the residuals had a magnitude of 10 mpg or less. This suggests that the fitted model makes a meaningful contribution to explaining the variability in mpg.
To build an intuitive understanding of the RFS plot, consider an extreme scenario where the model is reduced to just the mean fitted to the data.
```{r echo=FALSE, fig.height=3, fig.width=3, results='hide'}
M <- eda_lm(mtcars, hp, mpg, show.par = FALSE, poly = 0,
mean.l = FALSE, sd = FALSE)
```
Since the mean is a constant, it explains none of the variability in `mpg`. As a result, the residuals' quantile plot is expected to capture the entire range of variability in `hp` (approximately 25 units), as illustrated in the following RFS plot.
```{r echo=FALSE, fig.height=3, fig.width=4, results='hide'}
eda_rfs(M)
```
The greater the proportion of variability explained by the fitted model relative to the residual spread, the stronger the explanatory power of the model.
## Constructing an RFS plot with `eda_rfs`
The RFS plot can be constructed using the `eda_rfs` function from the `tukeyedar` package. It can take, as input, a model of class `lm`. For example;
```{r class.source="eda", fig.height=3, fig.width=4, results='hide'}
library(tukeyedar)
M <- lm(mpg ~ hp, mtcars)
eda_rfs(M)
```
## Constructing an RFS plot with `ggplot`
Some data manipulation is needed to construct an RFS plot in `ggplot`. Using the regression model `M` from the previous code chunk, the RFS dataset can be constructed as follows:
```{r fig.height=2.8, fig.width=4.5, results='hide'}
library(dplyr)
library(tidyr)
library(ggplot2)
df <- data.frame(Residuals = residuals(M),
`Fit minus mean` = predict(M) - mean(mtcars$mpg),
check.names = FALSE)
rf <- df %>%
pivot_longer(names_to = "type", values_to = "value", cols=everything()) %>%
group_by(type) %>%
arrange(value) %>%
mutate(fval = (row_number() - 0.5) / n())
ggplot(rf, aes(x = fval, y = value)) +
geom_point(alpha = 0.3, cex = 1.5) +
facet_wrap(~ type) +
xlab("f-value") +
ylab("mpg")
```