Question 1 5.5

CF=read.table("CH05PR05.txt")
x2=matrix(CF[,2],6,1)
x1=matrix(c(1,1,1,1,1,1),6,1)
x=cbind(x1,x2)
y=matrix(CF[,1],6,1)
n=nrow(CF)
tx=t(x)
ty=t(y)
#(1)
ty%*%y
##      [,1]
## [1,] 1259
#(2)
tx%*%x
##      [,1] [,2]
## [1,]    6   17
## [2,]   17   55
#(3)
tx%*%y
##      [,1]
## [1,]   81
## [2,]  261

5.13

library(MASS)
ginv(tx%*%x)
##            [,1]       [,2]
## [1,]  1.3414634 -0.4146341
## [2,] -0.4146341  0.1463415

5.24

(1) estimated regression coefficients

b=ginv(tx%*%x)%*%tx%*%y
b
##           [,1]
## [1,] 0.4390244
## [2,] 4.6097561

(2) residuals

Yhat=x%*%b
e=y-Yhat
e
##             [,1]
## [1,] -2.87804878
## [2,] -0.04878049
## [3,]  0.34146341
## [4,]  0.73170732
## [5,] -1.26829268
## [6,]  3.12195122

(3) SSR

n=nrow(CF)
J=matrix(1,6,6)
SSR=t(b)%*%tx%*%y-(1/n)*ty%*%J%*%y
SSR
##          [,1]
## [1,] 145.2073

(4) SSE

SSE=t(e)%*%e
SSE
##          [,1]
## [1,] 20.29268
SSE2=ty%*%y-t(b)%*%tx%*%y
SSE2#Equivalent Forms
##          [,1]
## [1,] 20.29268

(5) Variance-covariance matrix of b

MSE=(1/(n-2))*SSE
m=matrix(MSE,2,2)
V = m*ginv(tx%*%x)
V
##           [,1]       [,2]
## [1,]  6.805473 -2.1035098
## [2,] -2.103510  0.7424152

(6) point estimate of E{Yh} when Xh= 4

Xh=matrix(c(1,4),2,1)
t(Xh)%*%b
##          [,1]
## [1,] 18.87805

(7) Est variance of Yhath when Xh = 4

#a.
MSE*(1+t(Xh)%*%ginv(tx%*%x)%*%Xh)
##          [,1]
## [1,] 6.929209
#b. 
#s(bo,b1)
V[1,2]
## [1] -2.10351
#s(b0)^2
V[1,1]
## [1] 6.805473
#sb1
sqrt(V[2,2])
## [1] 0.8616352
#c 
H= x%*%ginv(t(x)%*%x)%*%t(x)
H
##             [,1]       [,2]       [,3]      [,4]      [,5]        [,6]
## [1,]  0.36585366 -0.1463415 0.02439024 0.1951220 0.1951220  0.36585366
## [2,] -0.14634146  0.6585366 0.39024390 0.1219512 0.1219512 -0.14634146
## [3,]  0.02439024  0.3902439 0.26829268 0.1463415 0.1463415  0.02439024
## [4,]  0.19512195  0.1219512 0.14634146 0.1707317 0.1707317  0.19512195
## [5,]  0.19512195  0.1219512 0.14634146 0.1707317 0.1707317  0.19512195
## [6,]  0.36585366 -0.1463415 0.02439024 0.1951220 0.1951220  0.36585366
#d
I=diag(n)
MSE=as.vector(MSE)
MSE*(I-H)
##            [,1]       [,2]       [,3]       [,4]       [,5]       [,6]
## [1,]  3.2171327  0.7424152 -0.1237359 -0.9898870 -0.9898870 -1.8560381
## [2,]  0.7424152  1.7323022 -1.9797739 -0.6186794 -0.6186794  0.7424152
## [3,] -0.1237359 -1.9797739  3.7120761 -0.7424152 -0.7424152 -0.1237359
## [4,] -0.9898870 -0.6186794 -0.7424152  4.2070196 -0.8661511 -0.9898870
## [5,] -0.9898870 -0.6186794 -0.7424152 -0.8661511  4.2070196 -0.9898870
## [6,] -1.8560381  0.7424152 -0.1237359 -0.9898870 -0.9898870  3.2171327

Question 2

(a)

Expy = matrix(c(1,2,3),3,1)
Expy
##      [,1]
## [1,]    1
## [2,]    2
## [3,]    3
var=matrix(c(6,0,-4,0,7,5,-4,5,8),3,3)
var
##      [,1] [,2] [,3]
## [1,]    6    0   -4
## [2,]    0    7    5
## [3,]   -4    5    8
a=matrix(c(10,20,30),3,1)
tA=matrix(c(1,2,3,4,5,6,7,8,9),3,3)
A=t(tA)

(b)

t(a)%*%Expy
##      [,1]
## [1,]  140
A%*%Expy
##      [,1]
## [1,]   14
## [2,]   32
## [3,]   50

(c)

t(a)%*%var%*%a
##       [,1]
## [1,] 14200
A%*%var%*%tA
##      [,1] [,2] [,3]
## [1,]  142  301  460
## [2,]  301  667 1033
## [3,]  460 1033 1606