2
$\begingroup$

I am having difficulty reproducing linear regression coefficients from Casella and Berger. On page 583 figure 12.2.2. He shows two regression lines I am interested in.

(a) a regression of y on x

(b) a regression of x on y

I can get the beta hat matrix for (a), but not for (b).

Here is my work to get (a)

The key is that I am using the following for my beta hat matrix

$$ {\hat{\beta}} = \left(X^\mathsf{T}X\right)^{-1} X^\mathsf{T}Y $$

df <- structure(list(Y = c(3.22, 4.87, 0.12, 2.31, 4.25, 2.24, 2.81, 
3.71, 3.11, 0.9, 4.39, 4.36, 1.26, 3.13, 4.05, 2.28, 3.6, 5.39, 
4.12, 3.16, 4.4, 1.18, 2.54, 4.89), B0 = c(1, 1, 1, 1, 1, 1, 
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1), B1 = c(3.74, 
3.66, 0.78, 2.4, 2.18, 1.96, 0.2, 2.5, 3.5, 1.35, 2.36, 3.13, 
1.22, 1, 1.29, 0.95, 1.05, 2.92, 1.76, 0.51, 2.17, 1.99, 1.53, 
2.6)), class = "data.frame", row.names = c(NA, -24L))

Y <- df$Y
B0 <- df$B0 # the intercept
B1 <- df$B1
X <- cbind(B0, B1)
solve(t(X)%*%X)%*%t(X)%*%Y

For (a) I correctly get $y = 1.86 + .68x$


I was assuming that a regression of x on y would just swap the $X$ and $Y$

This would give me the following beta hat matrix:

$$ {\hat{\beta}} = \left(Y^\mathsf{T}Y\right)^{-1} Y^\mathsf{T}X $$

The problem is that this does not match with the text

solve(t(Y)%*%Y)%*%t(Y)%*%X

My coefficients are

$y = .27 + .57x$

But the correct values are

$y = -2.31 + 2.82x$

enter image description here

My question is this:

What is the correct beta hat matrix when solving for horizontal distance from the point to the fitted line (x regressed on y)?


Solution

The answer is

$$ {\hat{\beta}} = \left(Y^\mathsf{T}Y\right)^{-1} Y^\mathsf{T}X $$

I had added an extra column of ones in my X matrix by mistake when I tried to solve (b).

$\endgroup$

1 Answer 1

1
$\begingroup$

Your dataframe has B1 where it should have X and an unnecessary B0

Try this instead

df <- structure(list(Y = c(3.22, 4.87, 0.12, 2.31, 4.25, 2.24, 2.81, 
3.71, 3.11, 0.9, 4.39, 4.36, 1.26, 3.13, 4.05, 2.28, 3.6, 5.39, 
4.12, 3.16, 4.4, 1.18, 2.54, 4.89),  X = c(3.74, 
3.66, 0.78, 2.4, 2.18, 1.96, 0.2, 2.5, 3.5, 1.35, 2.36, 3.13, 
1.22, 1, 1.29, 0.95, 1.05, 2.92, 1.76, 0.51, 2.17, 1.99, 1.53, 
2.6)), class = "data.frame", row.names = c(NA, -24L))

X1 <- cbind(1, df$X)
solve(t(X1) %*% X1) %*% t(X1) %*% df$Y
lm(df$Y ~ df$X)

which gives your regression line $y = 1.8614 + 0.6763x$ as you found

If you want the regression of $x$ on $y$ then the corresponding code would be

Y1 <- cbind(1, df$Y)
solve(t(Y1) %*% Y1) %*% t(Y1) %*% df$X
lm(df$X ~ df$Y) 

which would give a regression line of $x = 0.8203 + 0.3547y$ and rearranging this would give $y = -2.3124+2.8190 x$

$\endgroup$
1
  • $\begingroup$ Thanks. I did not realize that my X had an extra column of ones when I tried to solve (b). $\endgroup$
    – Alex
    Commented Oct 24, 2022 at 2:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.