A model is a claim about what value to predict.
gf_model() draws that claim on top of the data it was fit
to, so the claim and the evidence are in the same picture.
The empty model
The empty model predicts the same value, the mean, for every observation. On a histogram it is a single vertical line.
empty_model <- lm(body_mass_kg ~ NULL, data = penguins)
gf_histogram(~body_mass_kg, data = penguins, binwidth = 0.25) %>%
gf_model(empty_model)
A plot already implies a model, so you can leave the argument out and
let gf_model() work out which one. A distribution on its
own implies the empty model.
gf_histogram(~body_mass_kg, data = penguins, binwidth = 0.25) %>%
gf_model()
Group models
With a categorical explanatory variable the model predicts one value
per group, and gf_model() draws a mark at each group
mean.
species_model <- lm(body_mass_kg ~ species, data = penguins)
gf_jitter(body_mass_kg ~ species, data = penguins, width = .1) %>%
gf_model(species_model)
The same model layers onto faceted histograms, which makes the comparison between groups a comparison between panels.
gf_histogram(~body_mass_kg, data = penguins, binwidth = 0.25) %>%
gf_facet_grid(species ~ .) %>%
gf_model(species_model)
Regression models
With a quantitative explanatory variable the prediction is a line.
flipper_model <- lm(body_mass_kg ~ flipper_length_m, data = penguins)
gf_point(body_mass_kg ~ flipper_length_m, data = penguins) %>%
gf_model(flipper_model)
Layering two models in different colors puts the comparison on one plot: the empty model ignores flipper length, the regression model uses it.
gf_point(body_mass_kg ~ flipper_length_m, data = penguins) %>%
gf_model(empty_model, color = "dodgerblue") %>%
gf_model(flipper_model, color = "firebrick")
Residuals
A residual is the distance from what the model predicted to what
actually happened. gf_resid() draws them as segments. With
hundreds of points this is unreadable, so take a small sample first.
set.seed(1)
penguins_20 <- sample(penguins, 20)
sample_model <- lm(body_mass_kg ~ flipper_length_m, data = penguins_20)
gf_point(body_mass_kg ~ flipper_length_m, data = penguins_20) %>%
gf_model(sample_model) %>%
gf_resid(sample_model, color = "firebrick")
gf_square_resid() draws each residual as a square
instead of a segment, which makes squared error visible as
area: the quantity least squares is actually minimizing.
gf_point(body_mass_kg ~ flipper_length_m, data = penguins_20) %>%
gf_model(sample_model) %>%
gf_square_resid(sample_model, color = "firebrick")
Comparing the empty model’s squares to the regression model’s is the clearest picture of what the predictor bought you.
empty_20 <- lm(body_mass_kg ~ NULL, data = penguins_20)
gf_point(body_mass_kg ~ flipper_length_m, data = penguins_20) %>%
gf_model(empty_20) %>%
gf_square_resid(empty_20, color = "dodgerblue")
Where the error goes
The blue squares are the error the empty model leaves. The red ones
are what is left once flipper length is accounted for. The difference
between them is the part the predictor explains, and
gf_reduce() draws it: for each penguin, the distance from
the grand mean to the model’s prediction.
gf_point(body_mass_kg ~ flipper_length_m, data = penguins_20) %>%
gf_model(empty_20) %>%
gf_model(sample_model) %>%
gf_reduce(sample_model, color = "forestgreen")
Drawn as squares, the three areas complete each other.
gf_point(body_mass_kg ~ flipper_length_m, data = penguins_20) %>%
gf_model(sample_model) %>%
gf_square_resid(empty_20, color = "dodgerblue") %>%
gf_square_resid(sample_model, color = "firebrick") %>%
gf_square_reduce(sample_model, color = "forestgreen")
Blue is red plus green: the error you started with is the error you
are left with plus the error the predictor explained. The green share of
the blue is PRE, the same number supernova() puts in a
table.
supernova(sample_model)
#> Analysis of Variance Table (Type III SS)
#> Model: body_mass_kg ~ flipper_length_m
#>
#> SS df MS F PRE p
#> ----- --------------- | ------ -- ----- ------ ----- -----
#> Model (error reduced) | 6.131 1 6.131 24.841 .5798 .0001
#> Error (from model) | 4.442 18 0.247
#> ----- --------------- | ------ -- ----- ------ ----- -----
#> Total (empty model) | 10.573 19 0.556The standard deviation as a typical residual
gf_sd_ruler() draws one standard deviation, anchored at
the mean. Under the empty model the residuals are deviations
from the mean, so the ruler is a picture of a typical one.
gf_point(Thumb ~ Height, data = Fingers, alpha = .4) %>%
gf_model(lm(Thumb ~ NULL, data = Fingers)) %>%
gf_sd_ruler(where = "mean")
On a histogram the outcome is on the x-axis, so the ruler turns horizontal and runs along the baseline from the mean to one SD above it.
gf_histogram(~Thumb, data = Fingers, binwidth = 5) %>%
gf_model(lm(Thumb ~ NULL, data = Fingers)) %>%
gf_sd_ruler(color = "red", linewidth = 2)
Two groups with the same mean and different spread get rulers of
different lengths, which is the comparison the ruler exists to make.
Zooming both to the same x range keeps them comparable.
coord_cartesian() does that by changing what you see;
gf_lims() would do it by discarding anything outside the
range, which is not what you want when the whole point is to compare
distributions.
set.seed(154)
no_feedback <- data.frame(median_time = round(rnorm(100, 13, 6), 1))
set.seed(141)
feedback <- data.frame(median_time = round(rnorm(100, 13, 3), 1))
gf_histogram(~median_time, data = no_feedback, binwidth = 1) %>%
gf_refine(coord_cartesian(xlim = c(0, 31))) %>%
gf_model(lm(median_time ~ NULL, data = no_feedback)) %>%
gf_sd_ruler(color = "red", linewidth = 2)
gf_histogram(~median_time, data = feedback, binwidth = 1) %>%
gf_refine(coord_cartesian(xlim = c(0, 31))) %>%
gf_model(lm(median_time ~ NULL, data = feedback)) %>%
gf_sd_ruler(color = "red", linewidth = 2)