- Linear regression shows how a numeric outcome changes with one or more predictors.
- We’ll keep it simple: use the built-in mtcars data.
- Goal: fit \(\text{mpg} = \beta_0 + \beta_1\,\text{hp} + \epsilon\) and interpret.
- We’ll include two ggplot2 visuals, one plotly 3D plot, two math slides, and a code slide.
Why this topic?
Data: mtcars
- 32 cars from 1973–74 Motor Trend US magazine
- Variables: miles-per-gallon (mpg), horsepower (hp), weight (wt), cylinders (cyl), etc.
library(ggplot2) library(dplyr)
## ## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats': ## ## filter, lag
## The following objects are masked from 'package:base': ## ## intersect, setdiff, setequal, union
library(plotly)
## ## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2': ## ## last_plot
## The following object is masked from 'package:stats': ## ## filter
## The following object is masked from 'package:graphics': ## ## layout
data(mtcars)
df <- mtcars %>%
mutate(cyl = factor(cyl),
am = factor(am, labels = c("Automatic","Manual")))
str(df)
## 'data.frame': 32 obs. of 11 variables: ## $ mpg : num 21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ... ## $ cyl : Factor w/ 3 levels "4","6","8": 2 2 1 2 3 2 3 1 1 2 ... ## $ disp: num 160 160 108 258 360 ... ## $ hp : num 110 110 93 110 175 105 245 62 95 123 ... ## $ drat: num 3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ... ## $ wt : num 2.62 2.88 2.32 3.21 3.44 ... ## $ qsec: num 16.5 17 18.6 19.4 17 ... ## $ vs : num 0 0 1 1 0 1 0 1 1 1 ... ## $ am : Factor w/ 2 levels "Automatic","Manual": 2 2 2 1 1 1 1 1 1 1 ... ## $ gear: num 4 4 4 3 3 3 3 4 4 4 ... ## $ carb: num 4 4 1 1 2 1 4 2 2 4 ...
ggplot #1: mpg vs hp (with regression line)
- Negative slope is expected: higher horsepower often lowers fuel economy.
ggplot(df, aes(hp, mpg, color = cyl)) +
geom_point(size = 2, alpha = 0.8) +
geom_smooth(method = "lm", se = TRUE) +
labs(title = "Fuel Economy vs Horsepower",
x = "Horsepower (hp)",
y = "Miles per Gallon (mpg)",
color = "Cylinders") +
theme_minimal(base_size = 14)
## `geom_smooth()` using formula = 'y ~ x'
Fit the simple linear model
fit <- lm(mpg ~ hp, data = df) summary(fit)
## ## Call: ## lm(formula = mpg ~ hp, data = df) ## ## Residuals: ## Min 1Q Median 3Q Max ## -5.7121 -2.1122 -0.8854 1.5819 8.2360 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 30.09886 1.63392 18.421 < 2e-16 *** ## hp -0.06823 0.01012 -6.742 1.79e-07 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 3.863 on 30 degrees of freedom ## Multiple R-squared: 0.6024, Adjusted R-squared: 0.5892 ## F-statistic: 45.46 on 1 and 30 DF, p-value: 1.788e-07
- The coefficient for hp is \(\hat{\beta}_1\). A negative value means mpg drops as hp rises.
Math #1: Model & OLS estimator
- Model: \(Y_i = \beta_0 + \beta_1 x_i + \varepsilon_i\), with \(\mathbb{E}[\varepsilon_i]=0\), \(\operatorname{Var}(\varepsilon_i)=\sigma^2\), and iid errors.
- Least Squares chooses \((\hat{\beta}_0,\hat{\beta}_1)\) to minimize \[ S(\beta_0,\beta_1) = \sum_{i=1}^n (y_i - \beta_0 - \beta_1 x_i)^2. \]
- The slope estimator can be written as \[ \hat{\beta}_1 = \frac{\sum_i (x_i-\bar x)(y_i-\bar y)}{\sum_i (x_i-\bar x)^2}. \]
ggplot #2: Residuals vs Fitted
aug <- data.frame(
fitted = fitted(fit),
resid = resid(fit)
)
ggplot(aug, aes(fitted, resid)) +
geom_hline(yintercept = 0, linetype = "dashed") +
geom_point(alpha = 0.9) +
labs(title = "Residuals vs Fitted",
x = "Fitted values",
y = "Residuals") +
theme_minimal(base_size = 14)
plotly (3D): mpg ~ hp + wt
- 3D scatter to visualize mpg vs hp and wt.
plot_ly(
data = df,
x = ~hp, y = ~wt, z = ~mpg,
type = "scatter3d", mode = "markers",
marker = list(size = 4)
) %>%
layout(scene = list(
xaxis = list(title = "Horsepower"),
yaxis = list(title = "Weight (1000 lbs)"),
zaxis = list(title = "MPG")
))
Math #2: t-test & p-value
- Hypotheses for slope: \(H_0\!:\;\beta_1 = 0\) vs \(H_a\!:\;\beta_1 \neq 0\).
- Test statistic: \[ t = \frac{\hat{\beta}_1 - 0}{\operatorname{SE}(\hat{\beta}_1)} \sim t_{n-2}\quad \text{under } H_0. \]
- The p-value is \(2\cdot P(|T_{n-2}| \ge |t_{obs}|)\). Small p-values suggest the slope is nonzero.
coef(summary(fit))["hp", , drop = FALSE]
## Estimate Std. Error t value Pr(>|t|) ## hp -0.06822828 0.0101193 -6.742389 1.787835e-07
Code slide: end-to-end
# Packages
library(ggplot2); library(dplyr); library(plotly)
# Data
df <- mtcars %>%
mutate(cyl = factor(cyl),
am = factor(am, labels = c("Automatic","Manual")))
# Model
fit <- lm(mpg ~ hp, data = df)
# Plot 1 (ggplot): scatter + lm
p1 <- ggplot(df, aes(hp, mpg, color = cyl)) +
geom_point(size = 2, alpha = 0.8) +
geom_smooth(method = "lm", se = TRUE) +
labs(title = "Fuel Economy vs Horsepower", x = "Horsepower", y = "MPG", color = "Cylinders") +
theme_minimal(base_size = 14)
# Plot 2 (ggplot): residuals vs fitted
aug <- data.frame(fitted = fitted(fit), resid = resid(fit))
p2 <- ggplot(aug, aes(fitted, resid)) +
geom_hline(yintercept = 0, linetype = "dashed") +
geom_point(alpha = 0.9) +
labs(title = "Residuals vs Fitted", x = "Fitted", y = "Residuals") +
theme_minimal(base_size = 14)
# Plot 3 (plotly 3D)
p3 <- plot_ly(df, x = ~hp, y = ~wt, z = ~mpg, type = "scatter3d", mode = "markers")
# Show
p1; p2; p3
## `geom_smooth()` using formula = 'y ~ x'
