Principal Nested Spheres Analysis

• PNS is an intrinsic dimension reduction method that constructs a backward sequence of optimally fitting nested subspheres to capture curved data variations.
• It generalizes PCA on spheres by minimizing geodesic residuals and projecting data onto lower-dimensional embedded subspheres.
• Applications span functional data analysis, shape analysis, molecular dynamics, and high-dimensional omics, with fast implementations addressing computational challenges.

Principal Nested Spheres (PNS) is an intrinsic dimension reduction methodology for data that reside on a sphere or, more broadly, on Riemannian manifolds. It generalizes principal component analysis (PCA) by constructing a backward sequence of optimally fitting embedded subspheres, providing a nested, non-linear decomposition that is intrinsic to the spherical geometry. PNS is designed to capture modes of variation that standard linear techniques cannot, particularly in data exhibiting non-geodesic or "curved" structures, such as those in shape analysis, functional data registration, molecular dynamics, and high-dimensional genomics (Lu et al., 2013, Monem et al., 11 Nov 2025, Dryden et al., 2019).

1. Theoretical Foundation and Mathematical Structure

Let $S^d = \{x \in \mathbb{R}^{d+1} : \|x\| = 1\}$ denote the unit sphere in $(d+1)$-dimensional Euclidean space, equipped with the geodesic (great-circle) metric $\rho(x, y) = \arccos(x^T y)$. The essential object of PNS is a sequence of nested subspheres of decreasing dimension: $$S^d \supset S^{d-1}(c_1, r_1) \supset \cdots \supset S^0(c_d, r_d)$$ where each $S^{k}(c, r) = \{x \in S^d : \|x - c\|_2 = r\}$ defines a $k$-dimensional subsphere with center $c \in \mathbb{R}^{d+1}$, radius $r \in (0, \pi/2]$, and $k = d - j + 1$ at the $j$-th step (Lu et al., 2013, Monem et al., 11 Nov 2025, Dryden et al., 2019).

The PNS decomposition is "backward": at each stage, given data $x_1, \ldots, x_n \in S^k$, the algorithm solves

$$(c_j, r_j) = \arg \min_{c \in \mathbb{R}^{k+1}, r > 0} \sum_{i=1}^n d_{S^k}\left(x_i, S^{k-1}(c, r)\right)^2$$

where $$d_{S^k}(x, S^{k-1}(c, r)) = \min_{y \in S^{k-1}(c, r)} \arccos(\langle x, y \rangle)$$, thereby finding at each stage the best fitting $(k - 1)$-subsphere in the geodesic least-squares sense.

2. Algorithmic Workflow

PNS proceeds recursively:

1. Initialization: Begin with data $\{x_i\} \subset S^d$.
2. Subsphere Fit: At dimension $k$, fit a $(k-1)$-subsphere $S^{k-1}(c_j, r_j)$ minimizing the sum of squared geodesic residuals.
3. Orthogonal Projection: Project each $x_i$ to $S^{k-1}(c_j, r_j)$ along the geodesic from $c_j$, yielding $y_i^{(j)}$.
4. Recursion: Treat the projected data as points on a lower-dimensional sphere; repeat steps 2 and 3 until reaching $S^0$ (a pair of antipodal points).
5. Score Extraction: The signed geodesic distance from each data point to the subsphere at each step constitutes the PNS principal component score.

Explicitly, for projection, given $u = (x - c)/\|x - c\|$, the projected point is $c + r u$, subsequently normalized to unit length (Lu et al., 2013). In ambient coordinates, projections can also be given by

$$\pi_d(x_i) = \frac{\sin(r)}{\sin(\rho_i)} x_i + \frac{\sin(\rho_i - r)}{\sin(\rho_i)} c$$

with $\rho_i = \arccos(c^T x_i)$ (Eltzner et al., 2015, Dryden et al., 2019).

Backward nesting ensures that at each stage, the residual sum-of-squares geodesic variance is tightly controlled, partitioning total variation into contributions along each nested component (Lu et al., 2013, Dryden et al., 2019).

3. Model Selection, Overfitting Control, and Statistical Theory

Distinguishing between a great subsphere ($r = \pi/2$) and a small subsphere ($r < \pi/2$) is critical to avoid overfitting noise. Several formal hypothesis tests are deployed:

• Kolmogorov-Smirnov Test: Compares empirical distributions of geodesic residuals under great- and small-subsphere fits.
• Variance Test: Analyzes the ratio of sample variances of absolute residuals for the two fits.
• Likelihood Ratio Test: Evaluates the improvement in log-likelihood from allowing a small-subsphere parameterization, with the null model being the great subsphere (Monem et al., 11 Nov 2025, Eltzner et al., 2015).

PNS and its generalization to backward nested families of descriptors (BNFD) possess strong asymptotic guarantees. Huckemann and Eltzner show almost sure consistency and joint asymptotic normality for the nested descriptors, including CLTs for the last nested mean or geodesic, provided certain regularity, differentiability, and coercivity conditions hold. Factoring charts allow these joint CLTs to be projected to marginal CLTs for the final descriptor (Huckemann et al., 2016).

