Chapter 7, Moving Beyond Linearity - Question 11

(a)

set.seed(1)
x1 = rnorm(100)
x2 = rnorm(100)
err = rnorm(100)

(b)

b0 = 3
b1 = 5
b2 = 8
y = b0 + b1*x1 + b2*x2 + err

(c)

b1_hat = 6
m = y-b1_hat*x1
b2_hat = lm(m~x2)$coef[2]
b2_hat

x2: 7.94744587206254

(d)

m=y-b2_hat*x2
b1_hat = lm(m~x1)$coef[2]
b1_hat

x1: 5.02111113400588

(e)

b1_hat = 6
b0_hat_lst = c()
b1_hat_lst = c()
b2_hat_lst = c()

for(i in 1:1000){
    m1 = y-b1_hat*x1
    model = lm(m1~x2)
    b2_hat = model$coef[2]
    
    m2 = y-b2_hat*x2
    model = lm(m2~x1)
    b1_hat = model$coef[2]
    
    b0_hat_lst[i]=model$coef[1]
    b1_hat_lst[i]=b1_hat
    b2_hat_lst[i]=b2_hat
}

plot(NULL,NULL,xlab="Iteration",ylab="coef",xlim=c(1,1000),ylim=c(1,10))
lines(b0_hat_lst,col="red")
lines(b1_hat_lst,col="blue")
lines(b2_hat_lst,col="green")
legend(0,10,legend=c("b0_hat","b1_hat","b2_hat"),col=c("red","blue","green"),lty=1, cex=0.8)

(f)

lm.model = lm(y~x1+x2)
summary(lm.model)

Call:
lm(formula = y ~ x1 + x2)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.94359 -0.43645  0.00202  0.63692  2.63941 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   3.0254     0.1052   28.76   <2e-16 ***
x1            5.0211     0.1168   43.01   <2e-16 ***
x2            7.9465     0.1095   72.58   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.043 on 97 degrees of freedom
Multiple R-squared:  0.9865,	Adjusted R-squared:  0.9863 
F-statistic:  3556 on 2 and 97 DF,  p-value: < 2.2e-16

plot(NULL,NULL,xlab="Iteration",ylab="coef",xlim=c(1,1000),ylim=c(1,10))
lines(b0_hat_lst,col="red")
lines(b1_hat_lst,col="blue")
lines(b2_hat_lst,col="green")
abline(lm.model$coef[1],0,col="red",lty=2)
abline(lm.model$coef[2],0,col="blue",lty=2)
abline(lm.model$coef[3],0,col="green",lty=2)
legend(0,10,legend=c("b0_hat","b1_hat","b2_hat"),col=c("red","blue","green"),lty=1, cex=0.8)

(g)

One iteration is enough to get a good approximation of the coefficients