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
library(MASS)
ginv(tx%*%x)
## [,1] [,2]
## [1,] 1.3414634 -0.4146341
## [2,] -0.4146341 0.1463415
b=ginv(tx%*%x)%*%tx%*%y
b
## [,1]
## [1,] 0.4390244
## [2,] 4.6097561
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
n=nrow(CF)
J=matrix(1,6,6)
SSR=t(b)%*%tx%*%y-(1/n)*ty%*%J%*%y
SSR
## [,1]
## [1,] 145.2073
SSE=t(e)%*%e
SSE
## [,1]
## [1,] 20.29268
SSE2=ty%*%y-t(b)%*%tx%*%y
SSE2#Equivalent Forms
## [,1]
## [1,] 20.29268
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
Xh=matrix(c(1,4),2,1)
t(Xh)%*%b
## [,1]
## [1,] 18.87805
#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
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)
t(a)%*%Expy
## [,1]
## [1,] 140
A%*%Expy
## [,1]
## [1,] 14
## [2,] 32
## [3,] 50
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