Linear Algebra
Review & Recent Progress

Gene Cheung

York University
October 17, 20191 / 17

Outline

What is linear algebra?
Why should I care?
System of Linear Equations
Eigen-decomposition

2 / 17

What is Linear Algebra?

At the risk of over-simplification, linear algebra solves two
problems:

12Matrix A and vector y can also be complex: A ∈ Cn×n , y ∈ Cn .
One can also pose the generlized eigenvalue problem: Av = λBv.
3 / 17

What is Linear Algebra?

At the risk of over-simplification, linear algebra solves twoproblems:
1

System of linear equations:
given matrix1 A ∈ Rn×n and vector y ∈ Rn , find x
Ax = y12(1)Matrix A and vector y can also be complex: A ∈ Cn×n , y ∈ Cn .
One can also pose the generlized eigenvalue problem: Av = λBv.
3 / 17

What is Linear Algebra?

At the risk of over-simplification, linear algebra solves two
problems:
1

System of linear equations:
given matrix1 A ∈ Rn×n and vector y ∈ Rn , find x
Ax = y2Eigen-decomposition:
given matrix A ∈ Rn×n , find eigen-pair2 (λ, v), where λ ∈ R
and v ∈ Rn , such that
Av = λv

12(1)

(2)

Matrix A and vector y can also be complex: A ∈ Cn×n , y ∈ Cn .
One can also pose the generlized eigenvalue problem: Av = λBv.
3 / 17

Why should I care? (Part I)
Ax = y yields solution to the least-square problem3 :4 / 17

Why should I care? (Part I)
Ax = y yields solution to the least-square problem3 :1Overdetermined equation: For “tall” matrix H, solve
min ky − Hxk22x>=⇒ H H x∗ = H> y(3)(4)

4 / 17

Why should I care? (Part I)
Ax = y yields solution to the least-square problem3 :1Overdetermined equation: For “tall” matrix H, solve
min ky − Hxk22x>=⇒ H H x∗ = H> y2

Underdetermined equation: For “fat” matrix H, solve
min kxk22
s.t. y = Hxx


 

x
0
I − 21 H>
=⇒=µy
H0(3)
(4)

(5)
(6)

4 / 17

Why should I care? (Part I)
Ax = y yields solution to the least-square problem3 :1

Overdetermined equation: For “tall” matrix H, solve
min ky − Hxk22x>=⇒ H H x∗ = H> y2

3

Underdetermined equation: For “fat” matrix H, solve
min kxk22s.t. y = Hxx x0
I − 21 H>
=⇒=µy
H
0(3)
(4)

(5)
(6)

Regularization:=⇒

min ky − Hxk22 + λkAxk22x

H> H + λA> A x∗ = H> y(7)(8)

3 http://eeweb.poly.edu/iselesni/lecture notes/least squares/index.html
4 / 17Why should I care? (Part II)
Evaluate stability of Ax = y.
Condition number κ(A) of coefficient matrix A:|λmax (A)|κ(A) =|λmin (A)|(9)5 / 17

Why should I care? (Part II)
Evaluate stability of Ax = y.
Condition number κ(A) of coefficient matrix A:|λmax (A)|κ(A) =|λmin (A)|(9)

E-optimality criteria:max λmin (A)A∈S(10)5 / 17

Why should I care? (Part II)
Evaluate stability of Ax = y.
Condition number κ(A) of coefficient matrix A:|λmax (A)|κ(A) =|λmin (A)|(9)

E-optimality criteria:max λmin (A)
A∈S

(10)

Compute graph frequencies 4 :1Define variation operator Φ on graph, e.g., graph Laplacianmatrix L:
L=D−W
(11)

2

Compute Fourier modes for L via eigen-decomposition.
L = VΛV>(12)

where V contains eigenvectors v1 , v2 , . . . as columns, andΛ = diag(λ1 , λ2 , . . .).
4 A. Ortega et al., “Graph signal processing: Overview, challenges, and applications,” Proceedings of the IEEE, vol.
106, no. 5, pp. 808–828, 2018.
5 / 17

Graph Frequencies: Discrete Cosine Transform (GCT)
Cosines are eigenfunctions of differential operator:
12Construct graph Laplacian L for line graph with weights 1.
Compute eigenvectors for L.

Figure: Eigenvectors of line graph Laplacian matrix L = D − W

6 / 17

Graph Frequencies: Discrete Fourier Transform (GFT)

Figure: 1st eigenvector of graph Laplacian matrix L = D − W
7 / 17

Graph Frequencies: Discrete Fourier Transform (GFT)

Figure: 2nd eigenvector of graph Laplacian matrix L = D − W
8 / 17

Graph Frequencies: Discrete Fourier Transform (GFT)

