Random-Design Gram Matrices and Stability

Why the Empirical Gram Matrix Matters

A series regression fit uses the observed basis vectors to estimate its coefficients. Even when the basis is orthonormal under the population design distribution, its columns need not be orthogonal in a finite random sample. The matrix measuring this discrepancy is the empirical Gram matrix. If it is close to the identity, least squares behaves like the population orthogonal projection; if it is poorly conditioned or singular, small changes in the data can produce large coefficient changes, or no unique coefficient estimate at all.

The goal here is to quantify that distinction. The key quantities are the sample size $n$, the basis dimension $m$, and the largest squared length of a basis vector at any design point. The last quantity measures how concentrated the basis can be at a single point and controls the concentration bound.

Setup and Definitions

Definition: Let $X_1,\ldots,X_n$ be independent and identically distributed design points with distribution $P_X$. Let $p_m(x)=(\phi_1(x),\ldots,\phi_m(x))^\mathsf{T}$ be a vector of $m$ basis functions that are orthonormal under $P_X$, meaning $$ E[p_m(X)p_m(X)^\mathsf{T}]=I_m. $$ The empirical Gram matrix is $$ G_n=\frac{1}{n}\sum_{i=1}^n p_m(X_i)p_m(X_i)^\mathsf{T}. $$ The leverage bound is $K_m=\operatorname*{ess\,sup}_{x}\|p_m(x)\|_2^2$, where the essential supremum is taken under $P_X$.

The identity $E[G_n]=I_m$ follows by taking expectations term by term. Also, $K_m\geq m$, since $E[\|p_m(X)\|_2^2]=\sum_{j=1}^m E[\phi_j(X)^2]=m$. A finite $K_m$ is the boundedness condition used below. It holds, for example, for finite trigonometric systems and many finite-dimensional spline bases.

For a response $Y_i$, the series least-squares coefficient estimate minimizes $n^{-1}\sum_{i=1}^n(Y_i-p_m(X_i)^\mathsf{T}b)^2$ over $b\in\mathbb{R}^m$. Its quadratic term in $b$ is $b^\mathsf{T}G_nb$. Thus the eigenvalues of $G_n$ determine whether this objective has a unique minimizer and how sensitive the fit can be to perturbations.

A Concentration Bound for the Gram Matrix

We use the matrix Bernstein inequality in its two-sided self-adjoint form. If $A_1,\ldots,A_n$ are independent, mean-zero, self-adjoint $m\times m$ random matrices, $\|A_i\|_{\mathrm{op}}\leq R$ almost surely, and $\|\sum_i E[A_i^2]\|_{\mathrm{op}}\leq \sigma^2$, then for $t>0$, $$ P\left(\left\|\sum_{i=1}^n A_i\right\|_{\mathrm{op}}\geq t\right) \leq 2m\exp\left(-\frac{t^2}{2\sigma^2+(2/3)Rt}\right). $$ Here $\|\cdot\|_{\mathrm{op}}$ denotes the operator norm, equal for a symmetric matrix to the largest absolute eigenvalue.

To apply it, set $A_i=p_m(X_i)p_m(X_i)^\mathsf{T}-I_m$. These matrices are independent and mean zero. Write $q=p_m(X_i)$. The matrix $qq^\mathsf{T}$ has largest eigenvalue $\|q\|_2^2$; its other eigenvalues are zero. Since $\|q\|_2^2\leq K_m$ almost surely and $K_m\geq m\geq1$, it follows that $\|A_i\|_{\mathrm{op}}\leq K_m$.

For the variance term, $(qq^\mathsf{T})^2=\|q\|_2^2qq^\mathsf{T}$, so $$ E[A_i^2] =E[\|q\|_2^2qq^\mathsf{T}]-I_m \preceq K_m E[qq^\mathsf{T}] =K_m I_m. $$ The notation $B\preceq C$ means that $C-B$ is positive semidefinite. Since $E[A_i^2]$ is positive semidefinite, summing gives $\|\sum_i E[A_i^2]\|_{\mathrm{op}}\leq nK_m$. These bounds yield the following result.

