Inference for Estimated Reduction Subspaces

Why Infer a Subspace Rather Than a Direction?

Dimension-reduction methods such as SIR and SAVE often estimate a target subspace from the leading eigenspace of a symmetric matrix. The scientific target is the span of the relevant directions, not necessarily any particular ordered list of vectors spanning it. This distinction matters for inference: changing signs, or rotating a basis within the same subspace, can change the reported vectors without changing the reduction at all.

A convenient object for inference is therefore the projection onto the subspace. It is unchanged by a change of orthonormal basis, and its entries have a first-order perturbation formula that makes sampling uncertainty explicit. The formula also shows which features of the matrix estimate affect the subspace: to first order, only perturbations linking the selected eigenspace to its complement matter.

Definition: Let \(\mathcal S\) be a \(d\)-dimensional subspace of \(\mathbb R^p\), and let \(U\) be a \(p \times d\) matrix whose columns are an orthonormal basis for \(\mathcal S\). The orthogonal projection matrix onto \(\mathcal S\) is \(P=UU^T\). It satisfies \(P^T=P\), \(P^2=P\), and \(\operatorname{rank}(P)=d\). If \(U\) is replaced by \(UQ\), where \(Q\) is any \(d \times d\) orthogonal matrix, then \((UQ)(UQ)^T=P\).

In a reduction problem, suppose the population matrix \(M\) has a leading \(d\)-dimensional eigenspace equal to the target reduction subspace. Write its eigenvalues in decreasing order as \(\lambda_1\geq\cdots\geq\lambda_p\), and let \(\widehat M\) be a symmetric estimator. The estimated subspace is the leading \(d\)-dimensional eigenspace of \(\widehat M\), and \(\widehat P\) denotes its orthogonal projection. This interpretation requires that the population eigenspace genuinely identify the target; if the method's population operator does not recover the full central subspace, inference for its eigenspace does not by itself establish inference for the full central subspace.

First-Order Perturbation of a Spectral Projection

The key regularity condition is a positive gap between the selected eigenvalues and the unselected ones. Eigenvalues within either group may be repeated; individual eigenvectors may then be non-unique. The gap between the two groups, rather than distinctness of every eigenvalue, makes the selected subspace locally identifiable.

Theorem: Let \(M\) be a fixed real symmetric \(p \times p\) matrix with \(\lambda_d>\lambda_{d+1}\), and let \(P\) be its leading \(d\)-dimensional spectral projection. Suppose \(\widehat M\) is symmetric and \(\sqrt n(\widehat M-M)\xrightarrow{d}Z\), where \(Z\) is a finite random symmetric matrix. Let \(u_1,\ldots,u_p\) be an orthonormal eigenbasis of \(M\), with corresponding eigenvalues \(\lambda_1,\ldots,\lambda_p\). Then \[ \sqrt n(\widehat P-P)\xrightarrow{d} \sum_{i=1}^{d}\sum_{j=d+1}^{p} \frac{u_j^TZu_i}{\lambda_i-\lambda_j} \left(u_ju_i^T+u_iu_j^T\right). \] The limit is determined by the cross-eigenspace entries of \(Z\). This statement does not assert that the limit is normally distributed.

Here is a derivation. Express matrices in the eigenbasis of \(M\), with the selected coordinates first, and write \(M=\operatorname{diag}(A,D)\), where \(A=\operatorname{diag}(\lambda_1,\ldots,\lambda_d)\) and \(D=\operatorname{diag}(\lambda_{d+1},\ldots,\lambda_p)\). Write the perturbation \(E=\widehat M-M\) in conformable blocks \(E_{11},E_{12},E_{21},E_{22}\). For sufficiently small \(\|E\|\), the estimated subspace can be represented as the columns of \(\begin{pmatrix}I\\T\end{pmatrix}\) for a small \((p-d)\times d\) matrix \(T\). Invariance of that subspace under \(\widehat M\) gives

E_{21}+(D+E_{22})T=T(A+E_{11}+E_{12}T).

Keeping the terms linear in \(E\) and \(T\) yields \(E_{21}+DT-TA=0\). Thus, for a selected index \(i\leq d\) and unselected index \(j>d\),

T_{ji}=\frac{E_{ji}}{\lambda_i-\lambda_j}+\text{terms of order }\|E\|^2.

