Rank‑r Tensor Decomposition

• Rank‑r tensor decomposition is a method that expresses a tensor as a sum of rank‑1 outer products, extending matrix factorization to higher-order arrays.
• It leverages canonical polyadic decomposition (CPD) and maximally square unfolding to provide tight lower bounds on tensor rank and ensure uniqueness under generic conditions.
• Practical implementations use algorithms like randomized sampling, polynomial reformulation, and spectral methods, as illustrated by a 20×20×20×20 tensor having a 400×400 unfolding with rank bounds below 400.

A rank-$r$ tensor decomposition expresses a tensor as a sum of exactly $r$ rank-1 tensors, where each rank-1 tensor is an outer product of vectors—generalizing matrix rank factorizations to higher-order arrays. The minimal such $r$ is the tensor rank, a notion that grounds the canonical polyadic decomposition (CPD) and underpins both computational algorithms and complexity-theoretic aspects of tensor analysis. Determining or certifying this minimal $r$ is challenging, as tensor rank is NP-hard to compute, and both theoretical lower bounds and efficient algorithms depend intricately on algebraic and geometric structures of tensors.

1. Definitions and Characterization of Rank-$r$ Tensor Decomposition

Let $$X\in\mathbb{R}^{I_1\times I_2\times\cdots\times I_N}$$ be an order-$N$ tensor. $X$ is said to have rank-1 if it factors as $$X = a^{(1)}\otimes a^{(2)}\otimes \cdots\otimes a^{(N)}$$, with $a^{(n)}\in\mathbb{R}^{I_n}$ for all $n$. The tensor rank $R$ of $X$ is the smallest integer $r$ such that

$$X = \sum_{k=1}^r \lambda_k\, a_k^{(1)}\otimes a_k^{(2)}\otimes\cdots\otimes a_k^{(N)},$$

for some scalars $\lambda_k$ and factor vectors $a_k^{(n)}$. The CPD (also known as CANDECOMP/PARAFAC) of $X$ of rank $R$ is this minimal representation, often written in multilinear notation as $X = \Vert D; A^{(1)}, \dotsc, A^{(N)}\Vert$, with factor matrices $A^{(n)}\in\mathbb{R}^{I_n\times R}$ and a super-diagonal core $D$ containing the $\lambda_k$ coefficients (Calvi et al., 2019).

2. Lower Bounds on Tensor Rank: Maximally Square Unfolding

A central theoretical advance is the identification of tight lower bounds on tensor rank via matrix unfoldings—particularly, the maximally square unfolding. For an $N$-way tensor, the mode-$n$ unfolding $X_{\langle n\rangle}$ rearranges $X$ into a matrix with row index $(i_1,\ldots,i_n)$ and column index $(i_{n+1},\ldots,i_N)$. The maximally square unfolding partitions the modes into sets $S_1, S_2$ so that $|\prod_{i \in S_1} I_i - \prod_{i \in S_2} I_i|$ is minimized (i.e., the resulting matrix is as close to square as possible). The main theorem is:

$$\operatorname{rank}(X^{\rho^*}_{\langle k^*\rangle}) \leq R,$$

where $\rho^*$ is the optimal permutation and $k^*$ the corresponding split point. Under mild genericity, if

$$\operatorname{rank}(X^{\rho^*}_{\langle k^*\rangle}) < \min\left\{\prod_{i\in S_1} I_i,\, \prod_{i\in S_2} I_i\right\},$$

then $$R = \operatorname{rank}(X^{\rho^*}_{\langle k^*\rangle})$$ (Calvi et al., 2019). Thus, the matrix rank of the maximally square unfolding provides a computable tight lower bound on the CPD rank, which, under generic conditions, is exact.

For practical use, this lower bound $r_{\mathrm{lb}}$ can be computed using dynamic programming or heuristics, and prevents unnecessary computation for $R < r_{\mathrm{lb}}$. For example, a $20\times 20\times 20\times 20$ tensor has maximally square unfolding of size $400\times 400$, so any CPD of rank $R\le399$ must satisfy $R\ge399$. This substantially reduces the search space for rank-incremental decomposition methods.

3. Algorithms and Computational Complexity

The computation of a rank-$r$ tensor decomposition is inherently NP-hard. However, for generic or structured tensors, a variety of algorithmic approaches have been developed:

• Polynomial System Reformulation: The CPD can be reformulated as a system of polynomial equations, with factor matrices as unknowns, and solved using tools from computational algebraic geometry. The normal-form method exploits multigraded regularity to transform the overdetermined system into a well-conditioned eigenproblem. For many formats and ranks, effective bounds on regularity guarantee that the algorithm runs in polynomial time in the input size; in the unbalanced regime, this is conjectured to hold more generally (Telen et al., 2021).
• Algebraic and Eigenvalue Methods: For tensors with additional structure (e.g., multilinear rank-$(1, L_r, L_r)$ terms), uniqueness can be guaranteed and decomposition can be computed using explicit algebraic recipes (e.g., joint block-diagonalization or minors-based algorithms) (Domanov et al., 2018).
• Spectral and Polynomial-Time Algorithms in Special Regimes: For random or "generic" tensors, particularly in the overcomplete regime (rank > dimension), sum-of-squares (SOS) and fast spectral algorithms have enabled decomposition up to rank $r \sim n^{3/2}/\mathrm{polylog}(n)$ for order-3 tensors. These methods rely on injective norm certification, spectral gap amplification via tensor networks, and polynomial-time spectral rounding (Ge et al., 2015, Ding et al., 2022).
• Leverage-Based Sampling for Large-Scale Decomposition: For very large sparse or dense tensors, computational cost can be further controlled by randomized, leverage-based sampling in the ALS subproblems. Probabilistic guarantees for solution accuracy ensure that this "sketching" approach produces approximations with error at most $1+\varepsilon$ at a fraction of the cost and memory, with sample complexity independent of the full tensor size (Larsen et al., 2020).