Theorem: Suppose $X_1,\ldots,X_n$ are independent and identically distributed, the $m$ basis functions are orthonormal under $P_X$, and $K_m<\infty$. For $0<\delta<1$, $$ P\left(\|G_n-I_m\|_{\mathrm{op}}>\frac{\delta}{2}\right) \leq 2m\exp\left(-\frac{n\delta^2}{K_m(8+4\delta/3)}\right) \leq 2m\exp\left(-\frac{n\delta^2}{10K_m}\right). $$

For the first inequality, apply matrix Bernstein to $\sum_i A_i=n(G_n-I_m)$ at $t=n\delta/2$, using $R=K_m$ and $\sigma^2=nK_m$. Its exponent is $$ \frac{(n\delta/2)^2}{2nK_m+(2/3)K_m(n\delta/2)} =\frac{n\delta^2}{K_m(8+4\delta/3)}. $$ For the second inequality, $0<\delta<1$ implies $8+4\delta/3<28/3$. Consequently $1/(8+4\delta/3)>3/28>1/10$, which gives the stated, slightly weaker, exponential bound. This comparison is important: the exact denominator, rather than a looser intermediate estimate, justifies the coefficient $1/10$.

What Concentration Guarantees for Least Squares

When $\|G_n-I_m\|_{\mathrm{op}}\leq\delta/2$, every eigenvalue of $G_n$ lies between $1-\delta/2$ and $1+\delta/2$. This follows from the spectral theorem: for every unit vector $v$, the operator-norm bound gives $|v^\mathsf{T}(G_n-I_m)v|\leq\delta/2$, so $1-\delta/2\leq v^\mathsf{T}G_nv\leq1+\delta/2$. Taking the minimum and maximum over unit vectors gives the eigenvalue bounds.

Because $0<\delta<1$, the lower endpoint is positive. Hence $G_n$ is invertible on this event, the least-squares coefficient estimate is unique, and $$ \operatorname{cond}(G_n) =\frac{\lambda_{\max}(G_n)}{\lambda_{\min}(G_n)} \leq\frac{1+\delta/2}{1-\delta/2}<3. $$ Here $\operatorname{cond}$ is the spectral condition number of a positive definite matrix. More generally, if $\|G_n-I_m\|_{\mathrm{op}}\leq\eta<1$, then its condition number is at most $(1+\eta)/(1-\eta)$.

Thus a useful sufficient condition for stability with high probability is that $n\delta^2/K_m$ be large compared with $\log(2m)$. The logarithm arises from the matrix dimension factor $2m$ in the tail bound. This is a sufficient condition, not a universal sharp threshold: the bound may be conservative for a particular basis or design distribution.

Worked Example: A Two-Term Polynomial Basis

Let $X$ be uniform on $[0,1]$, and use $\phi_1(x)=1$ and $\phi_2(x)=\sqrt{3}(2x-1)$. The first function has squared integral $1$. The second has mean zero, and its squared integral is $3\int_0^1(2x-1)^2\,dx=3(1/3)=1$; hence the functions are orthonormal. Their squared vector length is $1+3(2x-1)^2$, so $K_2=4$.

Take $n=1000$ and $\delta=1/2$. The theorem gives $$ P\left(\|G_n-I_2\|_{\mathrm{op}}>1/4\right) \leq 4\exp\left(-\frac{1000(1/2)^2}{4(8+4(1/2)/3)}\right) =4\exp\left(-\frac{250}{104/3}\right) =4\exp(-7.2115\ldots) <0.0030. $$ On the complementary event, every eigenvalue lies in $[3/4,5/4]$, and the condition number is at most $(5/4)/(3/4)=5/3$. This illustrates how the probability guarantee translates directly into a bound on the curvature of the least-squares objective.

Examples of Leverage and Finite-Sample Instability

Worked Example: Histogram Basis and Sparse Bins

