In structure from motion, quadrifocal tensors capture more information than their pairwise counterparts (essential matrices), yet they have often been thought of as impractical and only of theoretical interest. In this work, we challenge such beliefs by providing a new framework to recover n cameras from the corresponding collection of quadrifocal tensors. We form the block quadrifocal tensor and show that it admits a Tucker decomposition whose factor matrices are the stacked camera matrices, and which thus has a multilinear rank of (4, 4, 4, 4) independent of n. We develop the first synchronization algorithm for quadrifocal tensors, using Tucker decomposition, alternating direction method of multipliers, and iteratively reweighted least squares. We further establish relationships between the block quadrifocal, trifocal, and bifocal tensors, and introduce an algorithm that jointly synchronizes these three entities. Numerical experiments demonstrate the effectiveness of our methods on modern datasets, indicating the potential and importance of using higher-order information in synchronization.
Group synchronization is the problem of determining reliable global estimates from noisy local measurements on networks. The typical task for group synchronization is to assign elements of a group to the nodes of a graph in a way that respects group elements given on the edges which encode information about local pairwise relationships between the nodes. In this paper, we introduce a novel higher-order group synchronization problem which is formulated on a hypergraph and seeks to synchronize higher-order local measurements along hyperedges in order to obtain global estimates on the nodes. Higher-order group synchronization is motivated by applications to computer vision, image processing, robotics, and community detection, among other domains. First, we define the problem of higher-order group synchronization and discuss its mathematical foundations. Specifically, we give necessary and sufficient conditions for synchronizability which establish the importance of cycle consistency in higher-order group synchronization. Next, we propose the first computational framework for general higher-order group synchronization; it acts globally and directly on higher-order measurements via a message passing algorithm. We derive theoretical guarantees for our framework, including convergence analyses in the presence of outliers and noise. Finally, we demonstrate potential advantages of our method through numerical experiments. In particular, we show that in certain cases our higher-order method applied to rotational and angular synchronization outperforms standard pairwise synchronization methods and is more robust to outliers. We also show that our method has comparable performance on simulated cryo-electron microscopy (cryo-EM) data compared to a standard cryo-EM reconstruction package.
Spectral embedding methods are widely used for dimensionality reduction and clustering of high-dimensional datasets with intrinsic low-dimensional structures. Although many datasets of practical interest exhibit invariance under symmetries such as rotations, standard spectral embedding methods do not account for this, treating symmetry-related data points as unrelated. Our approach to this problem is to incorporate the symmetries directly into the affinity kernels used for spectral embedding. We analyze the case of a Riemannian data manifold M with symmetries given by a compact Lie group G and prove that, under suitable conditions, graph Laplacians constructed from three types of invariant kernels converge pointwise to explicit second-order differential operators on the quotient space M/G. Our analysis implies improved convergence rates, as the effective dimension drops according to the dimension of the group. We validate our approach on datasets with SO(2) or SO(3) symmetry, and show that G-invariant spectral embedding recovers the intrinsic geometry of the data, in contrast to standard spectral embedding, which fails to do so even in the limit of infinite data.
In this paper, we introduce a general framework for analyzing the numerical conditioning of minimal problems in multiple view geometry, using tools from computational algebra and Riemannian geometry. Special motivation comes from the fact that relative pose estimation, based on standard 5-point or 7-point Random Sample Consensus (RANSAC) algorithms, can fail even when no outliers are present and there is enough data to support a hypothesis. We argue that these cases arise due to the intrinsic instability of the 5- and 7-point minimal problems. We apply our framework to characterize the instabilities, both in terms of the world scenes that lead to infinite condition number, and directly in terms of ill-conditioned image data. The approach produces computational tests for assessing the condition number before solving the minimal problem. Lastly, synthetic and real data experiments suggest that RANSAC serves not only to remove outliers, but in practice it also selects for well-conditioned image data, which is consistent with our theory.
We study shallow neural networks with monomial activations and output dimension one. The function space for these models can be identified with a set of symmetric tensors with bounded rank. We describe general features of these networks, focusing on the relationship between width and optimization. We then consider teacher-student problems, which can be viewed as problems of low-rank tensor approximation with respect to non-standard inner products that are induced by the data distribution. In this setting, we introduce a teacher-metric data discriminant which encodes the qualitative behavior of the optimization as a function of the training data distribution. Finally, we focus on networks with quadratic activations, presenting an in-depth analysis of the optimization landscape. In particular, we present a variation of the Eckart-Young Theorem characterizing all critical points and their Hessian signatures for teacher-student problems with quadratic networks and Gaussian training data.
Cryo-electron microscopy is a powerful imaging technique for reconstructing three-dimensional molecular structures from noisy tomographic projection images of randomly oriented particles. We introduce a new data fusion framework, termed the method of double moments, which reconstructs molecular structures from two instances of the second-order moment of projection images obtained under distinct orientation distributions: one uniform, the other non-uniform and unknown. We prove that these moments generically uniquely determine the underlying structure, up to a global rotation and reflection, and we develop a convex-relaxation-based algorithm that achieves accurate recovery using only second-order statistics. Our results demonstrate the advantage of collecting and modeling multiple datasets under different experimental conditions, illustrating that leveraging dataset diversity can substantially enhance reconstruction quality in computational imaging tasks.
Motivated by a flurry of recent work on efficient tensor decomposition algorithms, we show that the celebrated moment matrix extension algorithm of Brachat, Comon, Mourrain, and Tsigaridas for symmetric tensor canonical polyadic (CP) decomposition can be made efficient under the right conditions. We first show that the crucial property determining the algorithm's complexity is the regularity of a target decomposition. This observation allows us to reduce the overall complexity of the vanilla algorithm, while also unifying results from existing literature. We introduce the notion of efficient formats: tuples (n,d,r) such that generic tensors in SdCn+1 of rank r are efficiently decomposable via a modification of the extension algorithm – solving a linear system; moreover, this implies a new criterion for determining uniqueness of decompositions. We prove that (n,4,r) is an efficient format for r≤2n+1, implying that the modified extension algorithm is efficient for ranks exceeding those allowed by simultaneous diagonalization. We then formulate a conjecture that (n,4,r) is an efficient format for r=O(n2), implying that the modified extension algorithm is efficient for ranks asymptotically matching existing algorithms, with improved leading coefficient. Towards an affirmative answer, we give computer-assisted proofs verifying correctness for n=2,…,17. Next, we demonstrate that some classes of tensors with nonunique decompositions can be decomposed efficiently via the modified extension algorithm, a first for efficient algorithms. Of particular interest is the class of monomials, for which the extension algorithm is not only efficient but also improves on existing theory by explicitly parameterizing the space of decompositions. Code for implementations of the efficient algorithm for generic tensors and monomials is provided, along with several numerical examples.
Principal component analysis (PCA) is a fundamental technique for dimensionality reduction and denoising; however, its application to three-dimensional data with arbitrary orientations – common in structural biology – presents significant challenges. A naive approach requires augmenting the dataset with many rotated copies of each sample, incurring prohibitive computational costs. In this paper, we extend PCA to 3D volumetric datasets with unknown orientations by developing an efficient and principled framework for SO(3)-invariant PCA that implicitly accounts for all rotations without explicit data augmentation. By exploiting underlying algebraic structure, we demonstrate that the computation involves only the square root of the total number of covariance entries, resulting in a substantial reduction in complexity. We validate the method on real-world molecular datasets, demonstrating its effectiveness and opening up new possibilities for large-scale, high-dimensional reconstruction problems.
Principal component analysis (PCA) is a tool to capture factors that explain variation in data. Across domains, data are now collected across multiple contexts (for example, individuals with different diseases, cells of different types, or words across texts). While the factors explaining variation in data are undoubtedly shared across subsets of contexts, no tools currently exist to systematically recover such factors. We develop multi-context principal component analysis (MCPCA), a theoretical and algorithmic framework that decomposes data into factors shared across subsets of contexts. Applied to gene expression, MCPCA reveals axes of variation shared across subsets of cancer types and an axis whose variability in tumor cells, but not mean, is associated with lung cancer progression. Applied to contextualized word embeddings from language models, MCPCA maps stages of a debate on human nature, revealing a discussion between science and fiction over decades. These axes are not found by combining data across contexts or by restricting to individual contexts. MCPCA is a principled generalization of PCA to address the challenge of understanding factors underlying data across contexts.
We devise fast and provably accurate algorithms to transform between an N × N × N Cartesian voxel representation of a three-dimensional function and its expansion into the ball harmonics, that is, the eigenbasis of the Dirichlet Laplacian on the unit ball in R 3 . Given ε > 0 , our algorithms achieve relative ℓ 1 - ℓ ∞ accuracy ε in time 𝒪 ( N 3 ( l o g N ) 2 + N 3 | l o g ε | 2 ) , while the naive direct application of the expansion operators has time complexity 𝒪 N 6 . We illustrate our methods on numerical examples.
This educational article highlights the geometric and algebraic complexities that distinguish tensors from matrices, to supplement coverage in advanced courses on linear algebra, matrix analysis, and tensor decompositions. Using the case of real-valued 2 × 2 × 2 tensors, we show how tensors violate many well-known properties of matrices: (1) The rank of a matrix is bounded by its smallest dimension, but a 2 × 2 × 2 tensor can be rank 3. (2) Matrices have a single typical rank, but the rank of a generic 2 × 2 × 2 tensor can be 2 or 3 - it has two typical ranks. (3) Any limit point of a sequence of matrices of rank r is at most rank r, but a limit point of a sequence of 2 × 2 × 2 tensors of rank 2 can be rank 3 (a higher rank). (4) Matrices always have a best rank-r approximation, but no rank-3 tensor of size 2 × 2 × 2 has a best rank-2 approximation. We unify the analysis of the matrix and tensor cases using tools from algebraic geometry and optimization, providing derivations of these surprising facts. To build intuition for the geometry of rank-constrained sets, students and educators can explore the geometry of matrix and tensor ranks via our interactive visualization tool.
Covering numbers are a powerful tool used in the development of approximation algorithms, randomized dimension reduction methods, smoothed complexity analysis, and others. In this paper we prove upper bounds on the covering number of numerous sets in Euclidean space, namely real algebraic varieties, images of polynomial maps and semialgebraic sets in terms of the number of variables and degrees of the polynomials involved. The bounds remarkably improve the best known general bound by Yosef Yomdin and Georges Comte (Tame Geometry with Application in Smooth Analysis, 2004), and our proof is much more straightforward. In particular, our result gives new bounds on the volume of the tubular neighborhood of the image of a polynomial map and a semialgebraic set, where results for varieties by Martin Lotz (In: Proceedings of the American Mathematical Society, 143.5:1875–1889, 2015) and Saugata Basu and Antonio Lerario (In: Mathematische Annalen, 1–31, 2022) are not directly applicable. We illustrate the power of the result on three computational applications. Firstly, we derive a near-optimal bound on the covering number of tensors with low canonical polyadic (CP) rank, quantifying their approximation properties and filling in an important missing piece of theory for tensor dimension reduction and reconstruction. Secondly, we prove a bound on dimensionality reduction of images of polynomial maps via randomized sketching, which has direct applications to large scale polynomial optimization. Finally, we deduce generalization error bounds for deep neural networks with rational or ReLU activation functions, improving or matching the best known results in the machine learning literature, while helping to quantify the impact of architecture choice on generalization error.
We present an algorithm for low rank decomposition of tensors of any symmetry type, from fully asymmetric to fully symmetric. It recovers the decomposition one summand at a time via the higher-order power method. This approach is known to fail in general: there need not be a relationship between the summands of a decomposition and the (partially symmetric) singular vector tuples (pSVTs) of the tensor. Our approach overcomes this problem by transforming the input to a tensor with orthonormal slices, via orthogonalization of a flattening. The summands of the decomposition of the original tensor can be recovered from the pSVTs of this new transformed tensor. We introduce a shifted power method for computing pSVTs and prove its global convergence. Numerical experiments demonstrate that our algorithm achieves higher accuracy and faster runtime than existing methods.
We introduce the Subspace Power Method (SPM) for calculating the CP decomposition of low-rank real symmetric tensors. This algorithm calculates one new CP component at a time, alternating between applying the shifted symmetric higher-order power method (SS-HOPM) to a certain modified tensor, constructed from a matrix flattening of the original tensor; and using appropriate deflation steps. We obtain rigorous guarantees for SPM regarding convergence and global optima for input tensors of dimension d an order m of CP rank up to 𝒪(d^⌊m/2⌋) , via results in classical algebraic geometry and optimization theory. As a by-product of our analysis we prove that SS-HOPM converges unconditionally, settling a conjecture in Kolda and Mayo (SIAM J. Matrix Anal. Appl. 32(4), 1095–1124 2011). We present numerical experiments which demonstrate that SPM is efficient and robust to noise, being up to one order of magnitude faster than state-of-the-art CP decomposition algorithms in certain experiments. Furthermore, prior knowledge of the CP rank is not required by SPM.
This paper develops the sketching (i.e., randomized dimension reduction) theory for real algebraic varieties and images of polynomial maps, including, e.g., the set of low rank tensors and tensor networks. Through the lens of norming sets, we provide a framework for controlling the sketching dimension for any sketch operator used to embed said sets, including sub-Gaussian, fast Johnson-Lindenstrauss, and tensor structured sketch operators. Leveraging norming set theory, we propose a new sketching method called the median sketch. It embeds such a set V using only 𝒪( V) tensor structured or sparse linear measurements.
We study the best low-rank Tucker decomposition of symmetric tensors. The motivating application is decomposing higher-order multivariate moments. Moment tensors have special structure and are important to various data science problems. We advocate for projected gradient descent (PGD) method and higher-order eigenvalue decomposition (HOEVD) approximation as computation schemes. Most importantly, we develop scalable adaptations of the basic PGD and HOEVD methods to decompose sample moment tensors. With the help of implicit and streaming techniques, we evade the overhead cost of building and storing the moment tensor. Such reductions make computing the Tucker decomposition realizable for large data instances in high dimensions. Numerical experiments demonstrate the efficiency of the algorithms and the applicability of moment tensor decompositions to real-world datasets. Finally we study the convergence on the Grassmannian manifold, and prove that the update sequence derived by the PGD solver achieves first- and second-order criticality.
We revisit the topic of common lines between projection images in single-particle cryo-electron microscopy (cryo-EM). We derive a novel low-rank constraint on a certain 2n × n matrix storing properly scaled basis vectors for the common lines between n projection images of one molecular conformation. Using this algebraic constraint and others, we give optimization algorithms to denoise common lines and recover the unknown 3D rotations associated with the images. As an application, we develop a clustering algorithm to partition a set of noisy images into homogeneous communities using common lines, in the case of discrete heterogeneity in cryo-EM. We demonstrate the methods on synthetic and experimental datasets.
We develop new tools to study landscapes in nonconvex optimization. Given one optimization problem, we pair it with another by smoothly parametrizing the domain. This is either for practical purposes (e.g., to use smooth optimization algorithms with good guarantees) or for theoretical purposes (e.g., to reveal that the landscape satisfies a strict saddle property). In both cases, the central question is: how do the landscapes of the two problems relate? More precisely: how do desirable points such as local minima and critical points in one problem relate to those in the other problem? A key finding in this paper is that these relations are often determined by the parametrization itself, and are almost entirely independent of the cost function. Accordingly, we introduce a general framework to study parametrizations by their effect on landscapes. The framework enables us to obtain new guarantees for an array of problems, some of which were previously treated on a case-by-case basis in the literature. Applications include: optimizing low-rank matrices and tensors through factorizations; solving semidefinite programs via the Burer–Monteiro approach; training neural networks by optimizing their weights and biases; and quotienting out symmetries.
Single-particle cryogenic electron microscopy (cryo-EM) is an imaging technique capable of recovering the high-resolution three-dimensional (3D) structure of biological macromolecules from many noisy and randomly oriented projection images. One notable approach to 3D reconstruction, known as Kam's method, relies on the moments of the two-dimensional (2D) images. Inspired by Kam's method, we introduce a rotationally invariant metric between two molecular structures, which does not require 3D alignment. Further, we introduce a metric between a stack of projection images and a molecular structure, which is invariant to rotations and reflections and does not require performing 3D reconstruction. Additionally, the latter metric does not assume a uniform distribution of viewing angles. We demonstrate the uses of the new metrics on synthetic and experimental datasets, highlighting their ability to measure structural similarity.
Single particle cryogenic electron microscopy (cryo-EM) is an imaging technique capable of recovering the high-resolution 3-D structure of biological macromolecules from many noisy and randomly oriented projection images. One notable approach to 3-D reconstruction, known as Kam's method, relies on the moments of the 2-D images. Inspired by Kam's method, we introduce a rotationally invariant metric between two molecular structures, which does not require 3-D alignment. Further, we introduce a metric between a stack of projection images and a molecular structure, which is invariant to rotations and reflections and does not require performing 3-D reconstruction. Additionally, the latter metric does not assume a uniform distribution of viewing angles. We demonstrate uses of the new metrics on synthetic and experimental datasets, highlighting their ability to measure structural similarity.