##      cte X1 X2
## [1,]   1  1  0
## [2,]   1  2  0
## [3,]   1  3  0
## [4,]   1  4  1
## [5,]   1  5  1
## [6,]   1  6  1

##        y
## [1,]  90
## [2,] 100
## [3,] 110
## [4,] 135
## [5,] 145
## [6,] 165

Calcular producto escalar 1° y 2° columna de x (se transpuso la columna 1).

##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    1    1    1    1    1

##      [,1]
## [1,]    1
## [2,]    2
## [3,]    3
## [4,]    4
## [5,]    5
## [6,]    6

##      [,1]
## [1,]   21

Calculo de X'X

tx<- t(mat_x)
print(tx)

##     [,1] [,2] [,3] [,4] [,5] [,6]
## cte    1    1    1    1    1    1
## X1     1    2    3    4    5    6
## X2     0    0    0    1    1    1

txx<- tx %*% mat_x
print(txx)

##     cte X1 X2
## cte   6 21  3
## X1   21 91 15
## X2    3 15  3

Calculo de X'Y (como ya hemos determinado la transpuesta de X está de más volver a poner el proceso así que directamente multiplicamos la transpuesta de X por Y).

txy<- tx %*% mat_y
print(txy)

##        y
## cte  745
## X1  2875
## X2   445

1° Calcular la inversa de X'X

inv_txx <- solve(txx)
print(inv_txx,2)

##      cte    X1    X2
## cte  1.3 -0.50  1.17
## X1  -0.5  0.25 -0.75
## X2   1.2 -0.75  2.92

print(txy)

##        y
## cte  745
## X1  2875
## X2   445

Calculo del producto \([X'X]^{-1}\times (X'Y)\)

Beta<- inv_txx %*% txy
rownames(Beta)<- c("Bo", "B1", "B2")
round(Beta,2)

##        y
## Bo 75.00
## B1 12.50
## B2 10.83

Calculo de matriz de proyección P.

transpuesta_x<- t(mat_x)
matrizA<- inv_txx %*% transpuesta_x
matrizP<- mat_x %*% matrizA
round(matrizP,3)

##        [,1]  [,2]   [,3]   [,4]  [,5]   [,6]
## [1,]  0.583 0.333  0.083  0.250 0.000 -0.250
## [2,]  0.333 0.333  0.333  0.000 0.000  0.000
## [3,]  0.083 0.333  0.583 -0.250 0.000  0.250
## [4,]  0.250 0.000 -0.250  0.583 0.333  0.083
## [5,]  0.000 0.000  0.000  0.333 0.333  0.333
## [6,] -0.250 0.000  0.250  0.083 0.333  0.583

Comprabar si el producto \(P\times P=P\)

m_idemp <- matrizP %*% matrizP
round(m_idemp,3)

##        [,1]  [,2]   [,3]   [,4]  [,5]   [,6]
## [1,]  0.583 0.333  0.083  0.250 0.000 -0.250
## [2,]  0.333 0.333  0.333  0.000 0.000  0.000
## [3,]  0.083 0.333  0.583 -0.250 0.000  0.250
## [4,]  0.250 0.000 -0.250  0.583 0.333  0.083
## [5,]  0.000 0.000  0.000  0.333 0.333  0.333
## [6,] -0.250 0.000  0.250  0.083 0.333  0.583

Comprobamos que efectivamente la matriz P es una matriz idempotente.

Calculo de proyección Y sobre X, donde \(Y^ =PY\)

Y_<- matrizP %*% mat_y 
round(Y_,2)

##           y
## [1,]  87.50
## [2,] 100.00
## [3,] 112.50
## [4,] 135.83
## [5,] 148.33
## [6,] 160.83

Calculo del error \[ \varepsilon = Y - ¯Y \]

error<- mat_y - Y_
round(error,2)

##          y
## [1,]  2.50
## [2,]  0.00
## [3,] -2.50
## [4,] -0.83
## [5,] -3.33
## [6,]  4.17

Obtener autovalores de la matriz \[[X'X]\]

autovalorx<-eigen(txx)
round(autovalorx$values,3)

## [1] 98.357  1.378  0.266

Observamos que todos los autovalores son positivos, porque implica que la condicón es un mínimo.

Obtener autovalores de la matriz P. Comprobar que la traza P es igual a la suma de autovalores.

autovalorp<- eigen(matrizP)
round(autovalorp$values,2)

## [1] 1 1 1 0 0 0

print(autovalorp$values,5)

## [1]  1.0000e+00  1.0000e+00  1.0000e+00  1.9752e-15  1.7902e-16 -3.1122e-15

Autovalores de P aproximados y autovalores con 5 decimales.

¿La traza de P (\(\sum\) de la diagonal principal) es igual a la suma de autovalores?

trazaP<- sum(diag(matrizP))
print(trazaP)

## [1] 3

sum_autovp<- sum(autovalorp$values)
print(sum_autovp)

## [1] 3

Sí, la traza de p y la suma de sus autovalores son iguales.

Obtener autovalor de la matriz \(I-P\)

identidad<- diag(1,6,6)
id_p<- identidad - matrizP
round(id_p,3)

##        [,1]   [,2]   [,3]   [,4]   [,5]   [,6]
## [1,]  0.417 -0.333 -0.083 -0.250  0.000  0.250
## [2,] -0.333  0.667 -0.333  0.000  0.000  0.000
## [3,] -0.083 -0.333  0.417  0.250  0.000 -0.250
## [4,] -0.250  0.000  0.250  0.417 -0.333 -0.083
## [5,]  0.000  0.000  0.000 -0.333  0.667 -0.333
## [6,]  0.250  0.000 -0.250 -0.083 -0.333  0.417

auto_idp <- eigen(id_p)
round(auto_idp$values,2)

## [1] 1 1 1 0 0 0

Obtener autovalor de \(I+P\)

idmasp<- identidad + matrizP
round(idmasp,3)

##        [,1]  [,2]   [,3]   [,4]  [,5]   [,6]
## [1,]  1.583 0.333  0.083  0.250 0.000 -0.250
## [2,]  0.333 1.333  0.333  0.000 0.000  0.000
## [3,]  0.083 0.333  1.583 -0.250 0.000  0.250
## [4,]  0.250 0.000 -0.250  1.583 0.333  0.083
## [5,]  0.000 0.000  0.000  0.333 1.333  0.333
## [6,] -0.250 0.000  0.250  0.083 0.333  1.583

auto_idmasp<- eigen(idmasp)
round(auto_idmasp$values,2)

## [1] 2 2 2 1 1 1