Let $X$ be uniform on $[0,1]$, divide the interval into $m=4$ equal bins $B_1,\ldots,B_4$, and define $\phi_j(x)=2\,1\{x\in B_j\}$. Each function has squared integral $4(1/4)=1$, and different functions have disjoint supports, so this is an orthonormal basis. At each $x$, exactly one function is nonzero, giving $\|p_4(x)\|_2^2=4$ and $K_4=4$.

For a sample of $n=8$ observations with bin counts $(3,2,2,1)$, the empirical Gram matrix is diagonal. Its $j$th diagonal entry is $n^{-1}\sum_i4\,1\{X_i\in B_j\}=4N_j/8$. Therefore $$ G_8=\operatorname{diag}(1.5,1,1,0.5). $$ The eigenvalues are the diagonal entries, and the condition number is $1.5/0.5=3$. If instead one bin has count zero, the corresponding diagonal entry is zero and the matrix is singular. Orthonormality in the population does not prevent such a sample-level failure; enough observations must reach every relevant part of the design space.

Worked Example: A Discrete Basis with a Singular Realization

Suppose $X$ takes values $-1$ and $1$, each with probability $1/2$, and let $p_2(X)=(1,X)^\mathsf{T}$. Since $E[X]=0$ and $E[X^2]=1$, the two coordinates are orthonormal, and $K_2=2$.

For a sample of six design points containing four $+1$'s and two $-1$'s, the sample mean of the $X_i$'s is $1/3$. Direct calculation gives $$ G_6=\begin{pmatrix}1&1/3\\1/3&1\end{pmatrix}. $$ Its eigenvalues are $1+1/3=4/3$ and $1-1/3=2/3$, so the condition number is $2$. If all six values are $+1$, then $G_6=\begin{pmatrix}1&1\\1&1\end{pmatrix}$, whose eigenvalues are $2$ and $0$. The basis remains population-orthonormal, but this realization contains no variation in the second coordinate after accounting for the first. The probability of this singular event is $2(1/2)^6=1/32$, counting the all-$+1$ and all-$-1$ samples.

Why Excessive Dimension Causes Instability

There is a basic obstruction that does not require a concentration inequality. The matrix $G_n$ is a sum of $n$ matrices of rank at most one, so $\operatorname{rank}(G_n)\leq n$. If $m>n$, it cannot have full rank and must be singular. Even when $m\leq n$, a basis with high leverage can have directions that are poorly represented in the sample. The histogram example illustrates this: empty bins create zero empirical eigenvalues.

The concentration theorem gives a complementary warning. For many bounded systems, $K_m$ grows with $m$; then increasing the number of basis functions also weakens the sufficient concentration condition. A basis can be population-orthonormal at every dimension while its empirical Gram matrix becomes less stable as the dimension grows. Choosing $m$ therefore involves not only approximation quality, but also whether the sample supports reliable estimation of all $m$ coefficients.

One common mistake is to infer that $E[G_n]=I_m$ means $G_n$ is close to $I_m$ for a given sample. The expectation describes an average over repeated samples, not a guarantee for each realization. The correction is to quantify fluctuation using a probability bound and inspect its dependence on $n$, $m$, and $K_m$. A second mistake is to treat $m\leq n$ as sufficient for stability. It only removes the automatic rank obstruction; it does not rule out nearly dependent empirical columns or empty regions of the design space.

Key takeaway: Population orthonormality gives $E[G_n]=I_m$, while finite-sample stability requires the realized eigenvalues of $G_n$ to stay away from zero. Bounded leverage and a sufficiently large sample provide a useful concentration guarantee, but excessive dimension or sparse design coverage can still make least squares unstable.

Check Your Understanding

Use the definitions and results above to answer these questions.

Why is the expected empirical Gram matrix the identity when the basis is orthonormal under the design distribution?
For $0<\delta<1$, what eigenvalue interval follows from $\|G_n-I_m\|_{\mathrm{op}}\leq\delta/2$, and why does it imply a unique least-squares coefficient estimate?
For the four-bin histogram basis, what is $K_4$, and what happens to $G_n$ if one bin has no observations?
Why must $G_n$ be singular whenever $m>n$, regardless of population orthonormality?
What role does the leverage bound $K_m$ play in the concentration guarantee?