Figure: 3rd eigenvector of graph Laplacian matrix L = D − W

9 / 17

Graph Frequencies: Discrete Fourier Transform (GFT)

Figure: 9th eigenvector of graph Laplacian matrix L = D − W

10 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.
1

Gaussian elimination.

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.
1

Gaussian elimination.

2

If A is symmetric and positive semi-definite (PSD), thenCholesky decomposition.

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.
1

Gaussian elimination.

2

If A is symmetric and positive semi-definite (PSD), thenCholesky decomposition.

3

If A is diagonally dominant, Jacobi method.

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.
1

Gaussian elimination.

2

If A is symmetric and positive semi-definite (PSD), thenCholesky decomposition.

3

If A is diagonally dominant, Jacobi method.

4

If A is diagonally dominant, or symmetric and PD,
Gauss-Seidel method.

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.
1

Gaussian elimination.

2

If A is symmetric and positive semi-definite (PSD), thenCholesky decomposition.

3

If A is diagonally dominant, Jacobi method.

4

If A is diagonally dominant, or symmetric and PD,
Gauss-Seidel method.
Successive approximations using Krylov subspace methods:

5

5M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving a System of Linear Equations
Q: How to solve Ax = b?
A: Many methods.
1

Gaussian elimination.

2

If A is symmetric and positive semi-definite (PSD), thenCholesky decomposition.

3

If A is diagonally dominant, Jacobi method.

4

If A is diagonally dominant, or symmetric and PD,
Gauss-Seidel method.
Successive approximations using Krylov subspace methods:

5

If A symmetric and PD, Conjugate Gradient5 .

5
M. R. Hestenes and E. Stiefel, “Methods of conjugate gradients for solving linear systems,” Journal of
Research of the National Bureau of Standards, 1952, vol. 49, no. 1.
11 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.
1Power method.

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.
1Power method.

2

QR algorithm.

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.
1Power method.

2

QR algorithm.
Successive approximations using Krylov subspace methods:

3

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.
1Power method.

2

QR algorithm.
Successive approximations using Krylov subspace methods:

3

Arnoldi iteration.

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.
1Power method.

2

QR algorithm.
Successive approximations using Krylov subspace methods:

3

Arnoldi iteration.
Divide-and-conquer algorithm.

12 / 17

Solving the Eigen-decomposition Problem

Q: How to solve Av = λv?
1

How to compute eigenvalues λ’s?

2

How to compute eigenvectors v’s?

A: Many methods.
1Power method.

2

QR algorithm.
Successive approximations using Krylov subspace methods:

3

Arnoldi iteration.
Divide-and-conquer algorithm.
If A Hermetian, Lanczos algorithm.

12 / 17

Lower-bounding Smallest Eigenvalue λmin
Recall E-optimality criteria:max λmin (A)
A∈S

6(13)

R. S. Varga, Gersgorin and His Circles Springer, 2004
13 / 17

Lower-bounding Smallest Eigenvalue λmin
Recall E-optimality criteria:max λmin (A)A∈S(13)

Gershgorin Circle Theorem (GCT) 6

6

R. S. Varga, Gersgorin and His Circles Springer, 2004
13 / 17

Lower-bounding Smallest Eigenvalue λmin
Recall E-optimality criteria:max λmin (A)A∈S(13)

Gershgorin Circle Theorem (GCT) 6
Given matrix
PA, Gershgorin disc i, Φi , has centre ci = Aii andradius ri = j6=i |Aij |.

6

R. S. Varga, Gersgorin and His Circles Springer, 2004
13 / 17

Lower-bounding Smallest Eigenvalue λmin
Recall E-optimality criteria:max λmin (A)
A∈S

(13)

Gershgorin Circle Theorem (GCT) 6
Given matrix
PA, Gershgorin disc i, Φi , has centre ci = Aii andradius ri = j6=i |Aij |.
For each eigenvalue λ of A, there exists a disc Φi such that:ci − ri ≤ λ ≤ ci + ri

6(14)R. S. Varga, Gersgorin and His Circles Springer, 2004
13 / 17

Lower-bounding Smallest Eigenvalue λmin
Recall E-optimality criteria:max λmin (A)A∈S(13)

Gershgorin Circle Theorem (GCT) 6
Given matrix
PA, Gershgorin disc i, Φi , has centre ci = Aii andradius ri = j6=i |Aij |.
For each eigenvalue λ of A, there exists a disc Φi such that:ci − ri ≤ λ ≤ ci + ri(14)Compute a lower bound for λmin (A) without computingeigen-pairs:
λ−min (A) = min ci − ri ≤ λmin (A)i6(15)

