Fast Estimation of the Spectral Norm of Random Matrices

In the section "Approximate Estimation" of Higher-Order MuP: A Simpler yet Smarter Spectral-Condition Scaling], we borrowed a conclusion in advance without proof: "A random matrix of size $n\times m$ with entries drawn from the standard normal distribution has a spectral norm of roughly $\sqrt{n}+\sqrt{m}$."

In this post, we go back and fill in that gap, presenting a fast way to estimate the spectral norm of a random matrix.

Random Matrix Theory

Consider a random matrix $\boldsymbol{W}\in\mathbb{R}^{n\times m}$, each of whose entries is drawn independently and identically from the standard normal distribution $\mathcal{N}(0,1)$. We want to estimate the spectral norm of $\boldsymbol{W}$, i.e., its largest singular value, and we take $\mathbb{E}[\Vert\boldsymbol{W}\Vert_2]$ as our final estimate. more

First, it should be pointed out that the study of random matrix behavior has already grown into its own specialized branch. As far as estimating the spectral norm of normal random matrices is concerned, the relevant keywords include the "Marchenko–Pastur distribution]," the "Bai-Yin law," "Gordon's theorem," and so on, all of which come with detailed estimation results for the singular value distribution. However, we are not going to introduce any of that material here; instead, we will present a quick way to arrive at $\sqrt{n}+\sqrt{m}$.

The method below is a simplification, by the author, of Section 5.3.1 of Introduction to the Non-Asymptotic Analysis of Random Matrices]. It's really more of a "popular-science" explanation than a rigorous derivation—its purpose is to let readers quickly get a feel for where $\sqrt{n}+\sqrt{m}$ comes from. A rigorous proof would require filling in a great many tedious, unremarkable details, which we will skip here.

Estimating the Spectral Norm

Our starting point is the identity

\begin{equation}\Vert\boldsymbol{W}\Vert_2 = \max_{\Vert \boldsymbol{u}\Vert=1, \Vert \boldsymbol{v}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v}\label{eq:w-uv}\end{equation}

where $\boldsymbol{u}\in\mathbb{R}^n,\boldsymbol{v}\in\mathbb{R}^m$. Direct computation of $\Vert\boldsymbol{W}\Vert_2$ is generally not easy, so it's natural to look for some kind of approximation. Let's consider the following two "half-finished products":

\begin{equation}\max_{\Vert \boldsymbol{u}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v}\qquad\qquad \max_{\Vert \boldsymbol{v}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v}\end{equation}

Intuitively, compared with equation $\eqref{eq:w-uv}$, each of the two expressions above only accomplishes half the job. Let's make a bold assumption: their results are each also "half" of the full result. Adding them together then gives us what we take to be a good approximation of the final result:

\begin{equation}\Vert\boldsymbol{W}\Vert_2 = \max_{\Vert \boldsymbol{u}\Vert=1, \Vert \boldsymbol{v}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v} \approx \max_{\Vert \boldsymbol{u}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v} + \max_{\Vert \boldsymbol{v}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v} = \Vert \boldsymbol{W}\boldsymbol{v}\Vert + \Vert \boldsymbol{u}^{\top}\boldsymbol{W}\Vert \label{eq:core-approx}\end{equation}

In other words, we separately sample $\boldsymbol{u},\boldsymbol{v}$ from the unit hypersphere in $\mathbb{R}^n,\mathbb{R}^m$, and then use the above formula to obtain an approximation of the spectral norm of $\boldsymbol{W}$.

Computing the Expectation

With this approximation in hand, we can compute

\begin{equation}\mathbb{E}[\Vert\boldsymbol{W}\Vert_2]\approx\mathbb{E}[\Vert \boldsymbol{W}\boldsymbol{v}\Vert] + \mathbb{E}[\Vert \boldsymbol{u}^{\top}\boldsymbol{W}\Vert] \approx \sqrt{\mathbb{E}[\Vert \boldsymbol{W}\boldsymbol{v}\Vert^2]} + \sqrt{\mathbb{E}[\Vert \boldsymbol{u}^{\top}\boldsymbol{W}\Vert^2]}\end{equation}

