Full worked answers are in the solutions sheet.
Three exercises on Lectures 11 and 12. They illustrate the material from the lectures, and at times go a little beyond it. Exercises 1 and 3 are pen-and-paper: work them out with a calculator and the statistical tables first, as you would at the midterm, and only then check your numbers in R. Exercise 2 is done in R.
The dataset stackloss, which ships with R, records 21
days of operation of a chemical plant. The response,
stack.loss, measures how much raw material the plant lost
that day: the larger it is, the less efficient the plant. The three
predictors are settings of the plant:
Air.Flow: how fast the plant was running (larger means
faster);Water.Temp: the temperature of the water used to cool
it (°C);Acid.Conc.: the concentration of a liquid used in the
process.No knowledge of chemistry is needed: think of them as three settings that may affect the loss.
A regression of stack.loss on the three predictors gives
the following output, with some values removed.
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -39.9197 11.8960 -3.356 0.00375 **
Air.Flow 0.7156 (a) 5.307 5.8e-05 ***
Water.Temp 1.2953 0.3680 3.520 0.00263 **
Acid.Conc. -0.1521 0.1563 (b) 0.34405
---
Residual standard error: 3.243 on (c) degrees of freedom
Multiple R-squared: 0.9136, Adjusted R-squared: (e)
F-statistic: (d) on 3 and (c) DF, p-value: 3.016e-09Questions.
Water.Temp, and interpret it in context.Acid.Conc. is related to stack loss, at
the 5% level. Does your conclusion mean that acid concentration has
nothing to do with the efficiency of the plant?Use the same data and the same model as in Exercise 1.
fit <- lm(stack.loss ~ Air.Flow + Water.Temp + Acid.Conc., data = stackloss)Questions.
Consider a day with Air.Flow = 60,
Water.Temp = 20 and Acid.Conc. = 87. Using
predict, give a 95% confidence interval for the mean stack
loss of all days with these settings, and a 95% prediction interval for
the stack loss of one such day.
With
new <- data.frame(Air.Flow = 60, Water.Temp = 20, Acid.Conc. = 87),
predict(fit, new, se.fit = TRUE)$se.fit returns the
standard error \(s\{\hat Y_h\}\) that
the confidence interval uses. Check that you obtain the same value by
computing it yourself: \(s^2\{\hat Y_h\} =
\boldsymbol{X}_h'\,\boldsymbol{s^2}\{\hat{\boldsymbol{\beta}}\}\,\boldsymbol{X}_h\),
with \(\boldsymbol{X}_h = (1, 60, 20,
87)'\), where vcov(fit) gives \(\boldsymbol{s^2}\{\hat{\boldsymbol{\beta}}\} =
MSE\,(\boldsymbol{X}'\boldsymbol{X})^{-1}\).
Do the same at Air.Flow = 75,
Water.Temp = 18, Acid.Conc. = 87. Is each of
these values within the range of the data? Compare the widths of the
intervals with those of (a).
Plot Water.Temp against Air.Flow for
the 21 days, and add the two points of (a) and (b). Compute \(\boldsymbol{X}_h'(\boldsymbol{X}'\boldsymbol{X})^{-1}\boldsymbol{X}_h\)
at both points, and compare these two numbers with the values of \(\boldsymbol{X}_i'(\boldsymbol{X}'\boldsymbol{X})^{-1}\boldsymbol{X}_i\)
for the 21 observed days (hatvalues(fit) gives them). What
do you conclude about the intervals of (b)? You can use the function
lev below, which returns \(\boldsymbol{X}_h'(\boldsymbol{X}'\boldsymbol{X})^{-1}\boldsymbol{X}_h\)
for a vector \(\boldsymbol{X}_h\), for
example lev(c(1, 60, 20, 87)).
X <- model.matrix(fit) # the design matrix
XtX_inv <- solve(t(X) %*% X)
lev <- function(xh) drop(t(xh) %*% XtX_inv %*% xh)An agronomist records the yield of \(n = 25\) plots of wheat, for which \(SST = 400\). Model A regresses the yield on the rainfall alone, with \(SSE = 160\). Model B adds two soil measurements to the rainfall, with \(SSE = 150\). Both models have an intercept.
Questions.