13-03: Exercises — Regression Inference and Multiple Regression¶
Notes reference: 13-03: Regression Inference and Multiple Regression
Q1: Test the slope¶
A regression on n = 25 observations gives b₁ = 2.84 with SE(b₁) = 0.71.
Test H₀: β₁ = 0 at α = 0.05 and build a 95% CI for β₁.
Solution
STEP 1 H₀: β₁ = 0 H₁: β₁ ≠ 0 TWO-tailed, α = 0.05
STEP 2 df = n − 2 = 23
STEP 3 b₁ − 0 2.84
t = ──────── = ────── = 4.00000
SE(b₁) 0.71
STEP 4 critical t(0.025, 23) = ±2.06866
|4.00000| > 2.06866 → REJECT H₀
p = 2 × P(T₂₃ > 4.0) = 2(0.000276) = 0.000552
STEP 5 There is strong evidence at the 5% level that x has linear
predictive value for y.
95% CI FOR β₁
2.84 ± 2.06866(0.71) = 2.84 ± 1.4688 = (1.371, 4.309)
The interval EXCLUDES 0 — consistent with rejecting H₀ — and it says
a one-unit increase in x is associated with a rise in y of somewhere
between 1.37 and 4.31 units. That range is what a decision-maker
actually needs; the p-value alone would not supply it.
=2.84/0.71 ' t -> 4.0
=T.INV.2T(0.05, 23) ' 2.06866
=T.DIST.2T(4.0, 23) ' 0.000552
=2.84-2.06866*0.71 ' 1.3712
=2.84+2.06866*0.71 ' 4.3088
Q2: Three equivalent tests¶
For a simple regression with r = 0.72 and n = 30, show that testing the correlation, testing the slope, and the model F-test all give the same answer.
Solution
TEST 1 — CORRELATION
______
/n − 2 ______
t = r · √ ─────── = 0.72 · √ 28/(1 − 0.5184)
1 − r²
______
= 0.72 · √ 28/0.4816
= 0.72 × √58.1395
= 0.72 × 7.62493
= 5.48995 df = 28
TEST 2 — SLOPE
The slope t-statistic is ALGEBRAICALLY IDENTICAL:
t = b₁/SE(b₁) = 5.48995 df = 28
TEST 3 — MODEL F-TEST
F = t² = 5.48995² = 30.1396 df = (1, 28)
Equivalently: F = r²(n−2)/(1−r²) = 0.5184(28)/0.4816 = 30.1396 ✓
ALL THREE p-VALUES ARE IDENTICAL
2 × P(T₂₈ > 5.48995) = P(F₁,₂₈ > 30.1396) = 0.0000072
WHY: for a SIMPLE regression there is only one predictor, so "does x
predict y?", "is the slope non-zero?", and "is the model useful?" are
literally the same question. The three tests are three faces of one
computation.
THIS BREAKS DOWN with more than one predictor: the F-test then asks
whether ANY predictor helps, while each t-test asks about ONE predictor
GIVEN the others.
# with data x, y
m <- lm(y ~ x)
summary(m)$coefficients["x", "t value"] # slope t
cor.test(x, y)$statistic # correlation t — identical
summary(m)$fstatistic[1] # F = t²
Q3: Read a full regression output¶
Call: lm(formula = price ~ sqft)
Residuals:
Min 1Q Median 3Q Max
-42.318 -12.075 -0.891 11.284 48.902
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 42.15600 11.82400 3.565 0.00108 **
sqft 0.11842 0.00614 19.287 < 2e-16 ***
Residual standard error: 19.84 on 36 degrees of freedom
Multiple R-squared: 0.9118, Adjusted R-squared: 0.9094
F-statistic: 372.0 on 1 and 36 DF, p-value: < 2.2e-16
- Write the fitted equation.
- Interpret the slope.
- What is
n? - Test the slope at
α = 0.05. - Build a 95% CI for the slope.
- Verify
F = t². - Predict the price of a 2,000 sq ft house.
Solution
1. price-hat = 42.156 + 0.11842 × sqft (price in $000s)
2. Each additional SQUARE FOOT is associated with a $118.42 increase in
predicted price (0.11842 × $1,000). Per 100 sq ft that is $11,842.
3. df_residual = n − 2 = 36 → n = 38 houses
4. t = 19.287, p < 2e-16 ≤ 0.05 → REJECT H₀: β₁ = 0.
Floor area has overwhelming linear predictive value for price.
5. t(0.025, 36) = 2.02809
0.11842 ± 2.02809(0.00614) = 0.11842 ± 0.012452 = (0.10597, 0.13087)
→ between $106 and $131 per square foot, with 95% confidence.
6. t² = 19.287² = 371.99 ≈ F = 372.0 ✓
7. price-hat = 42.156 + 0.11842(2000) = 42.156 + 236.84 = 278.996
→ about $279,000
CAUTION: this is a POINT estimate. With s_e = 19.84, a 95%
PREDICTION interval for one house is roughly
279.0 ± 2.028 × 19.84 × √(1 + 1/38 + …) ≈ 279 ± 41
i.e. roughly $238,000 to $320,000 — a range no estate agent would
call precise. Report the interval, not just the point.
ALSO WORTH NOTING
R² = 0.9118 — floor area alone explains 91% of price variation.
The residual quartiles (−12.1, +11.3) are roughly symmetric about 0,
which is reassuring, but you would still plot the residuals.
Q4: Multiple regression coefficients¶
Estimate Std. Error t value Pr(>|t|)
(Intercept) 18.2400 6.1200 2.980 0.0048
experience 2.4100 0.3800 6.342 1.2e-07
education 1.8700 0.6500 2.877 0.0063
age -0.0900 0.1400 -0.643 0.5238
Residual standard error: 4.72 on 41 degrees of freedom
Multiple R-squared: 0.7823, Adjusted R-squared: 0.7664
F-statistic: 49.11 on 3 and 41 DF, p-value: 3.2e-14
- Write the fitted equation.
- Interpret the
experiencecoefficient. - Which predictors are significant at
α = 0.05? - What is
n? - Should
agestay in the model? - Interpret
R²and adjustedR².
Solution
1. ŷ = 18.24 + 2.41(experience) + 1.87(education) − 0.09(age)
2. HOLDING EDUCATION AND AGE CONSTANT, each additional year of
experience is associated with a 2.41-unit increase in salary.
The phrase "holding the others constant" is not optional — it is
what distinguishes a multiple-regression coefficient from a simple
one.
3. SIGNIFICANT at α = 0.05:
experience p = 1.2e-07 ✓
education p = 0.0063 ✓
intercept p = 0.0048 ✓ (rarely of interest)
NOT significant:
age p = 0.5238 ✗
4. df_residual = n − k − 1 = 41, with k = 3 predictors
n = 41 + 3 + 1 = 45
5. AGE: t = −0.643, p = 0.52 — no evidence it adds anything ONCE
experience and education are in the model.
Arguments for DROPPING it: simpler model, one more residual df,
and age is probably highly correlated with experience (check the VIF).
Arguments for KEEPING it: if theory demands it as a control variable,
or if the goal is prediction rather than explanation.
RECOMMENDED: check VIF(age); if it is high, the non-significance is
a multicollinearity symptom, not evidence that age is irrelevant.
Then refit without age and compare adjusted R² and AIC.
6. R² = 0.7823 78% of salary variation is explained by the model
Adj R² = 0.7664 the honest figure once the 3 predictors are
charged for. The 1.6-point gap is small, which
suggests the predictors are mostly earning their
place. Use ADJUSTED R² when comparing models with
different numbers of predictors.
Overall F(3, 41) = 49.11, p = 3.2e-14 → the model as a whole is
highly useful, even though one individual predictor is not.
Q5: Adjusted R²¶
Three models on the same n = 50 data set:
| Model | k |
R² |
|---|---|---|
| A | 2 | 0.680 |
| B | 4 | 0.702 |
| C | 7 | 0.715 |
Compute adjusted R² and choose.
Solution
n − 1
Adj R² = 1 − (1 − R²) · ───────────
n − k − 1
Model A: 1 − (1 − 0.680) × 49/47 = 1 − 0.320(1.042553) = 1 − 0.333617 = 0.6664
Model B: 1 − (1 − 0.702) × 49/45 = 1 − 0.298(1.088889) = 1 − 0.324489 = 0.6755
Model C: 1 − (1 − 0.715) × 49/42 = 1 − 0.285(1.166667) = 1 − 0.332500 = 0.6675
| Model | k | R² | Adj R² |
|-------|---|-------|--------|
| A | 2 | 0.680 | 0.6664 |
| B | 4 | 0.702 | 0.6755 | ← BEST
| C | 7 | 0.715 | 0.6675 |
CHOOSE MODEL B.
WHY: raw R² always rises as predictors are added — Model C looks best on
R² (0.715) but WORST-but-one on adjusted R². The three extra predictors
in C bought only 1.3 points of R², which is less than they cost in
degrees of freedom.
RULE: compare models with ADJUSTED R², AIC, or BIC — never with raw R².
Also confirm with a nested F-test (anova(mB, mC)) and check that the
extra predictors are theoretically justified, not just data-dredged.
adj_r2 <- function(r2, n, k) 1 - (1 - r2) * (n - 1) / (n - k - 1)
adj_r2(c(0.680, 0.702, 0.715), 50, c(2, 4, 7))
Q6: Multicollinearity¶
A model predicting house price includes sqft, bedrooms, and rooms_total. The overall F is highly significant (p < 0.001) but no individual predictor has p < 0.05. VIFs are 8.4, 12.7, and 15.2.
Diagnose and fix.
Solution
DIAGNOSIS: SEVERE MULTICOLLINEARITY.
THE TELL-TALE PATTERN
• significant overall F BUT no significant individual t
• VIFs above 10 (12.7 and 15.2 here)
• coefficients with implausible signs or huge standard errors
• coefficients that swing wildly when a predictor is added or removed
WHY IT HAPPENS
sqft, bedrooms, and rooms_total measure almost the same thing —
house size. Once sqft is in the model, bedrooms adds almost no NEW
information, so its individual effect cannot be estimated precisely.
The variance of each coefficient is inflated by the VIF factor:
VIF = 15.2 means SE(b) is √15.2 = 3.9× larger than it would be with
uncorrelated predictors.
WHAT IT DOES *NOT* BREAK
• PREDICTION is fine. If you only want ŷ, multicollinearity is not
a problem — R² and the fitted values are unaffected.
• The overall F-test is valid.
FIXES
1. DROP redundant predictors. Keep sqft; drop rooms_total (highest
VIF). Refit and re-check.
2. COMBINE them into an index (e.g. sqft per room), or use principal
components.
3. CENTER the predictors if the collinearity comes from an
interaction or polynomial term.
4. COLLECT MORE DATA — more observations shrink every standard error.
5. Use RIDGE REGRESSION, which is designed for this situation.
VIF BENCHMARKS
< 5 fine
5 – 10 worth investigating
> 10 serious — act on it
library(car)
m <- lm(price ~ sqft + bedrooms + rooms_total, data = houses)
vif(m)
round(cor(houses[, c("sqft","bedrooms","rooms_total")]), 3)
m2 <- lm(price ~ sqft + bedrooms, data = houses) # drop the worst offender
vif(m2)
anova(m2, m) # did dropping it hurt?
AIC(m, m2)
from statsmodels.stats.outliers_influence import variance_inflation_factor
X = sm.add_constant(houses[["sqft", "bedrooms", "rooms_total"]])
pd.Series([variance_inflation_factor(X.values, i) for i in range(X.shape[1])],
index=X.columns)
Q7: Categorical predictors and interaction¶
Predict salary from years (numeric) and dept (Sales / Engineering / Marketing).
- How many dummy variables?
- Write the model.
- Interpret each coefficient.
- What does adding an interaction change?
Solution
1. 3 levels → 2 dummy variables. One level is the REFERENCE.
Never create three — that is the dummy-variable trap.
2. Let Engineering be the reference (alphabetically first in R).
ŷ = b₀ + b₁(years) + b₂(Marketing) + b₃(Sales)
where Marketing = 1 if the employee is in Marketing, else 0
Sales = 1 if the employee is in Sales, else 0
Engineering is coded 0, 0
3. Suppose the fit is
ŷ = 45.2 + 2.8(years) − 6.4(Marketing) − 3.1(Sales)
b₀ = 45.2 Predicted salary for an ENGINEERING employee with
0 years of experience.
b₁ = 2.8 Each additional year adds 2.8 units — the SAME slope
for all three departments (parallel lines).
b₂ = −6.4 Marketing employees earn 6.4 units LESS than Engineering
employees WITH THE SAME YEARS of experience.
b₃ = −3.1 Sales employees earn 3.1 units less than Engineering,
holding years constant.
To compare Marketing with Sales: −6.4 − (−3.1) = −3.3, so Marketing
is 3.3 below Sales. (To test that difference formally, refit with
Sales as the reference.)
4. ADDING AN INTERACTION (salary ~ years * dept)
ŷ = b₀ + b₁(years) + b₂(Mkt) + b₃(Sales)
+ b₄(years × Mkt) + b₅(years × Sales)
This lets each department have its OWN SLOPE, not just its own
intercept — the lines are no longer parallel.
b₄ ≠ 0 would mean experience is rewarded at a DIFFERENT RATE in
Marketing than in Engineering.
RULE: always keep the main effects when you include their
interaction, and interpret the main effect of `years` as the slope
for the REFERENCE group only.
m1 <- lm(salary ~ years + dept, data = df) # parallel lines
m2 <- lm(salary ~ years * dept, data = df) # separate slopes
anova(m1, m2) # is the interaction worth it?
levels(df$dept) # the first level is the reference
df$dept <- relevel(df$dept, ref = "Sales") # change the reference
model.matrix(m1) # see the dummy columns
library(ggplot2)
ggplot(df, aes(years, salary, colour = dept)) +
geom_point() + geom_smooth(method = "lm", se = FALSE) + theme_minimal()
m1 = ols("salary ~ years + C(dept)", data=df).fit()
m2 = ols("salary ~ years * C(dept)", data=df).fit()
sm.stats.anova_lm(m1, m2)
# Change the reference: C(dept, Treatment(reference='Sales'))
' Build the dummies by hand — Engineering is the reference:
=IF($C2="Marketing",1,0) ' D2
=IF($C2="Sales",1,0) ' E2
' Interaction columns:
=$B2*D2 ' F2: years × Marketing
=$B2*E2 ' G2: years × Sales
' Then: Data ▸ Data Analysis ▸ Regression
' Input Y Range = salary, Input X Range = B:G (must be CONTIGUOUS)
Q8: Build and compare models¶
You have y plus predictors x1, x2, x3. Fit a sequence of models and choose one.
Solution
m0 <- lm(y ~ 1, data = df) # intercept only — the baseline
m1 <- lm(y ~ x1, data = df)
m2 <- lm(y ~ x1 + x2, data = df)
m3 <- lm(y ~ x1 + x2 + x3, data = df)
# ── Compare ────────────────────────────────────────────────────────
data.frame(
model = c("m1", "m2", "m3"),
r2 = sapply(list(m1, m2, m3), function(m) summary(m)$r.squared),
adj_r2 = sapply(list(m1, m2, m3), function(m) summary(m)$adj.r.squared),
aic = sapply(list(m1, m2, m3), AIC),
bic = sapply(list(m1, m2, m3), BIC)
)
anova(m1, m2, m3) # nested F-tests: does each addition help?
# ── Diagnose the chosen model ──────────────────────────────────────
best <- m2
par(mfrow = c(2, 2)); plot(best); par(mfrow = c(1, 1))
library(car); vif(best); ncvTest(best); durbinWatsonTest(best)
shapiro.test(residuals(best))
which(cooks.distance(best) > 4/nrow(df))
# ── Report ─────────────────────────────────────────────────────────
summary(best)
confint(best)
import statsmodels.api as sm
from statsmodels.formula.api import ols
m1 = ols("y ~ x1", df).fit()
m2 = ols("y ~ x1 + x2", df).fit()
m3 = ols("y ~ x1 + x2 + x3", df).fit()
pd.DataFrame({
"model": ["m1", "m2", "m3"],
"r2": [m.rsquared for m in (m1, m2, m3)],
"adj_r2": [m.rsquared_adj for m in (m1, m2, m3)],
"aic": [m.aic for m in (m1, m2, m3)],
"bic": [m.bic for m in (m1, m2, m3)],
})
sm.stats.anova_lm(m1, m2, m3)
print(m2.summary())
m2.conf_int()
' Data ▸ Data Analysis ▸ Regression, once per model.
' Record from each output block:
' Adjusted R Square, Standard Error, Significance F,
' and the P-value of every coefficient.
' Excel has no AIC/BIC — compare on Adjusted R Square and on whether
' each added predictor is individually significant.
A DEFENSIBLE MODEL-SELECTION PROCEDURE
1. Start from THEORY, not from the data. Decide which predictors
belong in the model before you look at any p-values.
2. Fit and compare with adjusted R², AIC/BIC, and nested F-tests.
3. Check multicollinearity (VIF) BEFORE dropping a non-significant
predictor — the two look identical in the output.
4. Run the residual diagnostics on the model you choose, not on all
of them.
5. Report the model you SELECTED and say how you selected it.
WHAT NOT TO DO
• Automated stepwise selection on many predictors. It inflates R²,
biases the p-values, and rarely replicates.
• Adding predictors until R² looks good. That is overfitting.
• Keeping n small relative to k. Aim for at least 10–20 observations
per predictor.
⬅️ Previous: 13-02: Exercises — Simple Linear Regression
🎉 That is every exercise set. Next: apply it all in the Projects, or test yourself in the Quiz Hub.