17. Rapidly Decreasing Singular Values
2.71M subscribers
SUBSCRIBED
MIT 18.065 Matrix Methods in Data Analysis, Signal Processing, and Machine Learning, Spring 2018 Instructor: Alex Townsend View the complete course: https://ocw.mit.edu/18-065S18YouTube Playlist: https://www.youtube.com/playlist?list...Professor Alex Townsend gives this guest lecture answering the question 'Why are there so many low rank matrices that appear in computational math?' Working effectively with low rank matrices is critical in image compression applications. License: Creative Commons BY-NC-SA More information at https://ocw.mit.edu/termsMore courses at https://ocw.mit.edu
SHOW LESS
00:00
00:22
GILBERT STRANG: So let me use the mic
00:26
to introduce Alex Townsend, who taught here at MIT--
00:31
taught Linear Algebra 18.06 very successfully.
00:35
And then now he's at Cornell on the faculty,
00:38
still teaching very successfully.
00:40
And he was invited here yesterday
00:43
for a big event over in Engineering.
00:47
And he agreed to give a talk about a section of the book--
00:54
section 4.3--
00:57
which, if you look at it, you'll see is all about his work.
01:00
And now you get to hear from the creator himself.
01:04
OK.
01:05
01:10
ALEX TOWNSEND: OK.
01:11
Thanks.
01:11
Thank you, Gil.
01:12
Thank you for inviting me here.
01:14
I hope you're enjoying the course.
01:15
Today I want to tell you a little about why
01:19
there so many matrices that are low rank in the world.
01:24
So as computational mathematicians--
01:26
Gil and myself-- we come across low-rank matrices all the time.
01:30
And we started wondering, as a community, why?
01:36
What is it about the problems that we are looking at?
01:41
What makes low-rank matrices appear?
01:44
And today I want to give you that story--
01:46
or at least an overview of that story.
01:49
So for this class, x is going to be n by n real matrix.
01:57
So nice and square.
01:59
And you already know, or are very comfortable with,
02:02
the singular values of a matrix.
02:05
So the singular values of a matrix, as you know,
02:09
are a sequence of numbers that are monotonically
02:15
non-increasing that tell us all kinds of things
02:21
about the matrix x.
02:22
02:26
For example, the number of nonzero singular values
02:30
tell us the rank of the matrix x.
02:33
And they also, you probably know, tell us how well a matrix
02:37
x can be approximated by a low-rank matrix.
02:42
So let me just write two facts down that you already
02:45
are familiar with.
02:47
So here's a fact--
02:50
that, if I look at the number of non-zero singular values in x--
02:57
so I'm imagining there's going to be k non-zero singular
03:01
values--
03:01
03:06
then we can say a few things about x.
03:09
For example, the rank of x, as we know, is k--
03:16
the number of non-zero singular values.
03:19
But we also know from the SVD that we can decompose x
03:25
into a sum of rank 1 matrices--
03:29
in fact, the sum of k of them.
03:33
So because x is rank k, we can write down
03:37
a low-rank representation for x, and it involves k terms,
03:45
like this.
03:47
Each one of these vectors here is a column vector.
03:52
So if I draw this pictorially, this guy
03:57
looks like this, right?
03:58
And we have k of them.
04:01
So because x is rank k, we can write x as a sum of k
04:06
rank 1 matrices.
04:08
And we also have an initial fact that we already know--
04:11
that the dimension of the column space of x
04:16
is equal to k, and the same with the row space.
04:20
So the column space of x equals the row space of x--
04:27
the dimension-- and they all equal k.
04:35
And so there are three facts we can determine
04:37
from looking at this sequence of singular values of a matrix x.
04:44
Of course, the singular value sequence is unique.
04:47
X defines its own singular values.
04:50
04:54
What we're interested in here is, what makes x?
04:58
What are the properties of x that make
05:00
sure that the singular values have a lot
05:02
of zeros in that sequence?
05:05
Can we try to understand what kind of x makes that happen?
05:09
05:12
And we really like matrices that have a lot of zeros
05:16
here, for the following reason--
05:18
05:22
we say x is low rank if the following holds, right?
05:30
Because if we wanted to send x to our friend--
05:33
we're imagining x as picture where each entry
05:36
is a pixel of that image.
05:40
If that matrix-- that image-- was low rank,
05:43
we could send the picture to our friend in two ways.
05:48
We could send one every single entry of x.
05:53
And for us to do that, we would have
05:55
to send n squared pieces of information,
05:58
because we'd have to send every entry.
06:01
But if x is sufficiently low rank,
06:03
we could also send our friend the vectors--
06:07
u, u1, v1, uk, up to vk.
06:12
And how much pieces of data would we
06:15
have to send our friend to get x to them
06:18
if we sent in the low-rank form?
06:20
Well, there's 2n here, 2n here numbers.
06:27
There's k of them.
06:28
So we'd have to send 2kn numbers.
06:33
And we strictly say a matrix is low
06:36
rank if it's more efficient to send x to our friend
06:41
in low-rank form then in full-rank form.
06:46
So this, of course, by a little calculation,
06:49
just shows us that, provided the rank is
06:53
less than half the size of the matrix,
06:56
we are calling the matrix low rank.
07:00
Now, often, in practice, we demand more.
07:05
07:09
We demand that k is much smaller than this number,
07:15
so that it's far more efficient to send our friend the matrix x
07:19
in low-rank form than in full-rank form.
07:23
So the colloquial use of the word low rank
07:26
is kind of this situation.
07:29
But this is the strict definition of it.
07:31
07:34
So what do low-rank matrices look like?
07:39
And to do that, I have some pictures for you.
07:43
I have some flags--
07:44
the world flags.
07:47
So these are all matrices x--
07:50
these examples-- because their flags happen to not be square.
07:54
I hope you can all see this.
07:56
But the top row here are all matrices
08:01
that are extremely low rank.
08:04
For example, the Austria flag--
08:06
if you want to send that to your friend,
08:08
that matrix is of rank 1.
08:11
So all you have to do is send your friend two vectors.
08:14
You have to tell your friend the column space and the row space.
08:18
And there's only the dimensions of one of both.
08:21
For the English flag, you need to send them two column
08:24
vectors and two row vectors--
08:27
u1, v1, u2 and v2.
08:31
And as we go down this row, they get slowly fuller and fuller
08:35
rank.
08:36
So the Japanese flag, for example,
08:38
is low rank but not that small.
08:43
The Scottish flag is essentially full rank.
08:45
So it's very inefficient to send your friend the Scottish flag
08:50
in low-rank form.
08:51
You're better off sending almost every single entry.
08:55
So what do low-rank matrices look like?
08:58
09:11
Well, if the matrix is extremely low rank,
09:15
like rank 1, then when you look at that matrix--
09:18
like here, like the flag--
09:19
it's highly aligned with the coordinates--
09:25
with the rows and columns.
09:27
So if it's rank 1, the matrix is highly aligned--
09:33
like the Austria flag.
09:34
09:41
And of course, as we add in more and more rank here,
09:43
the situation gets a bit blurry.
09:46
For example, once we get into the medium rank situation,
09:49
which is a circle, it's very hard
09:51
to see that the circle is actually, in fact, low rank.
09:56
But what I'm going to do was try to understand
09:58
why the Scottish flag or diagonal patterns--
10:03
particularly a bad example for low rank.
10:07
So I'm going to take the triangular flag
10:12
to examine that more carefully.
10:15
So the triangular flag looks like--
10:19
I'll take a square matrix and I'll color in the bottom half.
10:27
So this matrix is the matrix of ones below the diagonal.
10:31
10:38
And I'm interested in this matrix and, in particular,
10:41
its singular values, to try to understand
10:43
why diagonal patterns are not particularly
10:47
useful for low-rank compression.
10:51
And this matrix of all ones has a really nice property that,
10:57
if I take its inverse, it looks a lot like--
11:01
getting close to Gil's favorite matrix.
11:04
So if I take the inverse of this matrix--
11:08
it has an inverse because it's got ones on the diagonal--
11:12
then its inverse is the following matrix,
11:21
which people familiar with finite difference schemes
11:24
will notice the familiarity between that
11:27
and the first order finite difference approximation.
11:32
In particular, if I go a bit further and times
11:35
two of these together, and do this,
11:39
then this is essentially Gil's favorite matrix,
11:45
except one entry happens to be different--
11:49
11:52
ends up being this matrix, which is
11:55
very close to the second order, central, finite difference
11:59
matrix.
12:01
And people have very well studied
12:02
that matrix and know its eigenvalues,
12:06
its singular values--
12:07
they know everything about that matrix.
12:10
And you'll remember that if we know
12:12
the eigenvalues of a matrix, like x transpose x,
12:18
we know the singular values of x.
12:21
So this allows us to show, by the fact
12:25
that we know that, that the singular values of this matrix
12:32
are not very amenable to low rank.
12:34
They're all non-zero, and they don't even decay.
12:39
So I'm getting this from--
12:41
12:44
I rang up Gil, and Gil tells me these numbers.
12:47
12:54
That allows us to work out exactly what the singular
12:57
values of this matrix are, from the connection
12:59
to finite differences.
13:02
And so we can understand why this is not
13:04
good by looking at the singular values.
13:06
So the first singular value of x from this expression
13:10
is going to be approximately 2n over pi.
13:15
And from this expression, again, for the last guy--
13:19
the last singular value of x is going
13:23
to be approximately a half.
13:26
So these singular values are all large.
13:28
They're not getting close to zero.
13:31
If I plotted these singular values on a graph--
13:37
so here's the first singular value, the second,
13:40
and the n-th--
13:43
then what would the graph look like?
13:45
Well, plot these numbers.
13:47
13:51
Divide by this guy so that they all
13:54
are bounded between 1 and 0 because of the normalization,
13:59
because I divided by sigma 1 of x.
14:03
And so we can plot them, and they will
14:05
look like this kind of thing.
14:11
14:15
This number happens to be here where
14:16
they come to be pi over 4n, which
14:21
is me dividing this number by this number, approximately.
14:27
So triangular patterns are extremely bad for low rank.
14:32
We need things-- or we at least intuitively think
14:35
that we need things-- aligned with the rows and columns,
14:39
but the circle case happens to also be low rank.
14:45
And so what happened to the Japanese flag?
14:49
14:54
Why is the Japanese flag convenient for low rank?
14:59
Well it's the fact that it's a circle,
15:03
and there's lots of symmetry in a circle.
15:06
So if I try to look at the rank of a circle, the Japanese flag,
15:13
then I can bound this rank by decomposing the Japanese flag
15:21
into two things.
15:24
So this is going to be less than or equal to the rank of sum
15:29
of two matrices, and I'll do it so that the decomposition works
15:33
out.
15:34
I have the circle.
15:35
I'm going to cut out a rank one piece that lives
15:40
in the middle of this circle.
15:41
15:47
OK?
15:47
And I'm going to cut out a square from the interior
15:51
of that circle.
15:54
OK?
15:55
And I can figure out-- of course the rank is just bounded
15:58
by the sum of those two ranks.
16:00
This guy is bounded by rank one because it's highly
16:04
aligned with the grid.
16:05
16:09
So this guy is bounded by rank one.
16:11
So this thing here plus 1.
16:23
16:26
And now I have to try to understand
16:29
the rank of this piece.
16:32
Now this piece has lots of symmetry.
16:35
For example, we know that the rank of that matrix
16:39
is the dimension of the column space
16:43
and the dimension of the row space.
16:46
So when we look at this matrix, because of symmetry,
16:49
if I divide this matrix in half along the columns,
16:55
all the columns on the left appear on the right.
16:58
So for example, the rank of this matrix
17:01
is the same as the rank of that matrix
17:04
because I didn't change the column space.
17:07
OK?
17:08
Now I go again and divide along the rows,
17:13
and now the row dimension of this matrix
17:17
is the same as the top half, because as I wipe out those,
17:20
I didn't change the dimension of the row space
17:23
because the rows are the same top-bottom.
17:26
And so this becomes the rank of that tiny little matrix there.
17:30
And because it's small, it won't have too large a rank.
17:35
So this is definitely less than-- if I divide that up,
17:42
a little guy here looks like that plus the other guy that
17:50
looks like that plus 1.
18:00
And so of course the row space of this matrix cannot be very
18:07
high because this is a very thin matrix.
18:10
There's lots of zeros in that matrix, only a few ones.
18:13
And so you can go along and do a bit
18:15
of trig to try to figure out how many rows are
18:19
non-zero in this matrix.
18:22
And a bit of trig tells you--
18:25
well it depends on the radius of this original circle.
18:29
So if I make the original radius r of this Japanese flag,
18:34
then the bound that you end up getting will be,
18:38
for this matrix, r 1 minus square root 2 over 2
18:43
for this guy.
18:44
That's a bit of trig.
18:46
I've got to make sure that's an integer.
18:48
And then again, here it's the same but for the column space.
18:52
So this is me just doing trig.
18:53
18:56
OK?
18:56
And that's bound on the rank.
18:57
It happens to be extremely good.
18:59
And if you work out what that rank is and try to look back,
19:03
you will find it's extremely efficient
19:05
to send the Japanese flag to your friend in low rank form,
19:10
because it's not full rank because these numbers are
19:12
so small.
19:13
So this comes out to be, like, approximately 1/2 r plus 1.
19:21
So much smaller than what you would expect,
19:23
because remember, a circle is almost the anti-version version
19:27
of a line with the grid, but yet, it's still low rank.
19:32
OK.
19:35
Now most matrices that we come up
19:39
with in computational math are not exactly of finite rank.
19:44
They are of numerical rank.
19:48
And so I'll just define that.
19:51
So the numerical rank of a matrix
19:58
is very similar to the rank, except we allow ourselves
20:01
a little bit of wiggle room when we define it,
20:04
and so that amount of wiggle room will be of parameter
20:09
called tol called epsilon.
20:12
That's a tolerance.
20:13
I'm thinking of epsilon as a tolerance.
20:16
That's the amount of wiggle room I'm going to give myself.
20:21
OK.
20:22
And we say that the numerical rank--
20:27
I'll put an epsilon there to denote numerical rank--
20:31
is k.
20:34
k is the first singular value, or the last singular value,
20:37
above epsilon.
20:39
In the following sense, I'm copying the definition above
20:42
but with epsilons instead of zeros.
20:45
If this singular value is less than epsilon, relatively,
20:52
and the kth one was not below.
20:56
So k plus 1 is the first singular value below epsilon
21:00
in this relative sense.
21:03
So of course the rank of 0x, if that was defined,
21:10
is the same as the rank of x.
21:13
OK?
21:14
So this is just allowing ourselves some wiggle room.
21:17
But this is actually what we're interested more in practice.
21:20
All right?
21:21
I don't want to necessarily send my friend
21:23
the flag to exact precision.
21:26
I would actually be happy to send
21:27
my friend the flag up to 16 digits of precision,
21:31
for example.
21:32
They're not going to tell the difference between those two
21:34
flags.
21:35
And if I can get away with compressing the matrix
21:39
a lot more once I have a little bit of wiggle room,
21:42
that would be a good thing.
21:44
So we know from the Eckart and Young
21:58
that the singular values tell us how well we can approximate
22:02
x by a low-rank matrix.
22:05
In particular, we know that the k plus 1 singular value of x
22:13
tells us how well x can be approximated by a rank k
22:18
matrix.
22:19
OK?
22:20
For example, when the rank was exactly k, the sigma k plus 1
22:26
was 0, and then this came out to be 0
22:29
and we found that x was exactly a rank k matrix.
22:33
Here, because we have the wiggle room, the epsilon,
22:36
we get an approximation, not an exact.
22:39
So this is telling us how well we can approximate
22:44
x by a rank k matrix.
22:47
22:50
OK?
22:51
That's what the singular values are telling us.
22:54
And so this allows us to try our best to compress matrices
22:59
but use low-rank approximation rather
23:03
than doing things exactly.
23:05
23:07
And of course, on a computer, when we're using floating point
23:11
arithmetic, or on a computer because we always round numbers
23:16
to the nearest 16-digit number, if epsilon was 16 digits,
23:21
your computer wouldn't be able to tell
23:23
the difference between x or x the rank k
23:29
approximation if this number satisfied this expression.
23:35
Your computer would think of x and xk as the same matrix
23:38
because it would inevitably round both
23:42
to epsilon, within epsilon.
23:45
OK.
23:46
So what kind of matrices are numerically of low rank?
23:49
24:03
Of course all low-rank matrices are numerically of low rank
24:08
because the wiggle room can only help you
24:16
but it's far more than that.
24:19
There are many full-rank matrices--
24:21
matrices that don't have any singular values that are zero--
24:24
but the singular values decay rapidly to zero.
24:27
That are full-rank matrices with low numerical rank because
24:32
of the wiggle room.
24:33
So for example, here is the classic matrix
24:38
that fits this regime.
24:43
If I give you this, this is called the Hilbert matrix.
24:45
24:51
This is a matrix that happens to have
24:53
extremely low numerical rank but it's actually
24:57
full rank, which means that I can approximate H by a rank k
25:05
matrix where k is quite small very well,
25:08
provided you give me some wiggle room,
25:10
but it's not a low-rank matrix in the sense
25:13
that if epsilon was zero here, you didn't allow me
25:16
the wriggle room, all the singular values
25:18
of this matrix are positive.
25:20
So it's of low numerical rank but it's not a low-rank matrix.
25:28
The other classical example which
25:32
motivated a lot of the research in this area
25:35
was the Vandermonde matrix.
25:37
So here is the Vandermonde matrix.
25:39
25:48
An n by n version of it.
25:50
Think of the xi's as real.
25:52
25:55
And this is Vandermonde.
25:57
26:02
This is the matrix that comes up when
26:03
you try to do polynomial interpolation at real points.
26:08
This is an extremely bad matrix to deal with because it's
26:13
numerically low rank, and often, you actually
26:17
want to solve a linear system with this matrix.
26:21
And numerical low rank implies that it's extremely hard
26:24
to invert, so numerical low rank is not always good for you.
26:32
OK?
26:33
Often, we want the inverse, which exists,
26:42
but it's difficult because V has low numerical rank.
26:56
27:03
OK.
27:04
So people have been trying to understand
27:06
why these matrices are numerically
27:09
of low rank for a number of years,
27:12
and the classic reason why there are
27:16
so many low-rank matrices is because the world is smooth,
27:21
as people say.
27:22
They say, the world is smooth.
27:25
That's why matrices are of numerical low rank.
27:32
And to illustrate that point, I will do an example.
27:38
So this is classically understood
27:41
by a man called Reade in 1983, and this
27:50
is what his reason was.
27:51
I have a picture of John Reade.
27:54
He's not very famous, so I try to make
27:56
sure his picture gets around.
28:00
He's playing the piano.
28:01
It's, like, one of the only pictures I could find of him.
28:04
So what is in this reason?
28:06
Why do people say this?
28:08
Well here's an example that illustrates it.
28:12
If I take a polynomial in two variables and I--
28:19
for example, this is a polynomial of two variables--
28:23
and my x matrix comes from sampling
28:27
that polynomial integers--
28:30
for example, this matrix--
28:38
then that matrix happens to be of low rank--
28:41
28:44
mathematically of low rank, with epsilon equals zero.
28:50
Why is that?
28:50
Well if I write down x in terms of matrices,
28:54
you could easily see it.
28:56
So this is made up of a matrix of all ones
29:00
plus a matrix of j-- so that's 1, 2, up to n, 1, 2, up to n,
29:11
because every entry of that matrix
29:12
just depends on the row index.
29:15
And then this guy depends on both j and k.
29:18
So this is a multiplication table, right?
29:21
So this is n, 2, 4, up to 2n, n, 2n, n squared.
29:31
OK.
29:34
Clearly, the matrix of all ones is a rank one matrix.
29:38
29:42
The same with this guy.
29:43
The column space is just of dimension one.
29:47
And the last guy also happens to be of rank one
29:51
because I can write this matrix in rank one form, which
29:58
is a column vector times a row vector.
30:03
OK.
30:04
So this matrix x is of rank three.
30:07
30:13
I guess at lowest rank three is what I've actually shown.
30:16
OK.
30:17
30:20
Now of course this hasn't got to numerical low rank yet,
30:23
so let's get ourselves there.
30:24
30:28
So Reade knew this, and he said to himself, OK,
30:32
well if I can approximate--
30:35
if x is actually coming from sampling a function,
30:38
and I approximate that function by polynomial,
30:41
then I'm going to get myself a low-rank approximation
30:45
and get a bound on the numerical rank.
30:48
So in general, if I give you a polynomial of two variables,
30:56
which can be written down--
30:58
it's degree n in both x and y.
31:04
Let's just keep these indexes away from the matrix index.
31:07
31:10
I give you this such polynomial, and I go away
31:14
and I sample it and make a matrix X, then X,
31:22
by looking at each term individually like I
31:24
did there, will have low rank mathematically,
31:30
with epsilon equals zero.
31:31
This will have, at most, m squared rank,
31:35
and if m is 3 or 4 or 10, it possibly
31:39
could be low because this X could be a large matrix.
31:43
OK.
31:44
So what Reade did for the Hilbert matrix was said,
31:47
OK, well look at that guy.
31:49
That guy looks like it's sampling a function.
31:52
It looks like it's sampling the function
31:53
1 over x plus y minus 1.
31:57
So he said to himself, well, that x,
32:02
if I look at the Hilbert matrix, then that
32:07
is sampling a function.
32:09
It happens to not be a polynomial.
32:13
It happens to be this function.
32:16
But that's OK because sampling polynomials, integers,
32:20
gives me low rank exactly.
32:22
Maybe sampling smooth functions, functions like this,
32:27
can be well approximated by polynomials
32:29
and therefore have low numerical rank.
32:32
And that's what he did in this case.
32:34
So he tried to find a p, a polynomial approximation to f.
32:42
In particular, he looked at exactly this kind
32:45
of approximation.
32:46
32:50
So he has some numbers here so that things
32:54
get dissolved later.
32:55
And he tried to find a p that did this kind of approximation.
33:01
So this approximates f.
33:03
33:08
And then he would develop a low-rank approximation to X
33:14
by sampling p.
33:16
So he would say, OK, well if I let y be a sampling of p, then
33:26
from the fact that f is a good approximation to p,
33:29
y is a good approximation to X. And so this has finite rank.
33:34
33:38
He wrote down that this must hold.
33:43
33:46
And the epsilon comes out here because these factors
33:49
were chosen just right.
33:51
The divide by n was chosen so that the epsilon came out just
33:54
there.
33:55
OK?
33:56
So, for many years, that was kind of the canonical reason
33:59
that people would give, that, well,
34:01
if the matrix X is sampled from a smooth function,
34:05
then we can approximate our function by a polynomial
34:10
and get polynomial rank approximations.
34:15
And therefore, the matrix X will be of low numerical rank.
34:18
34:22
There's an issue with this reasoning,
34:26
especially for the Hilbert matrix,
34:28
that it doesn't actually work that well.
34:31
So for example, if I take the 1,000 by 1,000 Hilbert matrix
34:38
and I look at its rank--
34:41
OK, well I've already told you this is full rank.
34:45
You'll get 1,000.
34:46
All the singular values are positive.
34:50
If I look at the numerical rank of this 1,000
34:55
by 1,000 Hilbert matrix and I compute it, I compute the SVD
35:00
and I look at how many are above epsilon where epsilon is 10
35:06
to the minus 15, so that means I can
35:10
approximate the 1,000 by 1,000 Hilbert matrix
35:13
by a rank 28 matrix and only give up
35:18
15-- there will be exact 15 digits, which is a huge amount.
35:24
So this is what we get in practice,
35:27
but Reade's argument here shows that the rank of this matrix,
35:42
the numerical rank, is at most.
35:45
35:49
So it doesn't do a very good job on the Hilbert matrix
35:53
for bounding the rank, right?
35:57
So Reade comes along, takes this function.
36:00
He tries to find a polynomial that does this, where
36:02
epsilon is 10 to the minus 15.
36:04
He finds that the number of terms
36:07
that he needs in this expression here is around 719,
36:13
and therefore, that's the rank that he gets.
36:16
The bound on the numerical rank.
36:19
The trouble is that 719 tells us that this is not
36:25
of low numerical rank, but we know
36:27
it is, so it's an unsatisfactory reason.
36:32
So there's been several people trying
36:36
to come up with more appropriate reasons that
36:39
explain the 28 here.
36:44
And so one reason that I've started to use
36:50
is another slightly different way
36:52
of looking at things, which is to say the world is Sylvester.
36:57
37:03
Now Sylvester, what does that mean?
37:10
What does the word "Sylvester" mean in this case?
37:13
It means that the matrices satisfy
37:14
a certain type of equation called the Sylvester equation,
37:20
and so the reason is really, many of these matrices
37:25
satisfy a Sylvester equation, and that takes the form--
37:30
37:36
for sum A, B, and C.
37:44
OK.
37:44
So X is your matrix of interest.
37:46
You want to show X is of numerical low rank.
37:50
And the task at hand is to find an A, B, and C so
37:54
that X satisfies that equation.
37:58
OK.
38:00
For example, the two matrices I've had on the board
38:05
satisfy a Sylvester equation--
38:09
a Sylvester matrix equation.
38:10
There is an A, a B, and a C for which they do this.
38:14
For example, remember the Hilbert matrix,
38:17
which we have there still, but I'll write it down again.
38:20
38:24
Has these entries.
38:26
So all we need to do is to try to figure out
38:28
an A, a B, and then a C so that we can make
38:32
it fit a Sylvester equation.
38:34
There's many different ways of doing this.
38:36
The one that I like is the following,
38:41
where if I put 1/2 here and 3/2 here,
38:45
all the way down to n minus 1/2, times this matrix--
38:51
38:53
so this is timesing the top of this matrix by 1/2
38:59
and then 3/2 and then 5/2.
39:02
So we're basically timesing each entry of this matrix
39:05
by j minus 1/2.
39:07
39:10
And then I do something on the right
39:12
here, which I'm allowed to do because I've got the B freedom,
39:14
and I choose this to be the same up to a minus sign.
39:19
39:23
Then when you think about this, what is it doing?
39:26
It's timing the jk entry--
39:30
this is-- by j minus 1/2.
39:33
That's what this is doing.
39:34
And what's this doing is timesing the jk entry
39:37
by k minus 1/2.
39:40
So this is, in total, timesing the jk entry
39:44
by j plus k minus 1/2 minus 1/2, which is minus 1,
39:49
so this is timesing the jk entry by j plus k minus 1.
39:54
So it knocks out the denominator.
39:57
And what we get from this equation is a bunch of ones.
40:02
40:11
So in this case, A and B are diagonal,
40:14
and C is the matrix of all ones.
40:17
OK?
40:17
We can also do this for Vandermonde.
40:20
So Vandermonde, you'll remember, looks like this.
40:25
40:30
And then over here, we have this guy,
40:35
the matrix that appears with polynomial interpolation.
40:40
OK.
40:41
So if I think about this, I could also
40:44
come up with an A, B, and C, and for example, here's
40:50
one that works.
40:52
I can stick the x's on the diagonal.
40:55
40:59
So if you imagine what that matrix on the left is doing,
41:03
it's timesing each column by the vector x.
41:08
OK?
41:08
So the first column of this matrix becomes x, the vector x.
41:12
The second becomes the vector x squared,
41:16
where squared is done entry-wise.
41:18
And then the third entry is now x cubed,
41:21
and when we get to the last, it's x to the n.
41:24
OK?
41:24
So that's like, multiply each column by the vector x.
41:30
So if I want to try to come up with a matrix--
41:32
so what's left is of low rank, is like of this form.
41:36
What I can do is shift the columns.
41:40
So I've noticed that this product here,
41:43
this diagonal matrix, has made the first column x.
41:46
So if I want to kill off that column,
41:48
I can take the second column and permute it to the first column.
41:52
I could take the third column and permute it
41:54
to the second, the last column and permute it
41:56
to the penultimate column here.
41:58
And that will actually kill off a lot
42:00
of what I've created in this matrix right here.
42:03
So let me write that down.
42:05
This is a circumshift matrix.
42:08
This does that permutation.
42:10
42:17
I've put a minus 1 there.
42:18
I could have put any number there.
42:19
It doesn't make any difference.
42:22
But this is the one that works out extremely nicely.
42:25
Now this zeros out lots of things because of the way
42:29
I've done the multiplication by x and the circumshift
42:32
of the columns.
42:33
And so the first column is zero because this first column is x,
42:39
this first column is x, so I've got x minus x.
42:43
This column was x squared minus x squared, so I got zero,
42:47
and I just keep going along until that last column.
42:51
That last column is a problem because the last column
42:54
of this guy is x to the n, whereas I
42:57
don't have x to the n in V, so there are some numbers here.
43:02
OK.
43:02
43:05
You'll notice that C in both cases
43:08
happens to be a low-rank matrix.
43:11
In these cases, it happens to be of rank one.
43:14
And so people were wondering, maybe it's
43:19
something to do with satisfying these kind of equations that
43:23
makes these matrices that appear in practice
43:27
numerically of low rank.
43:29
And after a lot of work in this area,
43:33
people have come up with a bound that
43:37
demonstrates that these kind of equations
43:42
are key to understanding numerical low rank.
43:46
So if X satisfies a Sylvester equation, like this, and A
44:04
is normal, B is normal--
44:06
I don't really want to concentrate
44:08
on those two conditions.
44:12
It's a little bit academic.
44:17
Then-- people have found a bound on the singular
44:21
values of any matrix that satisfies
44:24
this kind of expression, and they
44:27
found this following bound.
44:30
44:41
OK, so here, the rank of C is r.
44:49
So that goes there.
44:50
So in our cases, the two examples we have,
44:52
r is 1, so we can forget about r.
44:55
This nasty guy here is called the Zolotarev number.
45:02
45:07
E is a set that contains the eigenvalues of A,
45:12
and F is a set that contains the eigenvalues of B. OK.
45:23
Now it looks like we have gained absolutely nothing
45:26
by this bound, because I've just told you singular values
45:30
are bound by Zolotarev numbers.
45:32
That doesn't mean anything to anyone.
45:35
It means a little bit to me but not that much.
45:38
So the key to this bound--
45:42
the reason this is useful--
45:43
is that so many people have worked out what these Zolotarev
45:49
numbers actually mean.
45:52
OK?
45:52
So these are two key people that worked out
45:57
what this bound means.
45:59
And we have gained a lot because people
46:02
have been studying this number.
46:04
This is, like, a number that people
46:06
cared about from 1870 onwards to the present day,
46:11
and people have studied this number extremely well.
46:14
So we've gained something by turning it
46:17
into a more abstract problem that people have thought
46:21
about previously, and now we can go
46:23
to the literature on Zolotarev numbers,
46:26
whatever they are, and discover this whole literature of work
46:30
on this Zolotarev number.
46:32
And the key part--
46:34
I'll just tell you the key--
46:35
46:41
is that the sets E and F are separated.
46:44
46:53
So for example, in the Hilbert matrix, the eigenvalues of A
46:57
can be read off the diagonal.
46:59
47:05
What are they?
47:06
They are between minus 1/2 and n minus 1/2.
47:13
And the eigenvalues of B lie in the set minus 1/2
47:20
minus n plus 1/2.
47:23
And the key reason why the Hilbert matrix
47:27
is of low numerical rank is the fact
47:30
that these two sets are separated,
47:33
and that makes this Zolotarev number gets small extremely
47:36
quickly with k.
47:38
Now you might wonder why there is a question
47:41
mark on Penzl's name.
47:44
There is an unofficial curse that's
47:49
been going on for a while.
47:51
Both these men died while working on the Zolotarev
47:54
problem.
47:55
They both died at the age of 31.
48:00
One died by being hit by a train, Zolotarev.
48:03
It's unclear whether he was suicidal or it was accidental.
48:08
Penzl died at the age of 31 in the Canadian mountains
48:13
by an avalanche.
48:16
I am currently not yet 31 but going to be 31 very soon,
48:21
and I'm scared that I may join this list.
48:23
48:27
OK.
48:29
But for the Hilbert matrix, what you get from this analysis,
48:32
based on these two peoples' work,
48:36
is a bound on the numerical rank.
48:39
And the rank that you get is, let's say,
48:43
a world record bound.
48:45
For the Hilbert matrix is 34, which is not quite 28, not yet,
48:56
but it's far more descriptive of 28 than 719.
49:02
And so this technique of bounding singular values
49:07
by using these Zolotarev numbers is starting to gain popularity
49:11
because we can finally answer to ourselves why there are so
49:15
many low-rank matrices that appear in computational math.
49:20
And it's all based on two 31-year-olds that died.
49:27
And so if you ever wonder when you're
49:30
doing computational science when a low rank appears
49:33
and the smoothness argument does not work for you,
49:36
you might like to think about Zolotarev and the curse.
49:40
OK, thank you very much.
49:42
[APPLAUSE]
49:44
49:47
GILBERT STRANG: Thank you [INAUDIBLE] Excellent.
49:49
ALEX TOWNSEND: How does it work now?
49:51
GILBERT STRANG: We're good.
49:53
Yeah.
49:54
ALEX TOWNSEND: I'm happy to take questions
49:55
if we have a minute, if you have any questions.
49:57
GILBERT STRANG: How near of 31 are you?
49:59
50:02
ALEX TOWNSEND: [INAUDIBLE] I get a spotlight.
50:03
I'm 31 in December.
50:05
GILBERT STRANG: Wow.
50:06
OK.
50:07
ALEX TOWNSEND: So they died at the age of 31, so you know,
50:10
next year is the scary year for me.
50:14
So I'm not driving anywhere.
50:16
I'm not leaving my house until I become 32.
50:20
GILBERT STRANG: Well, thank you [INAUDIBLE]
50:22
ALEX TOWNSEND: OK, thanks.
50:24
[APPLAUSE]
50:26
English - CC (English)
###