In this post, we will learn about the Cholesky decomposition, definitions of the terms related to this decomposition, the syntax for the function used to implement Cholesky decomposition, and a few examples.
What is Matrix Decomposition?
In linear algebra, matrix decomposition or matrix factorization refers to breaking the matrix down into two or more constituent parts, which, when multiplied together, yield the original matrix.
When the matrix is decomposed into smaller parts, it allows for efficiency and data analysis. Matrix decomposition also optimizes some algorithms, such as least square estimations.
Cholesky Decomposition is the type of decomposition of a matrix into its lower triangular matrix and conjugate transpose.
Mathematically, it is defined as
A=LL* where A is the original matrix, L is a lower triangular matrix, and L* is its conjugate transpose.
There are three conditions for a matrix to be decomposed by Cholesky.
- The matrix should be symmetric
- In the case of complex matrices, it should be Hermitian
- The elements of the matrix should be positive-definite
Let us define each of these terms in detail.
A matrix is considered a square matrix if the number of rows in the matrix is equal to the number of columns.
One example of such a matrix is given below.
Lower Triangular Matrix
A square matrix in which the elements above the principal diagonal are zero is called a lower triangular matrix.
A square matrix is symmetric if it is equal to its transpose. Interchanging rows and columns obtain a transpose of a matrix.
An example is given below.
Conjugate Transpose and Hermitian Matrix
Suppose there is a matrix with complex entries (a+bi) form, then the conjugate of that matrix is obtained by replacing the complex conjugate of all the elements.
The conjugate transpose of a matrix is obtained by interchanging the rows and columns of the conjugate matrix.
Now a Hermitian matrix is a square matrix equal to its conjugate transpose matrix.
A symmetric matrix is said to be positive-definite if all the eigen values of the matrix are positive.
The eigenvalues of a matrix can be computed using the
np.linalg.eigenvalsh function of the NumPy library.
Let us see an example to check if the given matrix is positive definite.
#to check if the matrix is positive definite import numpy as np A = np.array([[2, 1], [1, 2]]) def is_pos_definite(A): eig_vals = np.linalg.eigvalsh(A) return np.all(eig_vals > 0) print(is_pos_definite(A))
import numpy as np: We are importing the NumPy library as np to create arrays and use its functions.
In the next line, we create a 2×2 matrix and store it in a variable called A.
Next, we create a function called
is_pos_definite to check if the matrix is positive definite. The matrix A is passed as an array to this function. The
def keyword is used to define a new function.
Inside this function, we have a new variable eig_vals, which is used to store the eigenvalues of the matrix A. The
np.linalg.eigvalsh is used to compute the eigenvalues of matrix A.
In the next line, we are checking if the eigenvalues of matrix A are all positive(greater than zero) using
Lastly, we are printing the result of this function. If the eigenvalues of matrix A are positive, this line prints True. Else it prints False.
The output is given below.
Applications of Cholesky Decomposition
The main usage of Cholesky’s Decomposition is in Monte Carlo Simulation.
If you are unfamiliar with the Monte Carlo simulation, Check out this article on Monte Carlo in Python.
Apart from higher-level computations, Cholesky Decomposition is also used to solve systems of linear equations, find out matrix inversions, and also perform principal component analysis(PCA).
Check out this article on PCA if you want to know more about data reduction.
Exploring the syntax of
The NumPy implementation of Cholesky decomposition only takes a Symmetric matrix(real-valued) or Hermitian matrix(complex-valued), but in both cases, the matrix should be positive definite.
np.linalg.cholesky takes a square matrix A and returns the LL* of the matrix.
There is no facility to check if the matrix is Hermitian.
The lower triangular matrix of A is used to compute the Cholesky decomposition.
The syntax of this method is as follows.
There is only one parameter in the syntax.
a: The matrix passed as an argument should be array-like. Hermitian and positive-definite matrices are accepted.
Return Type: Returns a lower triangular matrix.
This function may raise an exception when the matrix is not a positive-definite matrix. This exception is called
Now that we have seen the method’s basic functionality let us dive into the examples!
Computing Cholesky Decomposition for Symmetric Matrix
A symmetric matrix is a square matrix equal to its transpose.
Being symmetric is not the only criterion for a matrix to be decomposed. The matrix should also be positive-definite.
Here is the code for the decomposition of a symmetric matrix which is also positive-definite.
# Example1 2x2 symmetric matrix import numpy as np B=np.array([[7,2],[2,1]]) print(B)
In this code snippet, we just import the NumPy library as np, create a two-dimensional array(2×2 matrix), and store it in the variable B.
print(B) is used to print the matrix residing in B.
Let us find the Cholesky decomposition of this matrix.
We are calling the
linalg.cholesky function to decompose matrix B.
This result is stored in L.
The output is a lower triangular matrix, as shown below.
Computing Cholesky Decomposition for Hermitian matrix
A matrix is said to Hermitian if it equals its own conjugate transpose.
Again, being a Hermitian matrix is not enough for decomposition. The matrix should also be a positive-definite matrix.
Here is the code for the Cholesky decomposition of a Hermitian matrix
#2x2 hermitian matrix import numpy as np C=np.array([[2,1-2j],[1+2j,3]]) print(C)
In the first line, we are importing the numpy library.
Next, we create a new variable called C for storing the 2×2 Hermitian matrix with complex numbers.
Next, we print the elements of the matrix to the screen.
The output is given below.
Here is the code for Cholesky’s decomposition.
The decomposed matrix is stored in the variable called L, which is obtained by the method
In the next line, we are printing the elements of the matrix L, which is a lower triangular matrix.
Computing Cholesky Decomposition for 3×3 Hermitian matrix
Let us do something different in this example. Let us create a hermitian matrix and define a new function to check if the matrix is positive definite with the help of
linalg.eigvalsh. Then pass the matrix as an argument to the Cholesky method.
The code to create a 3×3 Hermitian matrix is given below.
#3x3 hermitian matrix D=np.array([[10, -2 + 3j, 1], [-2 - 3j, 6, 2 + 2j], [1, 2 - 2j, 7]]) print(D) def is_pos_definite(D): eig_vals = np.linalg.eigvalsh(D) return np.all(eig_vals > 0) print(is_pos_definite(D))
All we have done in the above snippet is create a three-dimensional array and store it in a variable called D.
Then, we defined a function called is_pos_definite to check if the matrix is positive-definite.
Inside this function, we have also called an instance of
linalg.einvalsh to which the newly created matrix is passed as an argument.
The eig_vals is an object used to return the eigenvalues of matrix D.
np.all is a method of the Numpy library that is used to check if all the eigenvalues of the matrix are greater than zero.
Lastly, we are printing the result of the above computation. The last print statement returns True if the eigenvalues are positive else, it returns False.
What do you think the output will be?
Well, the output is true because the matrix passed as an argument is a positive-definite matrix.
Next, we compute the Cholesky Decomposition.
And the output is as follows.
The output is a bit unclear, so here is the formulated lower triangular matrix.
Illustration of LinAlgError Exception
In the previous examples, we have seen the decomposition of positive-definite matrices. But now let us see what happens if we try to use this method on a matrix that is not a positive-definite matrix.
#example of LinAlgError exception #3x3 hermitian matrix E=np.array([[1,2-3j,3+4j],[2+3j,0,4-5j],[3-4j,4+5j,2]]) print(E) def is_pos_definite(E): eig_vals = np.linalg.eigvalsh(E) return np.all(eig_vals > 0) print(is_pos_definite(E))
In this example, a non-positive-definite Hermitian matrix is stored in a variable called E.
Next, we are using the previously created function
is_pos_definite to check if E has positive eigenvalues.
The output is given below.
The output clearly shows that E is not a positive-definite matrix. What will happen if we pass E as an argument to the
As seen in the above image, the linalgerror is raised. This means that the matrix E is not positive-definite.
To summarize everything we have discussed in this post, we first started with what is a matrix decomposition, then Cholesky decomposition, some basic terms related to linear algebra of matrices such as symmetric matrix, transpose of a matrix, Hermitian matrix, Positive-definite matrix, and so on.
We have also seen the conditions a matrix must satisfy to be decomposed by the Cholesky approach, one of the main conditions being the matrix should be positive-definite.
We have also seen the Numpy approach to check if the matrix is a positive-definite matrix.
Next, we have seen the applications of Cholesky decomposition.
We have also seen the syntax of the method, its parameters, and the return object, which is a lower triangular matrix.
Coming to the examples, we have seen the decomposition of a symmetric matrix and also a Hermitian matrix.
Lastly, we have seen how the LinAlgError exception is raised if the matrix we are trying to decompose is not a positive-definite matrix.