A nested bootstrap procedure can be constructed for inference or two-sample testing on manifold-valued data, leveraging the full chain of nested descriptors (Huckemann et al., 2016).

4. Computational Methodology and High-Dimensional Scaling

Direct PNS on high-dimensional spheres is computationally demanding. Fast PNS approaches address this by first projecting the data into a lower-dimensional sphere using a preliminary tangent-space PCA at the Euclidean mean. Specifically:

• Compute the raw mean $\mu$ and principal subspace in the tangent plane via eigen-decomposition.
• Project observations into this subspace and map back to the sphere with

$X_i^* = \mu\,\cos(s_i) + (U_i/s_i)\,\sin(s_i)$

where $U_i$ is the low-rank representation and $s_i = \|U_i\|$ (Monem et al., 11 Nov 2025).

This reduces computational complexity from $O(n\,d^2)$ to $O(n\,p^2)$, where $p \ll d$, making PNS feasible for problems with thousands of variables, as encountered in cancer genomics or proteomics applications (Monem et al., 11 Nov 2025).

The fitting of each subsphere generally relies on alternating minimization: for fixed center, compute the optimal radius as the mean angular distance; for fixed radius, update the center via gradient descent on the sphere (Eltzner et al., 2015, Dryden et al., 2019).

5. Applications and Domain-Specific Variants

Functional Data Analysis

PNS is particularly effective for analyzing horizontal variation in time-warped functional data. After domain warping functions are mapped to Hilbert spheres via the SRVF, PNS finds modes of phase variation more efficiently and interpretably than tangent-based methods (Lu et al., 2013).

Shape Analysis and Molecular Dynamics

PNS generalizes to shape spaces (via Procrustes embedding onto spheres), replacing the forward, subspace-based structure of tangent-PCA with backward, nested subsphere decompositions. For molecular dynamics, the principal nested shape space (PNSS) variant reveals clustering and state transitions among protein conformations more distinctly than linear projections (Dryden et al., 2019).

High-dimensional Omics Data

In genomics and proteomics, PNS (with fast-PNS acceleration) extracts small sets of highly explanatory features and achieves strong classification performance using only a handful of low-dimensional PNS scores, with variable selection guided by the PNS biplot (Monem et al., 11 Nov 2025).

Torus and Manifold Data

Applications to torus-valued data are enabled by a combination of torus-to-sphere deformations, pre-clustering to avoid singularities, and PNS analysis, followed by post-mode hunting for clustering on principal arcs (Eltzner et al., 2015).

Table: Selected PNS Applications

Application Area	Sphere Dimension	Key Adaptations
Functional Data	$\infty$ (Hilbert)	SRVF mapping
Molecular Dynamics	High/3D	Procrustes/PNSS
Cancer Omics	500–12,478	Fast-PNS, biplot
RNA/Torus Data	Moderate	Deformation, clustering

6. Comparison With Related Methods

Tangent-PCA (Principal Geodesic Analysis/PGA):

• Linearizes the manifold at a mean, performs Euclidean PCA in the tangent space, then maps components geodesically back to the sphere (Lu et al., 2013, Eltzner et al., 2015).
• Captures only great-sphere (zero-curvature) modes, often requiring several components to approximate curved variance.

Geodesic-PCA:

• Seeks geodesics maximizing variance but may suffer from interpretive issues such as winding or dense coverage in non-simply connected data (Eltzner et al., 2015).

PNS:

• Allows arbitrary small or great subspheres with curved geometry, often explaining more variance with fewer components and providing more interpretable principal arcs in cases where data cluster around non-geodesic structure (Lu et al., 2013, Dryden et al., 2019).

7. Extensions, Visualization, and Interpretability

PNS Scores and Biplots

PNS scores are constructed as scaled signed geodesic gaps at each stage. These may be circular or real-valued, corresponding to angular coordinates on nested subspheres. For interpretability, the PNS biplot visualizes how changes in a given PNS principal coordinate affect each original variable, via the "back-fitting" mapping from score space to data space. Long curves in the biplot highlight variables most associated with specific principal modes (Monem et al., 11 Nov 2025).

BNFD Framework

PNS fits within the more general BNFD framework for nested submanifolds with strong theoretical properties, supporting formal inference and bootstrapped two-sample tests (Huckemann et al., 2016).

• Principal Nested Spheres for Time Warped Functional Data Analysis (Lu et al., 2013)
• Principal nested spheres for high-dimensional data (Monem et al., 11 Nov 2025)
• Backward Nested Descriptors Asymptotics with Inference on Stem Cell Differentiation (Huckemann et al., 2016)
• Torus Principal Component Analysis with an Application to RNA Structures (Eltzner et al., 2015)
• Principal nested shape space analysis of molecular dynamics data (Dryden et al., 2019)