Eigenvalues and Matrix Diagonalization

Introduction to Data Science Linear Algebra Eigenvalues and Matrix Diagonalization

Eigen Pairs

In this section, let us develop our understanding of linear transformation by breaking it down into simpler linear transformations.

Let $T$ be a linear transformation from $\mathbb {R}^n$ to $\mathbb {R}^n$ .

Suppose:

$B$ is a basis of $\mathbb {R}^n$
$V$ is the span of some of the vectors in $B$
$W$ is the span of the remaining vectors in $B$ .

Then:

Any vector in $\mathbb {R}^n$ can be written as the sum of a vector $v$ in $V$ and a vector $w$ in $W$ .
Since $T(v + w) = T(v) + T(w)$ , we can see how $T$ behaves on all of $\mathbb {R}^n$ if we know how it behaves on $V$ and on $W$ .
This decomposition is particularly helpful if $V$ and $W$ are chosen so that $T$ behaves in a simple way on $V$ and on $W$ .

Given such a decomposition of $\mathbb {R}^n$ into the vector spaces $V$ and $W$ , the same idea can be applied to split $V$ and $W$ into lower-dimensional vector spaces and repeat until no more splits are possible. The most optimistic outcome of this procedure would be that we get all the way down to $n$ one-dimensional subspaces and that $T$ acts on each of these subspaces by simply scaling each vector in that subspace by some factor. In other words, we would like to find $n$ vectors $v$ for which $T(v)$ is a scalar multiple of $v$ . From this, the following definition can be obtained.

Mathematical Definition:

An eigenvector $v$ of an $n \times n$ matrix $A$ has the property $Av = \lambda v$ for some $\lambda$ ∈ $\mathbb {R}$ .
$\lambda$ is the eigenvalue of $A$ .
The eigenvector together with its eigenvalue is called eigenpair.

Being more general, an eigenvector is a vector whose direction remains unchanged when a linear transformation is applied to it.

In the figure above, the red vectors are the eigenvectors because their directions didn’t change even after the application of a linear transformation(scaling in this case) while the magenta vector is not an eigenvector as its direction is changed after the linear transformation.

Q: Find the eigen pairs of a $3 \times 3$ matrix $A = \begin{bmatrix} 1.5 & 0 & 1\\-0.5 & 0.5 & -0.5 \\-0.5 & 0 & 0 \end{bmatrix}$

				
					%use s2

// Create a matrix
val A = DenseMatrix(arrayOf
        (doubleArrayOf(1.5, 0.0, 1.0), 
         doubleArrayOf(-0.5, 0.5, -0.5), 
         doubleArrayOf(-0.5, 0.0, 0.0)))
// eigen decomposition
val eigen = Eigen(A, Eigen.Method.QR, 0.0)
println("eigenpairs:")
for (i in 0 until eigen.size()) {
    // get the i-th eigenvalue
    val eigenvalue = eigen.getEigenvalue(i)
    // get the property associated with the eigenvalue
    val property = eigen.getProperty(eigenvalue)
    // get the eigen vector from property
    val eigenVector = property.eigenVector()
    println("eigenvalue: $eigenvalue; eigenvector: $eigenVector")
}

eigenpairs:
eigenvalue: 1.0; eigenvector: [-2.000000, 1.000000, 1.000000] 
eigenvalue: 0.5; eigenvector: [-0.000000, 1.000000, 0.000000] 
eigenvalue: 0.49999999999999994; eigenvector: [-0.000000, 1.000000, 0.000000]

Note: Every non-zero vector is an eigenvector of the identity matrix with eigenvalue 1.

Eigen Decomposition

Let’s say, $A$ is an $n \times n$ having linearly y independent eigenvectors, then one-dimensional subspaces spanned by each of these vectors can be thought of as (not necessarily orthogonal) axes along which $A$ acts by scaling.

Speaking in matrix terms, we can define $Q$ to be the matrix with the eigenvectors of $A$ as columns. Then from the definition of an eigenpair, we have:

$AQ = QD$

where $D$ is a diagonal matrix whose diagonal entries are the eigenvalues (in the order corresponding to the columns of $Q$ ) and whose other entries are zero.

Thus, we can conclude that $A = QDQ^{-1}$ , where $D$ is a diagonal matrix and $A$ is diagonalizable.

S2 makes it easier to find the diagonal matrix $D$ for a given diagonalizable matrix $A$ .

				
					%use s2

// Create a symmetric matrix
val A = DenseMatrix(SymmetricMatrix(
    arrayOf(doubleArrayOf(1.0),
            doubleArrayOf(2.0, 3.0), 
            doubleArrayOf(4.0, 5.0, 6.0), 
            doubleArrayOf(7.0, 8.0, 9.0, 10.0))))
// eigen decomposition
val precision = 1e-4
val decomp = EigenDecomposition(A, precision)
val D = decomp.D() // get the diagonal matrix which contains the eigen values
println("eigen decomposition result: $D")

// verification
val Q = decomp.Q();
val Qt = decomp.Qt();
val B = Q.multiply(D).multiply(Qt); // QDQ' = A = B
println("$A =\n$B")

