Scribd
Cargar
513 vistas3 páginas

Factorización de Cholesky en R

El documento describe un algoritmo para calcular la factorización de Cholesky de una matriz definida positiva A de tamaño n x n. La factorización de Cholesky descompone A como el producto de una matriz triangular inferior L y su transpuesta. El algoritmo verifica primero que A sea simétrica y definida positiva antes de calcular iterativamente los elementos de L de acuerdo con los pasos descritos. Se proveen tres ejemplos ilustrativos.

Cargado por

fermad
Derechos de autor
© All Rights Reserved
Nos tomamos en serio los derechos de los contenidos. Si sospechas que se trata de tu contenido, reclámalo aquí.
Formatos disponibles
Descarga como PDF, TXT o lee en línea desde Scribd
513 vistas3 páginas

Factorización de Cholesky en R

El documento describe un algoritmo para calcular la factorización de Cholesky de una matriz definida positiva A de tamaño n x n. La factorización de Cholesky descompone A como el producto de una matriz triangular inferior L y su transpuesta. El algoritmo verifica primero que A sea simétrica y definida positiva antes de calcular iterativamente los elementos de L de acuerdo con los pasos descritos. Se proveen tres ejemplos ilustrativos.

Cargado por

fermad
Derechos de autor
© All Rights Reserved
Nos tomamos en serio los derechos de los contenidos. Si sospechas que se trata de tu contenido, reclámalo aquí.
Formatos disponibles
Descarga como PDF, TXT o lee en línea desde Scribd

FACTORIZACIÓN DE CHOLESKY

Algoritmo en R para calcular la matriz 𝑳 (triangular inferior) de una matriz definida


positiva 𝑨 de tamaño 𝑛 × 𝑛 tal que 𝑨 = 𝑳𝑳𝒕 .

ALGORITMO
> CHOLESKY = function(n , a){
+ A = matrix(a , n)
+
+ if(n!=sqrt(length(a))){message("ERROR\nla matriz no es Cuadrada")}else{
+ if(sum(sum(abs(t(A)-A)))>0){message("ERROR\nla matriz no es Simétrica")}else{
+ if(det(A)<=0){message("ERROR\nla matriz no es definida positiva")}else{
+ L = matrix(0,n,n)
+ L[1,1] = (A[1,1])^(1/2) #Esto es lo del paso 1
+ for (j in 2:n) {L[j,1] = A[j,1]/L[1,1]} # Esto es lo del paso 2
+ i = 1
+ while (i < (n-1)) {
+ i = i+1
+ k = 1:(i-1)
+ m=i+1
+ L[i,i] = (A[i,i] - sum(L[i,k]^2) )^(1/2) # Esto es lo del paso 4
+ for (j in m:n) {
+ L[j,i] = (1/L[i,i])*(A[j,i]-sum(L[i,k]*L[j,k])) # Esto es lo del paso 5
+ }
+ }
+ L[n,n] = (A[n,n] - sum(L[n,1:(n-1)]^2) )^(1/2) # Esto es lo del paso 6
+ return(L)
+ }
+ }
+ }
+ }

EJEMPLO 1
Sea
4 −2 0 2
−2 2 3 −2
𝑨= ( )
0 3 18 0
2 −2 0 4

Cálculo de la matriz 𝑳
> a = c(4,-2,0,2,-2,2,3,-2,0,3,18,0,2,-2,0,4) ; n=4
> L <- CHOLESKY(n , a)
> L
[,1] [,2] [,3] [,4]
[1,] 2 0 0 0
[2,] -1 1 0 0
[3,] 0 3 3 0
[4,] 1 -1 1 1
Prueba de que 𝑨 = 𝑳𝑳𝒕 .
> L%*%t(L)
[,1] [,2] [,3] [,4]
[1,] 4 -2 0 2
[2,] -2 2 3 -2
[3,] 0 3 18 0
[4,] 2 -2 0 4

EJEMPLO 2
Sea
4 3 7 1 7
3 5 5 9 7
𝑨= 7 5 1 2 4
1 9 2 3 6
(7 7 4 6 7)
Cálculo de la matriz 𝑳
> a = c(4,3,7,1,7,3,5,5,9,7,7,5,1,2,4,1,9,2,3,6,7,7,4,6,7) ; n=5
> L = CHOLESKY(n , a)
ERROR
la matriz no es definida positiva

EJEMPLO 3
Sea
1 1 1 1 1
1 11 2 2 2
𝑨= 1 2 6 3 3
1 2 3 16 1
(1 2 3 1 21)

Cálculo de la matriz 𝑳
> a = c(1,1,1,1,1,1,11,2,2,2,1,2,6,3,3,1,2,3,16,1,1,2,3,1,21) ; n=5
> L = CHOLESKY(n , a) ; L
[,1] [,2] [,3] [,4] [,5]
[1,] 1 0.0000000 0.0000000 0.0000000 0.000000
[2,] 1 3.1622777 0.0000000 0.0000000 0.000000
[3,] 1 0.3162278 2.2135944 0.0000000 0.000000
[4,] 1 0.3162278 0.8583325 3.7634114 0.000000
[5,] 1 0.3162278 0.8583325 -0.2223341 4.371937

Prueba de que 𝑨 = 𝑳𝑳𝒕 .


> L%*%t(L)
[,1] [,2] [,3] [,4] [,5]
[1,] 1 1 1 1 1
[2,] 1 11 2 2 2
[3,] 1 2 6 3 3
[4,] 1 2 3 16 1
[5,] 1 2 3 1 21

También podría gustarte