The projection onto the columns of \(\begin{pmatrix}I\\T\end{pmatrix}\) is \(\begin{pmatrix}I\\T\end{pmatrix}(I+T^TT)^{-1}\begin{pmatrix}I&T^T\end{pmatrix}\). Since \(T^TT\) is of order \(\|E\|^2\), its first-order change from \(P\) has off-diagonal blocks \(T\) and \(T^T\), and zero diagonal blocks. Substituting the entries of \(T\) gives the displayed linear perturbation. The positive eigenvalue gap ensures the denominators are nonzero and the local subspace representation is valid. Finally, \(\sqrt n E\) converges in distribution and is therefore bounded in probability after scaling, so \(\sqrt n\|E\|^2\xrightarrow{p}0\). Applying the continuous mapping theorem to the linear term proves the stated limit. \(\square\)

This result gives a direct way to construct uncertainty measures. For example, if the entries of \(\sqrt n(\widehat M-M)\) are jointly asymptotically normal, then the displayed linear transformation is also jointly normal. Its covariance follows by applying that transformation to the covariance of the matrix entries. Without asymptotic normality of \(Z\), the perturbation result still gives the limit, but not a normal limit.

Worked Example: A Normal Limit for a Two-Dimensional Subspace Estimate

Let \(M=\begin{pmatrix}3&0\\0&1\end{pmatrix}\), select the leading one-dimensional eigenspace, and suppose \(\sqrt n(\widehat M-M)\xrightarrow{d}Z\), where \(Z\) is symmetric and \(Z_{12}=Z_{21}\sim N(0,4)\). The population projection is \(P=\begin{pmatrix}1&0\\0&0\end{pmatrix}\). The gap is \(3-1=2\), so the perturbation formula gives

\sqrt n(\widehat P-P)\xrightarrow{d} \frac{Z_{12}}{2} \begin{pmatrix} 0&1\\ 1&0 \end{pmatrix}.

The off-diagonal entry has limiting distribution \(Z_{12}/2\), which is \(N(0,4/4)=N(0,1)\). The diagonal entries have zero first-order limit. Thus a first-order 95% normal interval for the population value \(P_{12}=0\), based on the estimated entry, is approximately \(\widehat P_{12}\pm1.96/\sqrt n\). This normal approximation is justified here by the stated normal limit for \(Z_{12}\), together with the eigengap and the asymptotic expansion. It would not follow from a variance calculation alone.

What the Projection Limit Does—and Does Not—Say

The projection records the subspace as a whole. For a one-dimensional subspace, \(u\) and \(-u\) are two equally valid unit eigenvectors, but both give \(uu^T\). For a higher-dimensional subspace, every orthonormal basis \(UQ\) gives the same \(UU^T\). Consequently, uncertainty in \(\widehat P\) measures changes in the subspace, not arbitrary choices of signs or rotations used to report its basis.

A useful geometric summary is the set of principal angles between two subspaces: these angles measure how far the subspaces themselves have moved, regardless of how their bases are chosen. The projection difference is another basis-invariant summary. In the first-order limit above, its nonzero entries link selected and unselected eigenvectors; perturbations internal to the selected group can rotate its basis without changing the subspace to first order.

Worked Example: A Variance Does Not Establish a Normal Confidence Interval

Use \(M=\begin{pmatrix}4&0\\0&2\end{pmatrix}\), again selecting the leading one-dimensional eigenspace. Suppose the only nonzero entries of \(Z\) are \(Z_{12}=Z_{21}\), where \(P(Z_{12}=0.3)=P(Z_{12}=-0.3)=1/2\). Then \(E[Z_{12}]=0\) and \(\operatorname{Var}(Z_{12})=0.09\), but \(Z_{12}\) is not normally distributed. The gap is \(4-2=2\), so the first-order limit is

\sqrt n\,\widehat P_{12}\xrightarrow{d}\frac{Z_{12}}{2}, \qquad P\left(\frac{Z_{12}}{2}=0.15\right) =P\left(\frac{Z_{12}}{2}=-0.15\right)=\frac12.

The limiting mean is zero and its variance is \(0.09/2^2=0.0225\). Its standard deviation is \(0.15\), but the limit is a two-point distribution, not a normal distribution. In particular, the interval \(\widehat P_{12}\pm1.96(0.15)/\sqrt n\) is not justified as an asymptotic 95% interval by these assumptions. The limiting distribution instead places half its mass at each of \(-0.15\) and \(0.15\); its quantiles and coverage cannot be replaced by normal quantiles merely because its variance is known.

For a concrete sequence, take \(\widehat M=M+n^{-1/2}Z\). Its leading projection has off-diagonal entry \(\widehat P_{12}=t/(2\sqrt{1+t^2})\), where \(t=Z_{12}/\sqrt n\). Therefore \(\sqrt n\,\widehat P_{12}=Z_{12}/(2\sqrt{1+Z_{12}^2/n})\), which converges to \(Z_{12}/2\), confirming the first-order calculation.