where

\begin{equation}\mathbb{E}[\Vert \boldsymbol{W}\boldsymbol{v}\Vert^2] = \mathbb{E}[ \boldsymbol{v}^{\top}\boldsymbol{W}^{\top}\boldsymbol{W}\boldsymbol{v}] = \boldsymbol{v}^{\top}\mathbb{E}[\boldsymbol{W}^{\top}\boldsymbol{W}]\boldsymbol{v} = \boldsymbol{v}^{\top}(n\boldsymbol{I}_m)\boldsymbol{v} = n\end{equation}

By the same token, $\mathbb{E}[\Vert \boldsymbol{u}^{\top}\boldsymbol{W}\Vert^2]=m$, and so

\begin{equation}\mathbb{E}[\Vert\boldsymbol{W}\Vert_2]\approx\sqrt{n} + \sqrt{m}\end{equation}

This result turns out to be remarkably accurate (which can be verified through simulation). Specifically, if $n=ka,m=kb$, with $a,b$ both constants, then

\begin{equation}\lim_{k\to\infty} \frac{\Vert\boldsymbol{W}\Vert_2}{\sqrt{n} + \sqrt{m}} = 1,\qquad \boldsymbol{W}\sim\mathcal{N}(0,1)\end{equation}

The reason it's so accurate is that we cheated a bit—the crucial equation $\eqref{eq:core-approx}$ was essentially reverse-engineered from the already-known correct answer. Aside from equation $\eqref{eq:core-approx}$, the only other conditions we used are $\mathbb{E}[\boldsymbol{W}^{\top}\boldsymbol{W}]=n\boldsymbol{I}_m$ and $\mathbb{E}[\boldsymbol{W}\boldsymbol{W}^{\top}]=m\boldsymbol{I}_n$, so we can regard this approximation as holding for any distribution with zero mean and unit variance.

Smallest Singular Value

The spectral norm is the largest singular value; in fact, we can use the same line of reasoning to estimate the smallest singular value as well. Of course, "smallest" here needs to be defined more carefully. Let $n\geq m$; the smallest singular value we refer to here is the $m$-th largest singular value of $\boldsymbol{W}$ (counting from the top), and it is equal to

\begin{equation}\sigma_{\min}(\boldsymbol{W}) = \min_{\Vert \boldsymbol{v}\Vert=1}\max_{\Vert \boldsymbol{u}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v}\end{equation}

Note that in this expression, the positions and objects of $\min$ and $\max$ cannot be swapped. Similarly, we consider the sum of two "half-finished" approximations, one over $\min$ and one over $\max$, as an approximation to it:

\begin{equation}\sigma_{\min}(\boldsymbol{W}) = \min_{\Vert \boldsymbol{v}\Vert=1}\max_{\Vert \boldsymbol{u}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v}\approx \min_{\Vert \boldsymbol{v}\Vert=1}\boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v} + \max_{\Vert \boldsymbol{u}\Vert=1} \boldsymbol{u}^{\top}\boldsymbol{W} \boldsymbol{v} = -\Vert \boldsymbol{u}^{\top}\boldsymbol{W}\Vert + \Vert \boldsymbol{W}\boldsymbol{v}\Vert \end{equation}

Everything that follows proceeds exactly as before, and the result is $\mathbb{E}[\sigma_{\min}(\boldsymbol{W})]\approx\sqrt{n}-\sqrt{m}$.

A Final Note

This post has offered a quick way of thinking about estimating the spectral norm of a random matrix—or, to put it more precisely, a popular-science, heuristic explanation rather than a rigorous, exact derivation. It could in principle be made rigorous, but doing so would require filling in a great many theoretical details, all of which we have skipped over here.

English translation of a post from 科学空间 | Scientific Spaces by 苏剑林. Original: https://kexue.fm/archives/11335
Translated automatically with claude-sonnet-5; all equations are reproduced verbatim from the source. Copyright remains with the original author.