eigen decomposition result: 4x4
	[,1] [,2] [,3] [,4] 
[1,] 24.638176, 0.000000, 0.000000, 0.000000, 
[2,] 0.000000, -0.038687, 0.000000, 0.000000, 
[3,] 0.000000, 0.000000, -0.513527, 0.000000, 
[4,] 0.000000, 0.000000, 0.000000, -4.085962, 
4x4
	[,1] [,2] [,3] [,4] 
[1,] 1.000000, 2.000000, 4.000000, 7.000000, 
[2,] 2.000000, 3.000000, 5.000000, 8.000000, 
[3,] 4.000000, 5.000000, 6.000000, 9.000000, 
[4,] 7.000000, 8.000000, 9.000000, 10.000000,  =
4x4
	[,1] [,2] [,3] [,4] 
[1,] 1.000000, 2.000000, 4.000000, 7.000000, 
[2,] 2.000000, 3.000000, 5.000000, 8.000000, 
[3,] 4.000000, 5.000000, 6.000000, 9.000000, 
[4,] 7.000000, 8.000000, 9.000000, 10.000000,

Positive Definite Matrices

A symmetric matrix $A$ is called a positive definite matrix if its eigenvalues are all positive. A symmetric matric matrix $A$ is a positive semidefinite matrix if its eigenvalues are all non-negative.

Equivalently, a matrix $A$ is positive semidefinite if $x'Ax \geq 0$ for all $x$ .

Analogously, a symmetric matrix $A$ is called a negative definite matrix if its eigenvalues are all negative. A symmetric matric matrix $A$ is a negative semidefinite matrix if its eigenvalues are all non-positive.

In S2, we have a constructor “PositiveDefiniteMatrixByPositiveDiagonal” that converts a matrix into a symmetric, positive definite matrix, if it is not already so, by forcing the diagonal entries in the eigen decomposition to a small non-negative number.

				
					%use s2
// define a matrix
var A = DenseMatrix(arrayOf(
    doubleArrayOf(1.0, 2.0, 3.0),
    doubleArrayOf(4.0, 5.0, 6.0),
    doubleArrayOf(-1.0, -2.0, -3.0)))

// A_PDM = Positive Definite Matrix of A
val epsilon = 1e-15
val A_PDM = PositiveDefiniteMatrixByPositiveDiagonal(A,epsilon,1.0)
println("POSITIVE DEFINITE MATRIX OF\n\n $A \n\n IS \n\n $A_PDM\n")

POSITIVE DEFINITE MATRIX OF

 3x3
	[,1] [,2] [,3] 
[1,] 1.000000, 2.000000, 3.000000, 
[2,] 4.000000, 5.000000, 6.000000, 
[3,] -1.000000, -2.000000, -3.000000,  

 IS 

 3x3
	[,1] [,2] [,3] 
[1,] 2.286186, 2.413588, 0.605269, 
[2,] 2.413588, 5.529211, 1.135817, 
[3,] 0.605269, 1.135817, 1.284835,

We also have another constructor in S2 named “PositiveSemiDefiniteMatrixNonNegativeDiagonal” which converts a matrix into a symmetric, positive semi-definite matrix, if it is not already so, by forcing the negative diagonal entries in the eigen decomposition to 0.

				
					%use s2
// define a matrix
var A = DenseMatrix(arrayOf(
    doubleArrayOf(1.0, 2.0, 3.0),
    doubleArrayOf(4.0, 5.0, 6.0),
    doubleArrayOf(-1.0, -2.0, -3.0)))

// A_PSDM = Positive Semi-Definite Matrix of A
val epsilon = 100.0
val A_PSDM = PositiveSemiDefiniteMatrixNonNegativeDiagonal(A,epsilon)
println("POSITIVE SEMI-DEFINITE MATRIX OF\n\n $A \n\n IS \n\n $A_PSDM\n")

POSITIVE SEMI-DEFINITE MATRIX OF

 3x3
	[,1] [,2] [,3] 
[1,] 1.000000, 2.000000, 3.000000, 
[2,] 4.000000, 5.000000, 6.000000, 
[3,] -1.000000, -2.000000, -3.000000,  

 IS 

 3x3
	[,1] [,2] [,3] 
[1,] 7.100233, -0.000000, -0.000000, 
[2,] -0.000000, 7.100233, -0.000000, 
[3,] -0.000000, -0.000000, 7.100233,

Spectral Theorem:

If $A$ is an $n \times n$ symmetric matrix, then $A$ is orthogonally diagonalizable, which means $A$ has $n$ eigenvectors which are pairwise orthogonal.

Conversely, every orthogonally diagonalizable matrix is symmetric.

Polar Decomposition

Before proceeding with the definition of polar decomposition, let us make ourselves familiar with “Gram Matrix”:

Gram Matrix: If $A$ is an $m \times n$ matrix, then $A'A$ is its Gram matrix. The Gram matrix of $A$ is always positive semidefinite.