R. S. Varga, Gersgorin and His Circles Springer, 2004
13 / 17

Lower-bounding Smallest Eigenvalue λmin

Q: How tight is GCT lower bound λ−min (A) for λmin (A)?

7
Y. Bai et al., “Fast Graph Sampling Set Selection Using Gershgorin Disc Alignment,” submitted to IEEE
Transactions on Signal Processing, July 2019.
14 / 17

Lower-bounding Smallest Eigenvalue λmin

Q: How tight is GCT lower bound λ−min (A) for λmin (A)?
A: Consider instead similar transform7 :B = SAS−1(16)where S = diag(s1 , s2 , . . .).

7
Y. Bai et al., “Fast Graph Sampling Set Selection Using Gershgorin Disc Alignment,” submitted to IEEE
Transactions on Signal Processing, July 2019.
14 / 17

Lower-bounding Smallest Eigenvalue λmin

Q: How tight is GCT lower bound λ−min (A) for λmin (A)?
A: Consider instead similar transform7 :B = SAS−1(16)where S = diag(s1 , s2 , . . .).
B has same eigenvalues as A.

7
Y. Bai et al., “Fast Graph Sampling Set Selection Using Gershgorin Disc Alignment,” submitted to IEEE
Transactions on Signal Processing, July 2019.
14 / 17

Lower-bounding Smallest Eigenvalue λmin

Q: How tight is GCT lower bound λ−min (A) for λmin (A)?
A: Consider instead similar transform7 :B = SAS−1(16)where S = diag(s1 , s2 , . . .).
B has same eigenvalues as A.
If A is a generalized graph Laplacian with positive edges, ∃Ssuch that λ−min (B) = λmin (B) = λmin (A).

7Y. Bai et al., “Fast Graph Sampling Set Selection Using Gershgorin Disc Alignment,” submitted to IEEE
Transactions on Signal Processing, July 2019.
14 / 17

Lower-bounding Smallest Eigenvalue λmin
Partition matrix A:
A=A11 A12A>
12 A22

(17)

8G. Cheung et al., “Robust Semisupervised Graph Classifier Learning with Negative Edge Weights,” IEEE
Transactions on Signal and Information Processing over Networks, vol. 4, no.4, pp.712-726, December 2018.
15 / 17

Lower-bounding Smallest Eigenvalue λmin
Partition matrix A:
A=A11 A12A>
12 A22

(17)

Define Schur Complement:−1
A/A11 = A22 − A>
12 A11 A12(18)8G. Cheung et al., “Robust Semisupervised Graph Classifier Learning with Negative Edge Weights,” IEEE
Transactions on Signal and Information Processing over Networks, vol. 4, no.4, pp.712-726, December 2018.
15 / 17

Lower-bounding Smallest Eigenvalue λmin
Partition matrix A:
A=A11 A12A>
12 A22

(17)

Define Schur Complement:−1
A/A11 = A22 − A>
12 A11 A12(18)Haynsworth Inertia Additivity formula:
In(A) = In(A11 ) + In(A/A11 )(19)

8G. Cheung et al., “Robust Semisupervised Graph Classifier Learning with Negative Edge Weights,” IEEE
Transactions on Signal and Information Processing over Networks, vol. 4, no.4, pp.712-726, December 2018.
15 / 17

Lower-bounding Smallest Eigenvalue λmin
Partition matrix A:
A=A11 A12A>
12 A22

(17)

Define Schur Complement:−1
A/A11 = A22 − A>
12 A11 A12(18)Haynsworth Inertia Additivity formula:
In(A) = In(A11 ) + In(A/A11 )(19)

If A is real and symmetric, can compute tight lower bound forλmin using a recursive procedure plus shifting (+µI) 8 .

8G. Cheung et al., “Robust Semisupervised Graph Classifier Learning with Negative Edge Weights,” IEEE
Transactions on Signal and Information Processing over Networks, vol. 4, no.4, pp.712-726, December 2018.
15 / 17

Computing Extreme Eigenvectors

locally optimal block preconditioned conjugated gradient(LOBPCG) method 9 :Compute extreme eigen-pairs for symmetric, PD matrices.
Cost per iteration competitive with Lanczos.
Can directly take advantage of preconditioning.
Can benefit from warm start.

9A. Knyazev, “Toward the optimal preconditioned eigensolver: Locally optimal block preconditioned conjugate
gradient method,” SIAM journal on scientific computing, 23(2):517–541, 2001.
16 / 17

Summary

Linear algebra solves:
12Ax = yAv = λv.

Mature algorithms to solve Ax = y using Krylov methods.
Fast algorithms to find extreme eigen-pairs (λi , vi ).
New algorithms to find tight lower bounds for λmin .

17 / 17