Several Selected Directions and Basis Invariance

Worked Example: Which Perturbations Change a Two-Dimensional Subspace?

Consider \(M=\operatorname{diag}(5,3,1)\) and select its leading two-dimensional eigenspace. The population projection is \(\operatorname{diag}(1,1,0)\), and the gap separating the selected and unselected eigenvalues is \(3-1=2\). Suppose the only random entries of the symmetric limit matrix \(Z\) that matter are \(Z_{13}=Z_{31}\sim N(0,1)\) and \(Z_{23}=Z_{32}\sim N(0,4)\), independent. The first-order projection limit is

\sqrt n(\widehat P-P)\xrightarrow{d} \frac{Z_{13}}{5-1}(u_3u_1^T+u_1u_3^T) +\frac{Z_{23}}{3-1}(u_3u_2^T+u_2u_3^T),

where \(u_1,u_2,u_3\) are the coordinate vectors. Thus the limiting \((1,3)\) entry is \(Z_{13}/4\), with variance \(1/16\), while the limiting \((2,3)\) entry is \(Z_{23}/2\), with variance \(4/4=1\). A perturbation \(Z_{12}\) linking the two selected eigenvectors does not appear in this first-order projection limit: it can alter the basis within the selected plane, but not the plane itself to first order.

For instance, rotating the estimated basis by any orthogonal \(2\times2\) matrix \(Q\) changes its two columns from \(\widehat U\) to \(\widehat UQ\), but \((\widehat UQ)(\widehat UQ)^T=\widehat U\widehat U^T=\widehat P\). A reported change in basis columns is therefore not, by itself, evidence of increased subspace uncertainty.

Practical Inference and a Common Pitfall

When the limit matrix is asymptotically normal, the perturbation formula supplies a covariance for coordinates of \(\sqrt n(\widehat P-P)\). That covariance can support standard errors or joint inference on the projection, subject to the constraints that a projection matrix is symmetric, idempotent, and of fixed rank. In practice, one must estimate the covariance of the relevant entries of \(\sqrt n(\widehat M-M)\), account for the eigenvalue gaps in the denominators, and transform that covariance through the displayed linear map. If the matrix-estimator limit is nonnormal, inference should use its justified limiting distribution or another method with established coverage, rather than imposing normal quantiles.

Small gaps deserve special attention. As \(\lambda_i-\lambda_j\) shrinks for a selected-unselected pair, the corresponding first-order coefficient grows. The same matrix-estimation noise can then produce much greater subspace uncertainty. If the gap is zero, the selected subspace is not uniquely separated by this rule, and the theorem does not apply. A confidence procedure based on a well-separated eigenspace should not be interpreted as reliable when that separation is absent.

Common mistake: “The perturbation has mean zero and variance \(v\), so the estimated projection entry is asymptotically normal and I can use a normal 95% interval.”
Correction: The perturbation theorem transfers the distribution of the cross-eigenspace entries of \(Z\) through a linear map. Mean and variance alone do not determine that distribution. A normal limit and normal-based interval require asymptotic normality (or a separately justified normal approximation) and the relevant regularity conditions.
Another mistake: Treating a sign change or within-subspace rotation of estimated eigenvectors as uncertainty in the reduction itself. Compare projections or other basis-invariant subspace summaries instead.

Key takeaway: Under a positive eigenvalue gap, first-order uncertainty in an estimated reduction subspace is described by the projection perturbation, with each selected-to-unselected matrix entry scaled by its eigenvalue separation. This quantifies subspace movement, not arbitrary basis variation; the shape of the limiting distribution must be justified before choosing confidence quantiles.

Check Your Understanding

Use the eigengap and the distributional assumptions, not just a standard error, when answering.

Why is a positive gap between \(\lambda_d\) and \(\lambda_{d+1}\) needed for the leading \(d\)-dimensional projection to be locally identifiable?
In the perturbation formula, why do entries of \(Z\) linking two selected eigenvectors not appear in the first-order projection limit?
If a cross-eigenspace entry of \(Z\) has mean zero and variance \(0.09\), what additional information is needed to claim a normal limit?
Explain why replacing an orthonormal basis \(U\) by \(UQ\), for orthogonal \(Q\), does not change the estimated subspace.
What happens to the first-order uncertainty contribution when a selected eigenvalue approaches an unselected eigenvalue?