Q: Determine the Gram Matrix of $A = \begin{bmatrix}1 & 2 & 3\\4 & 5 & 6\\-1 & -2 & -3\end{bmatrix}$ .

				
					%use s2
// define a matrix
var A = DenseMatrix(arrayOf(
    doubleArrayOf(1.0, 2.0, 3.0),
    doubleArrayOf(4.0, 5.0, 6.0),
    doubleArrayOf(-1.0, -2.0, -3.0)))

// A_T = Tranpose of A
val A_T = A.t()

//A_T*A = Gram Matrix of A
val A_GM = A_T.multiply(A)
println("THE GRAM MATRIX OF \n $A \n\n IS \n\n $A_GM\n")

THE GRAM MATRIX OF 
 3x3
	[,1] [,2] [,3] 
[1,] 1.000000, 2.000000, 3.000000, 
[2,] 4.000000, 5.000000, 6.000000, 
[3,] -1.000000, -2.000000, -3.000000,  

 IS 

 3x3
	[,1] [,2] [,3] 
[1,] 18.000000, 24.000000, 30.000000, 
[2,] 24.000000, 33.000000, 42.000000, 
[3,] 30.000000, 42.000000, 54.000000,

Let’s define matrix $\sqrt{A'A}$ to be $QD^{1/2}Q'$ , where $QDQ'$ is the diagonalization of $A'A$ and $D^{1/2}$ is the matrix obtained by taking the diagonal entries of $D$ . Then $\sqrt{A'A}$ is symmetric and satisfies:

$\sqrt{A'A} \sqrt{A'A} = QD^{1/2}Q' QD^{1/2}Q' = A'A$

The matrix $\sqrt{A'A}$ is easier to understand than $A$ because it is symmetric and positive definite. However, it transforms space in nearly same way as $A$ . If $x \in \mathbb {R}^n$ , then:

$\left|Ax\right|^2 = x'A' Ax = x' \sqrt{A'A} \sqrt{A'A} x = \left| \sqrt{A'A}x \right|^2$

In other words, for all $x$ , the images of $x$ under $A$ and under $\sqrt{A'A}$ have equal norm. This means that for each $x \in \mathbb {R}^n$ , there is an orthogonal transformation from the range of $\sqrt{A'A}$ to the range of $A$ which sends
$Ax$ to $\sqrt{A'A}x$ . That means, this orthogonal transformation is the same for all $x$ .

Fig: From the above figure, we can observe that the grid-line images under $A$ and $\sqrt{A'A}$ are of the same shape.

$\implies$ They are related by orthogonal transformation.

This orthogonal transformation can be extended to an orthogonal transformation on $\mathbb {R}^n$ even if the range of $\sqrt{A'A}$ is not all of $\mathbb{R}^n$ . Hence, we arrive at a polar decomposition.

Polar Decomposition:

For any $n \times n$ matrix $A$ , there exists an orthogonal matrix $R$ such that

$A = R \sqrt{A'A}$

Let us implement the polar decomposition using S2.

				
					%use s2
//define matrix A
var A = DenseMatrix(arrayOf(
    doubleArrayOf(4.0, 5.0),
    doubleArrayOf(6.0, 7.0)))
 
//AT is transpose of A
var AT = A.t()
//ATA = AT*A
var ATA = AT.multiply(A)
println(ATA)

2x2
	[,1] [,2] 
[1,] 52.000000, 62.000000, 
[2,] 62.000000, 74.000000,

				
					%use s2
// define a matrix sqrt(A'A) from previous results
var RATA = DenseMatrix(arrayOf(
    doubleArrayOf((52.0).pow(0.5), (62.0).pow(0.5)),
    doubleArrayOf((62.0).pow(0.5), (74.0).pow(0.5))))

//Matrix A
var A = DenseMatrix(arrayOf(
    doubleArrayOf(4.0, 5.0),
    doubleArrayOf(6.0, 7.0)))

//calculate inverse of sqrt(A'A)
val INVRATA = Inverse(RATA)

// R = A*Inverse(sqrt(A'A))
val R = A.multiply(INVRATA)
println("R =\n\n $R\n")

//verification
var B = R.multiply(RATA)
println("\nR*sqrt(A'A) = \n\n $B \n\nSAME AS\n\n A = \n\n$A\n")

R =

 2x2
	[,1] [,2] 
[1,] -153.822883, 141.380679, 
[2,] -108.655461, 100.269860, 


R*sqrt(A'A) = 

 2x2
	[,1] [,2] 
[1,] 4.000000, 5.000000, 
[2,] 6.000000, 7.000000,  

SAME AS

 A = 

2x2
	[,1] [,2] 
[1,] 4.000000, 5.000000, 
[2,] 6.000000, 7.000000,

Previous Topic

Back to Lesson

Next Topic

Very helpful, I've needed this explanation for 15 years!

I believe there is a mistake in the above - most likely for the example of removing drift has been…

Great description! Thanks a lot!

I dont see the trail version of the license in 3.6.0 release. Can you please check? I also paid the…

[…] Transforming a QP Problem to an SOCP Problem [Numerical Method] […]