Notebook 16 – Math 2121, Fall 2020
In this notebook we explore more random matrices and their eigenvalues.
Last time we saw that the eigenvalues of a large orthogonal matrix chosen uniformly at random are, with high probability, evenly spaced complex numbers in the unit circle
Thus, the process of choosing a random eigenvalue of a very large, uniformly random orthogonal matrix is indistiguishable from choosing an element of
Today we will examine a similar phenomenon for random symmetric matrices.
xxxxxxxxxx
using LinearAlgebra
xxxxxxxxxx
using Plots
xxxxxxxxxx
using Distributions
First let's discuss a few more details about the standard normal distribution.
Here is an intuitive definition of this probability distribution:
Suppose some students are taking an exam with
different True/False questions.On each question, a student has equal chances of answering correctly or incorrectly.
Students receive
points for a correct answer, and points otherwise.The final score on the exam is the sum of these
points divided by .
Since we divide by
Definition: The standard normal distribution is the probability distribution that describes the random scores on this exam when the number of questions is very large.
10000
xxxxxxxxxx
questions = 10000
50000
xxxxxxxxxx
students = 50000
To observe this empirically, we can simulate the scores for 50000 students on a 10000 question exam.
50000×1 Array{Float64,2}:
-1.98
0.98
0.08
-0.1
-0.36
-0.34
0.26
⋮
-0.06
-2.22
0.1
-1.08
0.96
-0.24
xxxxxxxxxx
begin
points = rand([-1, 1], students, questions)
scores = sum(points, dims=2) / sqrt(questions)
end
Now we compare the histogram of scores to a histogram of 50000 samples of the standard normal distribution (computed using a Julia library function).
These plots should look almost the same if our parameters are large enough.
xxxxxxxxxx
begin
h1 = histogram(
scores,
legend=false, normalize=:pdf, xlim=(-4,4),
title="Histogram of $(students) scores on a $(questions) question exam")
h2 = histogram(
rand(Normal(0, 1), students, 1),
normalize=:pdf, legend=false, xlim=(-4,4),
title="Histogram of $(students) standard normal samples")
h3 = plot(
x -> exp(-x^2/2) / sqrt(2 * pi),
legend=false, xlim=(-4,4),
title="Standard normal density function")
plot(h1, h2, h3, layout=(3,1), size=(680, 560))
end
All matrices today have real entries.
A square matrix
3×3 Array{Int64,2}:
1 3 5
3 0 6
5 6 2
xxxxxxxxxx
[1 3 5; 3 0 6; 5 6 2]
All eigenvalues of a symmetric matrix are real numbers. For example:
-5.3258
-2.2393
10.5651
xxxxxxxxxx
eigen([1 3 5; 3 0 6; 5 6 2]).values
Compare with the eigenvalues of a random (non-symmetric) square matrix:
-0.716167+0.0im
-0.443853-0.600181im
-0.443853+0.600181im
-0.157667+0.0im
0.108804-0.339849im
0.108804+0.339849im
0.543573+0.0im
0.655452-0.530417im
0.655452+0.530417im
5.17306+0.0im
xxxxxxxxxx
eigen(rand(10, 10)).values
Here is a proof of this fact in general.
Assume
Suppose
The product
On the other hand we have both
and, because
These expressions must be equal, and the number
The set of symmetric
Here is a different way of sampling a random symmetric matrix.
First, we generate a random matrix with entries from the standard normal distribution.
5000
xxxxxxxxxx
n = 5000
xxxxxxxxxx
A = rand(Normal(0, 1), n, n);
Now to get a random symmetric matrix, we can just form the sum
To make later plots nicer, we also divide by the normalizing constant
5000×5000 Array{Float64,2}:
-0.00284525 0.0031909 -0.0115228 … 0.00735712 -0.00644351 -0.0120708
0.0031909 -0.0210216 0.00159086 -0.0109171 -0.00906114 0.0124589
-0.0115228 0.00159086 0.0177357 -0.00719807 -0.00122543 -0.00144259
-0.0088009 0.00934093 -0.00193897 0.017209 -0.00196763 0.0053541
0.010613 -0.000949646 -0.00718409 0.0034996 0.00706539 -0.0124317
-0.00219485 -0.00013567 -0.0111961 … 0.00365358 -0.000357063 -0.00754908
-0.00102919 0.00232089 -0.00990527 -0.0034211 -0.00695158 0.00374101
⋮ ⋱
0.010426 0.0027998 0.00216689 -0.0119353 -0.00485233 0.00267533
0.00492574 -0.00829167 0.0041719 … 0.0131883 -0.00866356 -0.00655834
-0.0016651 -0.000921308 0.00991947 -0.0120893 0.00111046 0.0184311
0.00735712 -0.0109171 -0.00719807 -0.00223987 0.00224638 -0.0108681
-0.00644351 -0.00906114 -0.00122543 0.00224638 -0.00971652 0.0104294
-0.0120708 0.0124589 -0.00144259 -0.0108681 0.0104294 -0.00556288
xxxxxxxxxx
S = (A + transpose(A)) / sqrt(8 * n)
Here are the (real) eigenvalues of this matrix, in sorted order:
-1.00232
-0.997142
-0.994154
-0.990228
-0.989757
-0.986384
-0.985493
-0.98449
-0.981978
-0.979637
-0.978236
-0.977056
-0.975703
-0.974927
-0.973047
-0.970075
-0.96974
-0.969321
-0.967861
-0.965469
0.979101
0.979994
0.980792
0.981251
0.985273
0.986169
0.987367
0.989432
0.991636
0.995266
The histogram of these eigenvalues, when
This histogram is the unit semicircle with height rescaled by factor
0.6366197723675814
xxxxxxxxxx
2/pi
xxxxxxxxxx
begin
bins = 100
xlim = (-1.8, 1.8)
q1 = histogram(
eig,
legend=false, normalize=:pdf, xlim=xlim, bins=bins,
title="Histogram of eigenvalues of random symmetric matrix")
q2 = plot(
x -> 2 * sqrt(1 - x^2) / pi,
xlim=xlim, legend=false,
title="Unit semicircle (height rescaled so area under curve is 1)")
plot(q1, q2, layout=(2,1), size=(680, 500))
end
Thus, the eigenvalues of a large random symmetric matrix with entries from the standard normal distribution obey yet another probability law that is neither the uniform distribution nor the normal distribution.
Instead, these random eigenvalues follow what is called the Wigner semicircle distribution.
We can again observe this empirically using Julia's built-in functions for the semicircle distribution:
xxxxxxxxxx
begin
q3 = histogram(
rand(Semicircle(1), 1000000, 1),
normalize=:pdf, legend=false, xlim=xlim, bins=bins,
title="Histogram of samples from Wigner semicircle distribution")
plot(q1, q3, layout=(2, 1))
end