Rotate, Round, and Debias
10 min read
Embedding search represents every stored item and every query by a vector. After the vectors are normalized, a larger inner product means a better match. The search system ranks items by that score; it does not need to reconstruct every stored vector perfectly.
Storage is the obstacle. A billion vectors in 768 dimensions occupy about three terabytes in 32-bit floating point. At one bit per coordinate, the raw codes occupy only 96 gigabytes. That saving is useful if the compressed vectors preserve the scores well enough to preserve the ranking.
Coordinatewise rounding is tempting, but it creates two separate problems. First, a fixed rounding grid treats some directions much better than others. Second, even small reconstruction errors can systematically shrink or shift inner products. Random rotation addresses the first problem; debiasing addresses the second.
RaBitQ [GL, SIGMOD ’24], its arbitrary-bit extension [GGXYLW, SIGMOD ’25], and TurboQuant [ZDHM, ICLR ’26] develop versions of this recipe. The four questions below build it one step at a time.
What’s Wrong with a Fixed Grid?
Why should the direction a vector points in matter to a rule that rounds one coordinate at a time? Start in two dimensions. The four one-bit codewords are the normalized corners of a square:
\[ \left(\pm\frac1{\sqrt2},\pm\frac1{\sqrt2}\right). \]
All four corners lie on the unit circle. A unit vector pointing toward a corner is represented exactly, while an equally long vector pointing along an axis must move to a corner. For the axis vector \((1,0)\), that squared gap is \(2-\sqrt2\approx0.586\).
The fixed grid therefore favors its diagonals and treats its axes badly. Higher-dimensional hypercubes have the same directional problem. An embedding distribution may happen to avoid the bad directions, but the grid itself gives no reason it should.
Now let \(\mathbf o\in\mathbb R^d\) be one normalized stored vector and let \(\mathbf q\in\mathbb R^d\) be a normalized query. Sample an orthogonal matrix \(\mathbf P\in\mathbb R^{d\times d}\) uniformly from Haar measure and rotate the stored vector before rounding:
\[ \mathbf u=\mathbf P^\top\mathbf o. \]
The same \(\mathbf P\) is shared across the database and the queries. Before rounding it changes no geometry, since:
\[ \langle\mathbf P^\top\mathbf o,\mathbf P^\top\mathbf q\rangle =\langle\mathbf o,\mathbf q\rangle. \]
Its role is to randomize how the fixed grid meets the data. For every fixed \(\mathbf o\), the rotated vector \(\mathbf u\) is uniform on the unit sphere. The diagonal and the axis now face the same distribution of rounding error. The rounding error remains, but the code’s alignment with the original vector is a number the encoder can compute, and dividing the code’s inner product with a query by that alignment gives an exactly unbiased score; what remains is its variance, and a control variate that lowers it.
What Does a Random Rotation Buy?
The rotated vector is uniform on the sphere, so what does one of its coordinates look like? A uniform point on the sphere can be generated from a standard Gaussian vector \(\mathbf g\sim\mathcal N(\mathbf0,\mathbf I_d)\):
\[ \mathbf u\ \stackrel{d}{=}\ \frac{\mathbf g}{\|\mathbf g\|_2}. \]
The norm in the denominator is about \(\sqrt d\), so a single coordinate is approximately \(g_i/\sqrt d\). The exact variance and concentration explain why a scalar quantizer can be reused across all coordinates.
Claim: For every fixed unit vector \(\mathbf o\) and Haar-random \(\mathbf P\), each coordinate of \(\mathbf u=\mathbf P^\top\mathbf o\) satisfies:
\[ \mathbb E[u_i]=0, \qquad \operatorname{Var}(u_i)=\frac1d. \]
For any failure probability \(0<\delta<1\), a union bound gives:
\[ \max_{i\in[d]}|u_i| =O\!\left(\sqrt{\frac{\log(d/\delta)}{d}}\right) \]
with probability at least \(1-\delta\).
Proof
The distribution of \(\mathbf u\) is unchanged by reflecting any coordinate, so \(u_i\) has mean zero. It is also unchanged by permuting coordinates, so all \(u_i^2\) have the same expectation. The coordinates have unit total squared length on every draw:
\[ \sum_{i=1}^d u_i^2=1. \]
Taking expectations and using exchangeability gives:
\[\begin{align} d\,\mathbb E[u_i^2] &=1, \\ \mathbb E[u_i^2] &=\frac1d. \end{align}\]
The standard spherical-cap inequality says that, for a numerical constant \(c>0\) and threshold \(0<t<1\), one coordinate obeys \(\Pr(|u_i|>t)\leq2e^{-cdt^2}\). Applying it to all \(d\) coordinates gives:
\[ \Pr\!\left(\max_i|u_i|>t\right) \leq 2d e^{-cdt^2}. \]
Choosing \(t\) proportional to \(\sqrt{\log(d/\delta)/d}\) proves the final statement.
The Gaussian approximation follows directly from the normalized-Gaussian construction:
\[ \sqrt d\,u_i = \frac{g_i}{\sqrt{d^{-1}\sum_{j=1}^d g_j^2}}. \]
The denominator approaches one by the law of large numbers. Slutsky’s theorem then gives \(\sqrt d\,u_i\Rightarrow\mathcal N(0,1)\) as \(d\) grows. The standardized coordinate \(\sqrt d\,u_i\) has variance one for every finite \(d\), not only in the limit. The coordinates are dependent because their squares must sum to one, but linearity of expectation still lets one shared scalar rule control their total mean-squared error.
For one bit, reconstruct the vector by rounding each rotated coordinate to its sign and rotating back:
\[ \bar{\mathbf o} =\mathbf P\frac{\operatorname{sign}(\mathbf u)}{\sqrt d}. \]
The reconstruction is a unit vector. Its realized alignment with the original vector is:
\[ a =\langle\bar{\mathbf o},\mathbf o\rangle =\frac{\|\mathbf u\|_1}{\sqrt d}. \]
The second equality writes \(\mathbf o=\mathbf P\mathbf u\), cancels \(\mathbf P\) against \(\mathbf P^\top\) because an orthogonal matrix preserves inner products, and uses \(\operatorname{sign}(u_i)\,u_i=|u_i|\) in each coordinate. The corrected score divides by this alignment, so the alignment has to be stable and stay away from zero.
Claim: The one-bit alignment approaches a constant in probability and in expectation:
\[ \mathbb E[a] \longrightarrow \sqrt{\frac2\pi}. \]
Moreover, deviations from its mean become exponentially unlikely as \(d\) grows.
Proof
Use the normalized-Gaussian representation of \(\mathbf u\) to write:
\[ a =\frac{d^{-1}\sum_{i=1}^d|g_i|} {\sqrt{d^{-1}\sum_{i=1}^d g_i^2}}. \]
The numerator approaches \(\mathbb E|g_i|=\sqrt{2/\pi}\) and the denominator approaches one, so the ratio converges in probability. Because \(0\leq a\leq1\), its expectation converges to the same limit.
Finally, the function \(\|\mathbf u\|_1/\sqrt d\) is \(1\)-Lipschitz on the sphere because Cauchy–Schwarz gives \(\|\mathbf x-\mathbf y\|_1\leq\sqrt d\|\mathbf x-\mathbf y\|_2\). Concentration of Lipschitz functions on the sphere gives the stated exponential concentration.
At moderate dimension, \(a\) is already close to \(\sqrt{2/\pi}\approx0.798\). The rotation has turned the fixed grid’s direction-dependent failure into a stable random alignment.
How Noisy Is the Debiased Inner Product?
With the alignment pinned near \(0.8\), how far does the corrected score land from the true inner product? The true score for the stored vector and query is their inner product. Call it \(\rho\) in this section because it appears throughout the error calculation:
\[ \rho=\langle\mathbf o,\mathbf q\rangle. \]
Both vectors are unit length, so \(-1\leq\rho\leq1\). Decompose the query into its component along \(\mathbf o\) and one perpendicular unit direction \(\mathbf e\):
\[ \mathbf q =\rho\mathbf o +\sqrt{1-\rho^2}\,\mathbf e, \qquad \langle\mathbf e,\mathbf o\rangle=0. \]
Taking the inner product with \(\bar{\mathbf o}\) separates the desired signal from the perpendicular rounding error:
\[ \langle\bar{\mathbf o},\mathbf q\rangle =a\rho +\sqrt{1-\rho^2}\,\langle\bar{\mathbf o},\mathbf e\rangle. \]
The stored code can therefore report the alignment-corrected estimate:
\[ \hat\rho =\frac{\langle\bar{\mathbf o},\mathbf q\rangle}{a}. \]
The correction is exact when the query is parallel to the stored vector. For every other query, Haar symmetry determines the full error distribution.
Claim: Let \(X\) be the first coordinate of a point drawn uniformly from the unit sphere in \(\mathbb R^{d-1}\). The variable \(X\) is independent of \(a\), has mean zero and variance \(1/(d-1)\), and:
\[ \hat\rho \ \stackrel{d}{=}\ \rho +\sqrt{1-\rho^2}\, \frac{\sqrt{1-a^2}}{a}\,X. \]
Consequently the estimator is unbiased and its conditional variance is:
\[ \operatorname{Var}(\hat\rho\mid a) =\frac{1-\rho^2}{d-1}\frac{1-a^2}{a^2}. \]
Its unconditional variance satisfies:
\[ \operatorname{Var}(\hat\rho) =\frac{1-\rho^2}{d-1} \mathbb E\!\left[\frac{1-a^2}{a^2}\right] \sim \left(\frac\pi2-1\right)\frac{1-\rho^2}{d}. \]
Proof
Condition on \(\mathbf u=\mathbf P^\top\mathbf o\). The rotated sign code \(\mathbf s=\operatorname{sign}(\mathbf u)/\sqrt d\) is then fixed. Split it into its component along \(\mathbf u\) and its perpendicular residual:
\[ \mathbf s =a\mathbf u+\mathbf w, \qquad \langle\mathbf w,\mathbf u\rangle=0. \]
Both \(\mathbf s\) and \(\mathbf u\) have unit norm, so Pythagoras gives:
\[ \|\mathbf w\|_2^2=1-a^2. \]
The vector \(\mathbf P^\top\mathbf e\) is a uniformly random unit direction in the \((d-1)\)-dimensional subspace perpendicular to \(\mathbf u\). Its inner product with the fixed residual \(\mathbf w\) therefore has distribution:
\[ \langle\bar{\mathbf o},\mathbf e\rangle =\langle\mathbf s,\mathbf P^\top\mathbf e\rangle \ \stackrel{d}{=}\ \sqrt{1-a^2}\,X. \]
Substituting this identity into the decomposition of \(\hat\rho\) proves the distributional formula. Spherical symmetry gives \(\mathbb E[X]=0\), while the same sum-of-squares argument used above gives \(\mathbb E[X^2]=1/(d-1)\). These two moments prove unbiasedness and the conditional variance.
The law of total variance gives the unconditional formula because the conditional mean is always \(\rho\). Finally, \(a\) concentrates at \(\sqrt{2/\pi}\), so its variance factor converges to:
\[ \frac{1-2/\pi}{2/\pi} =\frac\pi2-1. \]
The error is largest near \(\rho=0\), vanishes for parallel vectors, and shrinks as \(d^{-1/2}\). The plot shows both pieces of the calculation: alignment remains near \(0.8\), while the perpendicular component and the inner-product error contract with dimension.
The division argument is not tied to sign rounding. Any fixed rounding rule that maps \(\mathbf u\) to a reconstruction independently of the query inherits the reflection symmetry that makes the perpendicular term mean zero. Division by its realized alignment is unbiased whenever the ratio is integrable. Its variance can still be poor when that alignment is small or its perpendicular residual is large.
Can a Control Variate Push It Further?
One bit has fixed the bias; if the budget allows a few more bits per coordinate, where should they go? With a budget of \(b\) bits per coordinate, it is natural to choose the reconstruction that minimizes mean-squared error. Lloyd–Max quantization partitions the rotated-coordinate density into \(2^b\) intervals and represents each interval by its conditional mean [Lloyd, IEEE-IT ’82].
At one bit, the large-\(d\) levels are:
\[ \pm\sqrt{\frac{2}{\pi d}}. \]
They are the sign codewords shrunk by \(\sqrt{2/\pi}\). That shrinkage reduces reconstruction error, but for the Lloyd–Max reconstruction \(\widetilde{\mathbf o}\) the uncorrected expected inner product becomes approximately:
\[ \mathbb E\langle\widetilde{\mathbf o},\mathbf q\rangle =\frac2\pi\rho. \]
Mean-squared reconstruction error and unbiased inner products are different objectives.
Let \(\widetilde{\mathbf o}\) be any coarse reconstruction and keep its residual conceptually separate:
\[ \mathbf r=\mathbf o-\widetilde{\mathbf o}. \]
The target inner product splits without approximation:
\[ \rho =\langle\widetilde{\mathbf o},\mathbf q\rangle +\langle\mathbf r,\mathbf q\rangle. \]
Compute the first term from the coarse code and spend randomness only on the second. This split is a control-variate construction: a cheap approximation carries the predictable part, and the randomized estimator sees only the smaller remainder.
Quantized Johnson–Lindenstrauss (QJL) supplies a one-bit residual estimator [ZDH, AAAI ’25]. Draw \(m\) independent standard Gaussian directions \(\mathbf s_1,\ldots,\mathbf s_m\). Store \(\|\mathbf r\|_2\) and the \(m\) signs of their inner products with \(\mathbf r\). The query uses:
\[ \hat\rho^{\mathrm{cv}} =\langle\widetilde{\mathbf o},\mathbf q\rangle +\sqrt{\frac\pi2}\frac{\|\mathbf r\|_2}{m} \sum_{j=1}^m (\mathbf s_j^\top\mathbf q) \operatorname{sign}(\mathbf s_j^\top\mathbf r). \]
Claim: Conditional on the coarse reconstruction, the residual estimator is unbiased and has exact variance:
\[ \mathbb E[\hat\rho^{\mathrm{cv}}] =\rho, \]
\[ \operatorname{Var}(\hat\rho^{\mathrm{cv}}) =\frac1m\left( \frac\pi2\|\mathbf r\|_2^2 -\langle\mathbf r,\mathbf q\rangle^2 \right) \leq \frac\pi{2m}\|\mathbf r\|_2^2. \]
Proof
If \(\mathbf r=\mathbf0\), the claim is immediate. Suppose \(\mathbf r\neq\mathbf0\). For one Gaussian row \(\mathbf s\), rotational symmetry gives the identity:
\[ \mathbb E\!\left[ (\mathbf s^\top\mathbf q) \operatorname{sign}(\mathbf s^\top\mathbf r) \right] =\sqrt{\frac2\pi}\, \frac{\langle\mathbf r,\mathbf q\rangle}{\|\mathbf r\|_2}. \]
Multiplication by \(\sqrt{\pi/2}\|\mathbf r\|_2\) makes the row contribution unbiased for \(\langle\mathbf r,\mathbf q\rangle\). Adding the exact coarse term proves \(\mathbb E[\hat\rho^{\mathrm{cv}}]=\rho\).
Squaring removes the sign, so the row contribution has second moment:
\[\begin{align} \mathbb E\!\left[ \frac\pi2\|\mathbf r\|_2^2 (\mathbf s^\top\mathbf q)^2 \right] &=\frac\pi2\|\mathbf r\|_2^2 \mathbb E[(\mathbf s^\top\mathbf q)^2] \\ &=\frac\pi2\|\mathbf r\|_2^2. \end{align}\]
The last line uses \(\|\mathbf q\|_2=1\).
Subtracting the squared mean gives the single-row variance. The \(m\) rows are independent, so averaging divides that variance by \(m\).
The variance now charges only for \(\|\mathbf r\|_2^2\), the error left by the coarse quantizer. TurboQuant spends \(b-1\) bits per coordinate on the Lloyd–Max reconstruction and one bit on a QJL residual with \(m=d\). Its mean-squared residual falls as \(4^{-(b-1)}\), so the inner-product variance falls as \(4^{-b}/d\). This matches the optimal dependence on both bit width and dimension; the paper’s worst-case upper bound is \(\sqrt3\pi^2/(d4^b)\) for a unit query, compared with a lower bound of \(1/(d4^b)\).
Hutch++ uses the same architecture for matrix traces: compute the dominant low-rank part directly, then apply Hutchinson’s random probes only to the spectral tail. In both constructions, approximation alone may be biased, while an unbiased estimator of the residual restores the target and pays variance only for what the approximation missed.
Nearest-neighbor search uses alignment division when one-bit codes and sharp per-vector confidence bounds are the priority. Transformer key-value cache compression uses the residual estimator because attention reads keys through inner products. EDEN rotates and quantizes client updates before federated averaging [VBPMBM, ICML ’22], using scaling to preserve the mean update. Across these settings, the reusable design is the same: randomize the direction, spend the code on a small reconstruction residual, and correct the statistic the application actually reads.