The Remarkable Johnson–Lindenstrauss Lemma: Theory
Today we're going to study the Johnson–Lindenstrauss lemma; since the name is a bit long, we'll abbreviate it as the "JL lemma" from here on.
Personally, I think the JL lemma is one of those magical results that every computer science student should know. It's a famous result about dimensionality reduction, and also one of the classic examples of the many counterintuitive "curse of dimensionality" phenomena in high-dimensional spaces. One could say that the JL lemma is the theoretical foundation of various dimensionality reduction and hashing techniques in machine learning. Moreover, in modern machine learning, the JL lemma also provides important theoretical support for understanding and tuning model dimension-related hyperparameters.
Logarithmic Dimensions
The JL lemma can be stated in a very colloquial way as follows:
Colloquial version of the JL lemma: To fit $N$ vectors, you only need a $\mathcal{O}(\log N)$-dimensional space.
Specifically, the JL lemma says that regardless of how many dimensions these $N$ vectors originally have, we can reduce them down to $\mathcal{O}(\log N)$ dimensions while keeping the error in their relative distances within a controlled range. As you can imagine, this is an extremely strong, extremely counterintuitive, and extremely practical result. For example, suppose we want to do vector retrieval: the original vector dimension might be very large, making a full-scale retrieval very costly. The JL lemma tells us that we can transform these vectors into $\mathcal{O}(\log N)$ dimensions while the retrieval performance stays approximately the same — this really is a case of "a pie falling from the sky"!
Readers might wonder: for such a strong result, shouldn't the corresponding dimension reduction method be extremely complicated? The answer is quite the opposite — the reduction process only requires a random linear projection! Some have even remarked that the JL lemma is "a result that's easier to prove than to understand," meaning that proving it mathematically isn't actually all that hard, but intuitively grasping this counterintuitive result is not so easy.
Coincidentally, we've actually introduced two counterintuitive results before: in the article The angle distribution between two random vectors in n-dimensional space, we discussed how "any two random vectors in high-dimensional space are almost always perpendicular," which is clearly quite different from what happens in 2D or 3D space. In the article Understanding parameter initialization strategies for models from a geometric perspective, this result was further upgraded to "a $n\times n$ matrix sampled from $\mathcal{N}(0,1/n)$ is almost an orthogonal matrix," which conflicts sharply with our usual understanding that "orthogonality is a very strict condition (requiring the transpose to equal the inverse)."
But in fact, both of these results are not only correct, but also directly related to the JL lemma. One could say the JL lemma can be viewed as a refinement and application of these two results. So, we need to first describe these two results in more quantitative language — for "almost perpendicular," what exactly is the probability of being perpendicular, and for "approximately orthogonal," what exactly is the magnitude of the error?
Probability Inequalities
To this end, we need some probability theory, the most important of which is "Markov's inequality":
Markov's inequality: If $x$ is a nonnegative random variable and $a > 0$, then
\begin{equation}P(x\geq a)\leq \frac{\mathbb{E}[x]}{a}\end{equation}
Note that this inequality doesn't impose any other particular restriction on the distribution of $x$, other than requiring that the random variable takes values in the nonnegative range (or equivalently, that the probability of $x$ being negative is always 0). The proof is actually very simple:
\begin{equation}\mathbb{E}[x]=\int_0^{\infty}x p(x) \geq \int_a^{\infty}x p(x) \geq \int_a^{\infty} a p(x) = a P(x\geq a)\end{equation}
Markov's inequality requires the random variable to be nonnegative, but the random variables we usually deal with aren't necessarily nonnegative, so we typically need some transformation before we can use it. For instance, $x - \mathbb{E}[x]$ is not nonnegative, but $|x - \mathbb{E}[x]|$ is, so applying Markov's inequality gives us:
\begin{equation}P(|x - \mathbb{E}[x]|\geq a) = P((x - \mathbb{E}[x])^2\geq a^2) \leq \frac{\mathbb{E}[(x - \mathbb{E}[x])^2]}{a^2}=\frac{\mathbb{V}ar[x]}{a^2}\end{equation}
This is "Chebyshev's inequality."
Another classic technique is called the "Cramér–Chernoff method," which is also the main method we'll rely on later. It turns the random variable nonnegative via an exponential function: for any $\lambda > 0$, we have
\begin{equation}x \geq a \quad\Leftrightarrow\quad \lambda x \geq \lambda a \quad\Leftrightarrow\quad e^{\lambda x} \geq e^{\lambda a}\end{equation}
So applying Markov's inequality gives
\begin{equation}P(x \geq a) = P(e^{\lambda x} \geq e^{\lambda a})\leq e^{-\lambda a}\mathbb{E}[e^{\lambda x}]\end{equation}
The leftmost side is independent of $\lambda$, but the rightmost side contains a $\lambda$, and this inequality holds for any $\lambda > 0$. So in theory, we can find the $\lambda$ that minimizes the right-hand side, to obtain the tightest bound:
\begin{equation}P(x \geq a) \leq \min_{\lambda > 0} e^{-\lambda a}\mathbb{E}[e^{\lambda x}]\end{equation}
A Lemma Behind the Lemma
Now we can introduce the following result, which is a lemma underlying the JL lemma — indeed, one could say it's the theoretical foundation of everything in this article:
Unit-norm lemma: Let $u\in\mathbb{R}^{n}$ be a vector sampled i.i.d. from $\mathcal{N}(0,1/n)$, and let $\varepsilon \in (0, 1)$ be a given constant. Then we have
\begin{equation}P(|\Vert u\Vert^2 - 1| \geq \varepsilon) \leq 2\exp\left(-\frac{\varepsilon^2 n}{8}\right)\end{equation}
This lemma tells us that when $n$ is sufficiently large, the probability that the norm of $u$ deviates significantly from 1 is very small (once $\varepsilon$ is fixed, it decays to 0 exponentially in $n$). So a $n$-dimensional vector sampled from $\mathcal{N}(0,1/n)$ will be very close to a unit vector.
The proof relies precisely on the "Cramér–Chernoff method": first, $|\Vert u\Vert^2 - 1| \geq \varepsilon$ implies either $\Vert u\Vert^2 - 1 \geq \varepsilon$ or $1 - \Vert u\Vert^2\geq \varepsilon$, and we need to handle these two cases separately. Without loss of generality, let's first derive the probability for $\Vert u\Vert^2 - 1 \geq \varepsilon$. By the Cramér–Chernoff method, we have
\begin{equation}P(\Vert u\Vert^2 - 1 \geq \varepsilon) \leq \min_{\lambda > 0} e^{-\lambda \varepsilon}\mathbb{E}\big[e^{\lambda (\Vert u\Vert^2 - 1)}\big] = \min_{\lambda > 0} e^{-\lambda (\varepsilon + 1)}\mathbb{E}\big[e^{\lambda \Vert u\Vert^2}\big]\end{equation}
Writing $u$ in component form as $(u_1, u_2, \cdots, u_n)$, where each component is independent and identically distributed as $\mathcal{N}(0,1/n)$, we get
\begin{equation}\mathbb{E}\big[e^{\lambda \Vert u\Vert^2}\big] = \mathbb{E}\big[e^{\lambda\sum\limits_i u_i^2}\big] = \mathbb{E}\big[\prod_i e^{\lambda u_i^2}\big]=\prod_i \mathbb{E}\big[ e^{\lambda u_i^2}\big]\end{equation}
and since $\mathbb{E}\big[ e^{\lambda u_i^2}\big]=\int_{-\infty}^{\infty} \frac{1}{\sqrt{2\pi}}e^{-u_i^2/2}e^{\lambda u_i^2/n} du_i=\sqrt{n/(n-2\lambda)}$, we have
\begin{equation}P(\Vert u\Vert^2 - 1 \geq \varepsilon) \leq \min_{\lambda > 0} e^{-\lambda (\varepsilon + 1)}\left(\frac{n}{n-2\lambda}\right)^{n/2}\end{equation}
The minimum of the right-hand side is achieved at $\lambda = \frac{n\varepsilon}{2(1+\varepsilon)}$; we leave the derivation to the reader. Substituting back in, we get
\begin{equation}P(\Vert u\Vert^2 - 1 \geq \varepsilon) \leq e^{n(\log(1+\varepsilon) - \varepsilon)/2}\leq e^{-n\varepsilon^2/8}\end{equation}
where the proof of $\log(1+\varepsilon) - \varepsilon \leq -\varepsilon^2/4$ is also left to the reader. Similarly, we can derive the probability for $1 - \Vert u\Vert^2\geq \varepsilon$, obtaining:
\begin{equation}P(1 - \Vert u\Vert^2\geq \varepsilon) \leq e^{n(\log(1-\varepsilon) + \varepsilon)/2}\leq e^{-n\varepsilon^2/8}\end{equation}
where one can show $\log(1-\varepsilon) + \varepsilon \leq \log(1+\varepsilon) - \varepsilon$, so the above follows the same inequality relation as $\log(1+\varepsilon) - \varepsilon$. Now adding the two together, we get $P(|\Vert u\Vert^2 - 1| \geq \varepsilon)\leq 2e^{-n\varepsilon^2/8}$. This completes the proof.
Starting from the "unit-norm lemma," we can prove the "orthogonality lemma":
Orthogonality lemma: Let $u,v\in\mathbb{R}^{n}$ be two vectors sampled i.i.d. from $\mathcal{N}(0,1/n)$, and let $\varepsilon \in (0, 1)$ be a given constant. Then we have
\begin{equation}P(|\langle u, v\rangle| \geq \varepsilon) \leq 4\exp\left(-\frac{\varepsilon^2 n}{8}\right)\end{equation}
This lemma tells us that when $n$ is sufficiently large, the probability that the inner product of $u,v$ deviates significantly from 0 is very small (once $\varepsilon$ is fixed, it decays to 0 exponentially in $n$). So two $n$-dimensional vectors sampled from $\mathcal{N}(0,1/n)$ will be very nearly orthogonal. Combined with the "unit-norm lemma," we then arrive at the conclusion that "a matrix sampled from $\mathcal{N}(0,1/n)$ is almost an orthogonal matrix ($n\times n$)."
With the "unit-norm lemma" already established, this proof isn't too hard. We know that if $u,v\sim \mathcal{N}(0,1/n)$, then $\frac{u\pm v}{\sqrt{2}}\sim \mathcal{N}(0,1/n)$, so following the proof of the "unit-norm lemma," we have
\begin{equation}P\left(\left\Vert \frac{u+v}{\sqrt{2}}\right\Vert^2 - 1 \geq \varepsilon\right) \leq e^{-n\varepsilon^2/8},\quad P\left(1-\left\Vert \frac{u-v}{\sqrt{2}}\right\Vert^2 \geq \varepsilon\right) \leq e^{-n\varepsilon^2/8}\end{equation}
Note that adding equations $\left\Vert \frac{u+v}{\sqrt{2}}\right\Vert^2 - 1 \geq \varepsilon$ and $1-\left\Vert \frac{u-v}{\sqrt{2}}\right\Vert^2 \geq \varepsilon$ together gives $\langle u, v\rangle\geq \varepsilon$, so
\begin{equation}P(\langle u, v\rangle\geq \varepsilon) \leq P\left(\left\Vert \frac{u+v}{\sqrt{2}}\right\Vert^2 - 1 \geq \varepsilon\right) + P\left(1-\left\Vert \frac{u-v}{\sqrt{2}}\right\Vert^2 \geq \varepsilon\right) \leq 2e^{-n\varepsilon^2/8}\end{equation}
Similarly, we can prove $P(-\langle u, v\rangle\geq \varepsilon) \leq 2e^{-n\varepsilon^2/8}$, and combining the two gives us the "orthogonality lemma."
The Proof Process
Now we can go ahead and prove the JL lemma. Below is its mathematical statement:
Mathematical version of the JL lemma: Given $N$ vectors $v_1,v_2,\cdots,v_N\in\mathbb{R}^m$ and $n > \frac{24\log N}{\varepsilon^2}$, and a random matrix $A\in\mathbb{R}^{n\times m}$ sampled i.i.d. from $\mathcal{N}(0,1/n)$, with $\varepsilon \in (0, 1)$ a given constant, then with probability at least $\frac{N-1}{N}$, we have for all $i\neq j$:
\begin{equation}(1-\varepsilon)\Vert v_i - v_j\Vert^2 \leq \Vert Av_i - A v_j\Vert^2 \leq (1+\varepsilon)\Vert v_i - v_j\Vert^2\label{eq:bound}\end{equation}
The lemma tells us that regardless of the original vector dimension $m$, we only need a dimension of $n > \frac{24\log N}{\varepsilon^2}$ in order to accommodate $N$ vectors, such that the deviation in their relative distances doesn't exceed $\varepsilon$. Moreover, the JL lemma also gives us the dimension-reduction method: we simply need to randomly sample a $n\times m$ matrix $A$ from $\mathcal{N}(0,1/n)$, and then apply the transformation $v\to Av$, which achieves the goal with probability $\frac{N-1}{N}$. Truly simple and practical!
The proof is also a direct application of the "unit-norm lemma." First, if $u\in\mathbb{R}^m$ is a given unit vector, and $A\in\mathbb{R}^{n\times m}$ is sampled i.i.d. from $\mathcal{N}(0,1/n)$, then each component of $Au$ is independently distributed as $\mathcal{N}(0,1/n)$. This isn't hard to prove either: by definition, each component is $(Au)_i = \sum\limits_j A_{i,j}u_j$. Since $A_{i,j}$ are mutually independent, $(Au)_i$ are clearly mutually independent as well, and since $A_{i,j}\sim\mathcal{N}(0,1/n)$, the sum of normal random variables is still normally distributed, so $(Au)_i$ follows a normal distribution with mean $\sum\limits_j u_j\times 0=0$ and variance $\sum\limits_j u_j^2\times \frac{1}{n} = \frac{1}{n}$.
So, in plain terms, $Au$ is equivalent to a $n$-dimensional vector sampled i.i.d. from $\mathcal{N}(0,1/n)$. Now substituting in $u=\frac{v_i - v_j}{\Vert v_i - v_j\Vert}$ and using the "unit-norm lemma," we get
\begin{equation}P\left(\left|\left\Vert \frac{A(v_i - v_j)}{\Vert v_i - v_j\Vert}\right\Vert^2 - 1\right| \geq \varepsilon\right) \leq 2\exp\left(-\frac{\varepsilon^2 n}{8}\right)\end{equation}
This result holds for any $i\neq j$, so taking a union bound over all combinations of $i\neq j$, we find that the probability that at least one term $\geq \varepsilon$ fails is no more than
\begin{equation}P\left(\exists (i,j):\,\left|\left\Vert \frac{A(v_i - v_j)}{\Vert v_i - v_j\Vert}\right\Vert^2 - 1\right| \geq \varepsilon\right) \leq 2 {N\choose 2} \exp\left(-\frac{\varepsilon^2 n}{8}\right)\end{equation}
Or conversely, for any $i\neq j$, the probability that $\left|\left\Vert \frac{A(v_i - v_j)}{\Vert v_i - v_j\Vert}\right\Vert^2 - 1\right| \leq \varepsilon$ (equivalently $\eqref{eq:bound}$) holds simultaneously for all of them is at least
\begin{equation}1 - 2 {N\choose 2} \exp\left(-\frac{\varepsilon^2 n}{8}\right) = 1 - N(N-1)\exp\left(-\frac{\varepsilon^2 n}{8}\right)\end{equation}
Substituting in $n > \frac{24\log N}{\varepsilon^2}$, we get
\begin{equation}1 - N(N-1)\exp\left(-\frac{\varepsilon^2 n}{8}\right)\geq 1 - N(N-1)N^{-3}\geq 1-N^{-1}\end{equation}
This completes the proof.
The JL lemma above preserves Euclidean distance approximately. In many cases, however, retrieval relies on inner products (e.g., cosine similarity) rather than Euclidean distance. For this case, we have:
Inner-product version of the JL lemma: Given $N$ unit vectors $v_1,v_2,\cdots,v_N\in\mathbb{R}^m$ and $n > \frac{24\log N}{\varepsilon^2}$, and a random matrix $A\in\mathbb{R}^{n\times m}$ sampled i.i.d. from $\mathcal{N}(0,1/n)$, with $\varepsilon \in (0, 1)$ a given constant, then with probability at least $\frac{N-2}{N}$, we have for all $i\neq j$:
\begin{equation}\left|\langle Av_i, Av_j\rangle - \langle v_i, v_j\rangle\right|\leq\varepsilon\end{equation}
The proof is simple; we just mimic the proof of the "orthogonality lemma." According to the proof of the JL lemma, under the same conditions, with probability at least $\frac{N-2}{N}$, the following hold simultaneously for any $i\neq j$:
\begin{equation}\begin{aligned} (1-\varepsilon)\Vert v_i - v_j\Vert^2 \leq \Vert Av_i - A v_j\Vert^2 \leq (1+\varepsilon)\Vert v_i - v_j\Vert^2 \\ (1-\varepsilon)\Vert v_i + v_j\Vert^2 \leq \Vert Av_i + A v_j\Vert^2 \leq (1+\varepsilon)\Vert v_i + v_j\Vert^2 \end{aligned}\end{equation}
Multiplying the first by $-1$ gives $-(1+\varepsilon)\Vert v_i - v_j\Vert^2 \leq -\Vert Av_i - A v_j\Vert^2 \leq -(1-\varepsilon)\Vert v_i - v_j\Vert^2$, and adding it to the second gives
\begin{equation}4\langle v_i, v_j\rangle-2\varepsilon(\Vert v_i\Vert^2 + \Vert v_j\Vert)\leq 4\langle Av_i, Av_j\rangle \leq 4\langle v_i, v_j\rangle + 2\varepsilon(\Vert v_i\Vert^2 + \Vert v_j\Vert)\end{equation}
Noting that $v_i,v_j$ is a unit vector, the above is equivalent to $\left|\langle Av_i, Av_j\rangle - \langle v_i, v_j\rangle\right|\leq\varepsilon$.
Extreme Sufficiency
Anyone who has worked through the proof of the JL lemma by hand should be able to sense that the appearance of $\log N$ in the conclusion of the JL lemma is essentially because the probability term $2\exp\left(-\frac{\varepsilon^2 n}{8}\right)$ in the "unit-norm lemma" decays exponentially, and we could actually relax this to a polynomial decay rate instead, which is where the $\log N$ comes from.
Overall, the JL lemma tells us that to fit $N$ vectors with error $\varepsilon$, we only need a space of dimension $\mathcal{O}\left(\frac{\log N}{\varepsilon^2}\right)$, and the exact constant in front of $\frac{\log N}{\varepsilon^2}$ isn't actually that important. This is because the JL lemma is in fact a very sufficient condition, and in practice the actual requirements are often much more relaxed. For example, in the proof of the JL lemma, if we change the condition to $n > \frac{16\log N}{\varepsilon^2}$, then the probability that $\eqref{eq:bound}$ holds becomes no less than
\begin{equation}1 - N(N-1)\exp\left(-\frac{\varepsilon^2 n}{8}\right)\geq 1 - N(N-1)N^{-2}=1/N\end{equation}
Note that although $1/N$ is small, it is nonetheless greater than 0, so there still exists a $A$ that makes $\eqref{eq:bound}$ hold — it's just that the cost of finding such a $A$ is higher (since each attempt only succeeds with probability $1/N$). But if we only care about existence, that's already enough.
Moreover, the JL lemma only considers dimensionality reduction under random linear projections, and it already achieves $n > \frac{16\log N}{\varepsilon^2}$. If we use more refined reduction methods, such as SVD-based reduction, it's possible to get even better results (with a smaller leading coefficient); and if nonlinear reduction methods are considered as well, the results can become even better still. So there's no need to worry too much about the exact constant in front of $\frac{\log N}{\varepsilon^2}$ — we just need to know the order of magnitude of $\mathcal{O}\left(\frac{\log N}{\varepsilon^2}\right)$. If we actually need to apply it in practice, we usually need to determine the constant based on the specific situation, rather than directly invoking the theoretical result.
To Be Continued
In this article, we introduced the Johnson–Lindenstrauss lemma (JL lemma), an important and remarkable result about dimensionality reduction, and one of the key manifestations of how unusual high-dimensional spaces really are. It tells us that "we only need a $\mathcal{O}(\log N)$-dimensional space to fit $N$ vectors," allowing retrieval problems that were originally posed in a high-dimensional space to be reduced to a $\mathcal{O}(\log N)$-dimensional space instead.
This article mainly discussed the theoretical proof details of the JL lemma. In the next article, we'll try to apply it to understanding some machine learning problems — stay tuned!
Translated automatically with claude-sonnet-5; all equations are reproduced verbatim from the source. Copyright remains with the original author.