4. Uniqueness, Rank Detectability, and Algorithmic Barriers

A variety of uniqueness theorems delineate when the CPD of a tensor is identifiable and the corresponding decomposition is algorithmically feasible:

• Kruskal-Type and Structural Theorems: Sufficient conditions for uniqueness involve the Kruskal ranks of the factor matrices. For more general decompositions (e.g., multilinear rank-$(1, L_r, L_r)$), explicit block partitioning conditions and combinatorial sums of dimensions guarantee uniqueness and algebraic computability (Domanov et al., 2018).
• Maximally Square Unfolding in Generic Settings: Under genericity, if the maximally square unfolding is rank-deficient, then the CP rank is exactly equal to its rank. Detectability then reduces to efficient matrix rank computation, as described above (Calvi et al., 2019).
• Statistical–Computational Gaps in Average-Case Models: There is a pronounced separation between regimes where unique decomposition is theoretically possible (statistical limit $r \lesssim n^2$ for order-3 symmetric tensors) and where polynomial-time algorithms are known to succeed ($r \lesssim n^{3/2}$). Low-degree polynomial (LDP) and SOS frameworks rigorously establish that no polynomial-time algorithm can succeed beyond this threshold, unless new algorithmic paradigms are discovered (Wein, 2022).

5. Practical Implementation and Computational Recommendations

For practitioners seeking to compute rank-$r$ tensor decompositions or determine rank, the following recommendations and caveats are substantiated:

• Utilize Maximally Square Unfolding: Always compute $$r_{\mathrm{lb}} = \operatorname{rank}(X^{\rho^*}_{\langle k^*\rangle})$$ and do not attempt decomposition for $r < r_{\mathrm{lb}}$. When $$r_{\mathrm{lb}} < \min\{\prod_{i\in S_1} I_i, \prod_{i\in S_2} I_i\}$$, the exact rank is immediately determined and decomposition can proceed without further rank search (Calvi et al., 2019).
• Embed Lower Bound Tests in Iterative Routines: In alternating least squares (ALS) or other rank-incremental schemes, the maximally square unfolding test can prune fruitless low-rank trials and accelerate convergence.
• Numerical Stability and Conditioning: Algorithms that reduce the tensor to two frontal slices (matrix pencils) for eigendecomposition, such as Jennrich’s algorithm, can be numerically unstable and lose significant accuracy, especially when $n_3=2$ or the slices are nearly rank-deficient. It is preferable to avoid such reductions or to always supplement them with stable Newton or Gauss-Newton refinement (Beltrán et al., 2018). Iterative methods that operate directly on the full tensor are generally better conditioned.
• Randomized and Scalable Methods: For massive tensors, leverage-based sampling in ALS subproblems and matrix sketching techniques provide scalable, communication-efficient approximations, enabling decomposition with near-deterministic accuracy and substantial speedup (Larsen et al., 2020).

6. Limitations, Extensions, and Open Problems

Several fundamental limitations and open directions are structurally explicit:

• Noise and Approximate Decomposition: Existing exact lower bounds and uniqueness guarantees are for noiseless, exact, and generic tensors. In practice, noise is ubiquitous, requiring adaptation via singular value thresholding or robust approximate decomposition, which is not yet fully characterized (Calvi et al., 2019).
• Structured Tensors: Maximally square unfolding and related bounds do not directly exploit additional algebraic structure, such as sparsity, nonnegativity, or symmetry. Integrating such constraints remains an active research problem.
• Algorithmic Barriers at High Rank: No polynomial-time algorithm is known for decomposing random or generic tensors above the $n^{3/2}$-rank threshold for order-3, and LDP/SOS lower bounds indicate this threshold is algorithmically sharp unless new paradigms are developed (Wein, 2022). For some higher-order or structured tensors, normal-form and algebraic algorithms can provide polynomial time guarantees, but generality is limited by open conjectures in computational algebraic geometry (Telen et al., 2021).
• Selection of Unfoldings: While the optimal maximally square partitioning is an NP-hard problem, it is tractable for moderate $N$, and good heuristics perform well empirically (Calvi et al., 2019).
• Extension to Non-generic, Symmetric, and High-Order Cases: Many theoretical guarantees critically rely on genericity (drawn from continuous distributions). For deterministic or algebraically constrained tensors, practical and algorithmic guarantees remain less well understood.

In summary, the theory and practice of rank-$r$ tensor decomposition are shaped by intricate interactions between algebraic geometry, matrix analysis, complexity theory, and numerical linear algebra. Advances in lower bounds, uniqueness, and scalable algorithms continue to reduce computational cost and clarify fundamental limits; however, challenges at high rank, under noise, and for structured data remain open for future research (Calvi et al., 2019, Wein, 2022, Telen et al., 2021).