Abstract 1 Introduction 2 Proof Overview 3 Preliminaries 4 The decomposition algorithm 5 Decomposing high dimensional tensors References Appendix A Reshaping into 3-tensors

An Algorithmic Proof of Kruskal’s Tensor Decomposition Theorem

Vishwas Bhargava ORCID Department of Computer Science and Engineering, IIT Bombay, Mumbai, India    Leonard J. Schulman ORCID Department of Computing and Mathematical Sciences, Caltech, Pasadena, CA, USA    Shiri Sivan ORCID Department of Computing and Mathematical Sciences, Caltech, Pasadena, CA, USA
Abstract

A famous theorem of Kruskal gives the simplest and arguably most fundamental criterion under which a tensor is guaranteed a unique minimum-rank decomposition. Kruskal’s condition requires that the sum of the Kruskal ranks {ki}i=1m of the components satisfies ∑i∈[m]ki≥2⁢r+m−1, where r denotes the rank and m the order of the tensor. However, Kruskal’s original proof and subsequent simplifications/generalizations have remained non-constructive. With the sole exception of the case (k1=r,k2=r,k3=2), attributed to Jennrich – no algorithm has been established for decomposing tensors under the Kruskal condition without additional assumptions. In fact, whether there exists an efficient algorithm for decomposing a tensor under the Kruskal condition was explicitly posed as an open problem in the work of Bhaskara et al. (COLT 2014). Even slight variations of the Jennrich special case, such as the (r,r−1,3) case, have remained algorithmically open; specifically, no sub-exponential time bound was known.

In this work, we make progress on this problem by giving an elementary, constructive proof of Kruskal’s Theorem for general m-way tensors. Concretely, we give a randomized algorithm that decomposes any tensor satisfying the Kruskal condition by utilizing random projections to map the problem into a geometry of intersecting hyperplanes via a MinRank instance. Specifically for 3-way tensors satisfying k1+k2+k3=2⁢r+2, the algorithm achieves a runtime of nO⁢(k) where k=min⁡(k1,k2,k3). Thus, we extend smoothly beyond the Jennrich special case, achieving polynomial-time complexity for any family of tensors that satisfies the Kruskal condition, provided the least Kruskal rank is bounded.

Keywords and phrases:
Tensor decomposition, Kruskal’s theorem, tensor rank, MinRank, algebraic algorithms
Category:
Track A: Algorithms, Complexity and Games
Funding:
Vishwas Bhargava: Part of this work was done when the author was a postdoc at Caltech, supported by Prof. Christopher Umans’ Simons Foundation award.
Leonard J. Schulman: Supported in part by NSF CCF-2321079.
Shiri Sivan: Supported in part by NSF CCF-2321079.
Copyright and License:
[Uncaptioned image] © Vishwas Bhargava, Leonard J. Schulman, and Shiri Sivan; licensed under Creative Commons License CC-BY 4.0
2012 ACM Subject Classification:
Theory of computation → Algebraic complexity theory
; Theory of computation → Design and analysis of algorithms
Related Version:
Full Version: https://eccc.weizmann.ac.il/report/2026/073 [3]
Editors:
Sayan Bhattacharya, Danupon Nanongkai, Michael Benedikt, and Gabriele Puppis

1 Introduction

Tensor decomposition is a fundamental tool for modern science, with applications ranging from machine learning and statistics to signal processing and computational complexity. In a nutshell, tensors are multi-dimensional arrays with entries from a field 𝔽. For instance, a 3-dimensional (or 3-way or 3-mode) tensor can be written as 𝒯=(αi,j,k)∈𝔽n1×n2×n3.

A tensor T∈𝔽n1×⋯×nm can always be written as the sum of rank 1 tensors; such an expression is known as a tensor decomposition.111Also known as the CANDECOMP/PARAFAC (CP) decomposition. Here, a rank 1 tensor over a field 𝔽 is the outer product of m nonzero vectors, written as

v1⊗v2⊗…⊗vm∈𝔽n1×n2×⋯×nm.

If the decomposition consists of r summands, it is called an r-decomposition. The minimal r for which there exists an r-decomposition is defined as the tensor rank of T, and such a decomposition is then referred to as a rank decomposition.

When T has a decomposition

T=∑t=1rvt(1)⊗vt(2)⊗⋯⊗vt(m),

we define the components as, for each i∈[m], the matrix with columns vt(i) (t=1,…,r):

V(i)=[v1(i),v2(i),…,vr(i)]∈ℝni×r

The vector space containing vt(i) will be referred to as the i’th mode of the tensor. We also employ the shorthand notation

T=⟦V(1),V(2),…,V(m)⟧.

A primary reason for the wide applicability of tensors is that their rank decomposition can be unique (also known as identifiable), a property that matrices (unless of rank 1) lack. The motivations for identifiability in data analysis (beginning at least with work of the psychologist Spearman on tests of intelligence [28]) are abundantly discussed in many references e.g., [6, 12, 24]. We refer the interested reader to these works and the references therein.

The most celebrated criterion for uniqueness of tensor decomposition was given by Kruskal in 1977 [19] (and see [17, 18]). The Kruskal rank (or k-rank) of a matrix is the maximum integer k such that every subset of k columns is linearly independent. Clearly, k-rank ≤ rank. Kruskal’s uniqueness theorem states that if a tensor T has an r-decomposition ⟦V(1),…,V(m)⟧, and the k-ranks of the components satisfy

∑i=1mkV(i)≥2⁢r+m−1, (1)

then rank⁡(T)=r and the decomposition is unique up to permutation and scaling of the columns (going forward we omit the “up to scaling and permutation”). The inequality in (1) is the Kruskal Condition.

While the theorem guarantees that a unique decomposition exists, the original proof and subsequent simplifications (e.g., [29, 25, 21]) are existential rather than constructive. This limitation is widely recognized. The question of whether there is a hardness result for tensor decomposition, or conversely, whether there exists an efficient (polynomial-time) algorithm for decomposing a tensor under the general conditions of Kruskal’s uniqueness theorem – and whether there is an algorithmic proof of the theorem itself – was explicitly posed as a major open problem in COLT 2014 by Bhaskara et al. [5], carrying a $100 reward.

It is important to note that while the tensor rank problem (and therefore also the problem of finding a decomposition of specified rank) is NP-hard, it remains unknown whether this hardness persists when the input tensor is promised to satisfy the Kruskal condition. The only significant special case of the Kruskal condition in which an efficient algorithm is known is the case that is due (prior to Kruskal’s work) to Jennrich, in which two of the components are invertible (k1=k2=r) while the third has Kruskal rank at least 2. Even this restricted case has found surprisingly many applications in computer science and statistics. However, even a slight variation, such as the (r,r−1,3) case, has, to the best of our knowledge, remained algorithmically open. This was the original motivation for our work.

1.1 Our contribution

We give the first algorithmic proof of the Kruskal theorem. Rather than relying on an existential argument, we use structural implications of the Kruskal condition to recover the components directly. For ease of presentation, we state our results over the real field ℝ. Our algorithm and its runtime guarantee is stated in the Real RAM model, where arithmetic over the reals is exact and unit cost. However, the algorithm can be easily implemented over the rationals using standard Turing machines. That is, if the unknown decomposition is in fact rational, then our algorithm can work entirely over the rationals. There is a known polynomial-time reduction (in the input size ∏ni) that converts an m-way tensor satisfying the Kruskal condition into a 3-way tensor that also satisfies the Kruskal condition (see Theorem 3 in [27] or Lemma 24 in the Appendix).222Lemma 24 proves that every Kruskal m-mode tensor has a Kruskal 3-reshaping. Polynomial-time follows from the fact that the number of ways to reshape the tensor is polynomial in the input size. Consequently, we focus our analysis on the 3-way case.

As our objective is to design a decomposition algorithm for Kruskal tensors, we assume T satisfies the Kruskal condition:

k1+k2+k3≥2⁢r+2.

We set k:=2⁢r+2−max⁡{k1+k2,k1+k3,k2+k3}. This value governs both the time complexity and the core technical arguments of our algorithm.

Theorem 1 (Algorithmic Kruskal for 3-tensors).

Let 𝒯∈ℝn1×n2×n3 be a tensor with a decomposition 𝒯=∑i=1rui⊗vi⊗wi such that the components have Kruskal ranks respectively k1,k2,k3 satisfying the Kruskal condition k1+k2+k3≥2⁢r+2. Then there exists a randomized algorithm that recovers the components in (n1+n2+n3)O⁢(k) time with high probability.

By definition, k≤min⁡{k1,k2,k3}. If the Kruskal condition holds with some slack (i.e., ∑ki>2⁢r+2), we can effectively replace the smallest Kruskal rank with this smaller value k that saturates the bound. For example, if all Kruskal ranks are r, then k=2.

Our work settles the request of Bhaskara et al. [5] for an algorithmic proof of the Kruskal theorem, and makes substantial progress on the efficiency question without fully settling it: we provide a time complexity upper bound for decomposing Kruskal tensors that scales exponentially with respect to the minimum Kruskal rank of the modes. Notably, for (previously open) instances “close” to the Jennrich case – such as the (r,r−p,p+2) regime for fixed p – our algorithm is polynomial in the input size. By contrast, the best existing method, the brute-force approach of reducing the decomposition to an algebraic system, will typically require exp⁡(n⁢r) time for n=max⁡ni due to the n⁢r unknowns involved.

We should mention that in the m-to-3 reduction alluded to above, the initial m-way tensor does not strictly need to satisfy the Kruskal condition; it is sufficient that it should reduce to a 3-way instance which does. The reduction is achieved by “clubbing” various modes together through Khatri-Rao products – a standard reshaping technique in the literature. (For the m-to-3 reduction see [27, 22, 8], or Sections 5 and A.) By partitioning the modes into sets I,J, and ⁢L such that the resulting reshaped tensor satisfies the 3-way Kruskal criteria, we effectively provide a constructive proof for the more general Reshaped Kruskal Theorem.

Theorem 2 (Reshaped Kruskal Theorem).

Consider an m-way tensor 𝒯 with an r-rank decomposition ⟦V(1),…,V(m)⟧. If there exists a partition of the modes I⊔J⊔L=[m] such that the Kruskal ranks of the partitioned (Khatri-Rao) matrices satisfy:

kI+kJ+kL≥2⁢r+2

then rank⁡(𝒯)=r and the decomposition is unique. Furthermore, there exists a randomized algorithm that recovers the components in NO⁢(k) time with high probability, where N is the input size, and k:=2⁢r+2−max⁡{kI+kJ,kI+kL,kJ+kL}.

1.2 Related Work

Tensor decomposition is a highly active area of research with many fundamental open problems and landmark results. Here, we mention some of these results with a focus on the identifiability and computational complexity of the problem (and omitting entirely the extensive literature in algebraic complexity theory).

On identifiability, the work of Lovitz and Petrov [22] provides an elegant generalization of Kruskal’s theorem using ideas from matroid theory. It is not clear whether our method can be extended to their generalization.

On the complexity and algorithmic side, determining tensor rank is famously NP-hard [13, 14]. Consequently, one cannot expect efficient algorithms to determine rank or find decompositions for general input tensors. The Kruskal framework is the least restrictive in which hardness is not known. The Kruskal condition is also tight [10]. From this point of view our result is a worst-case tensor decomposition algorithm, unlike “generic” uniqueness or decomposition theorems. Indeed, most research on tensor decomposition algorithms has focused on tensors with additional non-degeneracy conditions; tensors drawn from specific distributions; or heuristic approaches. We discuss some of these results now, moving from no assumptions to progressively adding more structure; as we will see, the assumptions have a tremendous effect on the rank parameters we can efficiently decompose.

In the worst-case setting for 3-dimensional tensors, one essentially cannot outperform brute force, which involves solving systems over the base field. Formally, Schaefer and Stefankovic [26] showed that for any field 𝔽, given a system S of algebraic equations over 𝔽, we can in polynomial time construct a 3-dimensional tensor 𝒯S (with a linear blowup) and an integer k such that S has a solution in 𝔽 if and only if 𝒯S has rank at most k over 𝔽. However, for higher m-dimensional tensors, there exists an FPT-style worst-case algorithm with 2rO⁢(1)⁢poly⁢(n,m) runtime [4]. The next level is the generic setting, where the input tensor can be anything barring a (result-specific) lower-dimensional algebraic variety. For generic tensors, Chiantini and Ottaviani [7] prove the uniqueness of decomposition up to rank (n1−1)⁢(n2−1)⁢(n3−1)+1n1+n2+n3−2. Again in the generic setting, Kothari et al. [16] recently provided an efficient decomposition algorithm for families of tensors with rank approaching 2⁢n. Adding even more structure (and moving towards applications), one has the average-case setting; for instance, Ma et al. [23] deployed techniques from the Sum-of-Squares hierarchy to reach rank O~⁢(n1.5) for tensors drawn from Gaussian decomposition.

2 Proof Overview

We focus our analysis on 3-way tensors, as decomposition of m-way tensors can be reduced in a standard fashion to the 3-way case through the reshaping procedure alluded to in Section 1.1. Let T∈ℝn1×n2×n3 be a tensor that satisfies the Kruskal condition k1+k2+k3≥2⁢r+2, and has a decomposition T=∑j=1ruj⊗vj⊗wj. Since part of our goal is to reprove the Kruskal theorem, we shall not assume that this decomposition is unique. As per our notation, we define components U=[u1,…,ur]∈ℝn1×r, V=[v1,…,vr]∈ℝn2×r, and W=[w1,…,wr]∈ℝn3×r.

We start by recovering a single component, say V. Towards this end, we turn to a geometric point of view. Namely, we complete an equivalent task: recovering the r hyperplanes orthogonal to the columns of V. This task is completed in two steps. First we identify all the intersection points of the hyperplanes (collectively these form a zero-dimensional variety); then we reconstruct the hyperplane arrangement, given the intersection points. We refer to the combination of these two steps as the identifiability test, as this test establishes the identifiability of the (projected) component V.

Both steps require that the hyperplane arrangement lie in general position, meaning the intersection of any t hyperplanes has dimension max⁡{d−t,0}. Conveniently, this can be accomplished with probability 1 via a random projection (see Lemma 13).

The identifiability test

To identify the intersection points, we construct a MinRank instance333The MinRank search problem is the task of finding a linear combination of given matrices M1,…,Mn that minimizes the rank of the resulting matrix. See e.g., [11, 1] whose solutions correspond precisely to the aforementioned intersection points. Let T′=T×2A=⟦U,A⁢V,W⟧ be the projected tensor, where A∈ℝk2×n2 is a random matrix (×iA means that the i-th mode of the tensor is contracted against the rows of matrix A). For variables z=(z1,…,zd) where d=k2−1, and [z,1] is the vector (z1,…,zd,1), we define the V-mode contraction:

M⁢(z)=T′×2[z,1]=U⁢Diag⁢([z,1]⁢A⁢V)⁢WT. (2)

In words, M⁢(z) is the linear combination of 2nd mode slices associated with coefficients [z,1]. Under the Kruskal condition, we show that for all z∈ℝd, rank⁢(M⁢(z))≥r−d. Crucially, we prove that the matrix M⁢(z) achieves this minimum rank r−d if and only if [z,1]⁢A⁢V has exactly d zeros.

We only look at slice combinations of the form [z,1] because fixing the last entry ensures that d-intersections are 0-dimensional (points). A more subtle reason for this choice is the projective nature of tensor decomposition; since we can scale components freely, we have the freedom to look for one fewer variable in the search for these hyperplanes.

Geometrically, this implies that the minimizers of rank ⁢M⁢(z) are precisely the d-wise intersection points of the hyperplanes defined by the columns of the projected component A⁢V (see Theorem 14). This already is a significant step toward identifiability, as the matrix M⁢(z) depends solely on the observable tensor rather than the hidden decomposition itself. We are essentially reducing tensor decomposition to the search problem associated with MinRank. We denote

MinRank⁢(M⁢(z),ℓ)=minz∈ℝℓ−1⁡rank⁢(M⁢(z)), (3)

where M⁢(y):=W⁢Diag⁢([y,1]⁢A⁢V)⁢WT for a random A∈ℝℓ×nv.444The parameter ℓ is the number of matrices in the MinRank instance, since M⁢(z) is a sum of the matrices {W⁢Diag⁢(ai∗⁢V)⁢WT}i∈[ℓ], where ai∗ denotes the i-th row of A.

It is a geometric fact that the set of d-wise intersection points of an arrangement of hyperplanes in general position in ℝd, determine those hyperplanes; and in addition, the hyperplanes can be efficiently computed from the intersection points. This allows us to reconstruct the matrix A⁢V (see Algorithm HP-reconstruction and Lemma 17).

The identifiability test is the heart of our algorithmic proof of Kruskal’s Theorem. Establishing the uniqueness of the projected component A⁢V via this test effectively reduces the problem to the Jennrich setting, where uniqueness is well-understood and efficient decomposition applies (see Lemma 18). With U and W at our disposal, the identifiability of the original V follows.

The tensor decomposition algorithm

Our decomposition algorithm, sketched below, is directly based on the ideas we discussed above for identifiability. For convenience, at each step we link to the text for more details.

Algorithm 1 Rank decomposition algorithm for 3-way tensors (outline).

Input: T∈ℝn1×n2×n3.

For k′=2,3,…,min⁡{n1,n2,n3}:

  • For i=1,2,3:

    1. 1.

      Apply a random projection A∈ℝk′×ni on the i-th mode (Step 1).

    2. 2.

      Generate M⁢(z)= the formal linear combination of mode-i slices with coefficients (z1,…,zk′−1,1) (Step 2).

    3. 3.

      Solve MinRank⁢(M⁢(z),k′). (Step 3)

      1. (a)

        Denote the solution m. (see Equation 3 for the definition of MinRank)

      2. (b)

        Calculate the minimizers 𝒮={z∈ℝk′−1|rank⁡M⁢(z)=m}.

      3. (c)

        Set r′=m+k′−1.

    4. 4.

      Reconstruct A⁢V(i) (Step 4).

    5. 5.

      Restore the two other components V(j),V(ℓ) (Step 5).

    6. 6.

      Recover V(i) (Step 6).

    7. 7.

      If T=⟦V(1),V(2),V(3)⟧, output the decomposition.

Output fail.

We employ randomness in Steps 1, 3 and 4. The successful event in Step 1 occurs with probability 1 (see Lemma 25), whereas the randomness in Steps 3 and 4 is utilized solely in sake of computational efficiency and does not affect the uniqueness argument. Concretely, each randomized step admits a deterministic alternative: the random projection can be replaced by going over a sufficiently large grid (by Schwartz-Zippel); the system solving step (and thus MinRank) can be performed deterministically (for instance using the algorithm of Ierardi [15]); and the hyperplane reconstruction can be performed via the deterministic brute-force strategy described in Section 4.2 (specifically, Lemma 16 and the subsequent discussion).555The deterministic version can be implemented in rO⁢(d2) time, compared to rO⁢(d) for the randomized version. Thus the algorithm can in principle be made fully deterministic, and its output depends only on the input tensor – not on any random choices.

Thus, the algorithm’s output is uniquely determined for each i and k′. It only remains to show that different parameter values cannot result in distinct decompositions. Towards this end, we prove the existence of fixed parameters k and i, determined solely by T, for which the recovery algorithm is guaranteed to succeed. Specifically, if T admits a rank decomposition satisfying the Kruskal condition, the algorithm – called on these parameters – is guaranteed to recover the constituent components. Because the algorithm is deterministic, it maps the tensor T to a unique set of factors (up to scaling and permutation). If two distinct decompositions were to exist, both would necessarily be recovered as the same output by this deterministic procedure, which is a contradiction. Consequently, the decomposition is unique.

This parameter k is defined by

k:=2⁢r+2−max⁡{k1+k2,k1+k3,k2+k3},

and the corresponding i is the mode of minimal Kruskal rank. That is, k≤ki=min⁡{k1,k2,k3}. 666(k,i) are not the only parameters on which the algorithm is guaranteed to succeed. However, for uniqueness we only need one set of parameters that succeed on all decompositions. Computationally, the choice (k,i) delivers the best time complexity.

Since we do not know the tensor rank and Kruskal ranks we iterate over different k guesses in the outermost loop of the algorithm. The correctness of the algorithm, as well as the uniqueness of the decomposition, rely on the following facts:

  1. 1.

    If an r-decomposition satisfies the Kruskal condition, then the tensor has rank r and every rank decomposition has the same Kruskal ranks and thus k (see Lemma 19).

  2. 2.

    At iterations k′<k, the decomposition algorithm does not output a “redundant” decomposition, i.e., an r′-decomposition with r′>r (see Lemma 20).

  3. 3.

    At iteration k, the algorithm outputs a (unique) r-decomposition (see Theorem 21).

The proofs of facts 1 and 2 are provided in Section 4, while an overview of the 3rd fact follows below. For technical details regarding the subroutines, time complexity, and success probability, see Section 4.

For fact 3, we now “walk” through iteration (k,i) of the algorithm, assuming the success event occurs at every step of the algorithm.

Step 1. Apply a random projection

A random projection will ensure that A⁢V(i) (for i∈[3]) has full Kruskal rank with high probability. The success event implies that the hyperplanes perpendicular to the columns of Q:=A⁢V(i) lie in general position. More formally, denote the columns of Q by q1,…,qr. For all i∈[r] we define the hyperplane

Hi={z∈ℝd|⟨[z,1],qi⟩=0}.

We denote the full set of hyperplanes by ℋ={H1,…,Hr}. By Lemma 13, ℋ lies in general position with probability 1.

Step 2. Generate a MinRank instance

The MinRank instance generated at this step is the linear combination of i-mode slices associated with coefficients (z1,…,zd,1), defined by

M⁢(z)=T×i[z,1]⁢A=V(j)⁢Diag⁢([z,1]⁢A⁢V(i))⁢(V(ℓ))T.

Step 3. Solve the MinRank instance

The above specific factorization of M⁢(z) coupled with the Kruskal structure of the factor matrices V(j) and V(ℓ) ensures that the rank drops precisely when z lies on the intersection of the hyperplanes (see Theorem 14). Thus we have,

MinRank⁢(M⁢(z),k)=r−d

so at Step 3c we set r′=r, and then at 3b we find all the minimizers

𝒮={z∈ℝd|rank⁡M⁢(z)=r−d}.

By Theorem 14 this collection of points is exactly the set of d-wise intersection points of hyperplanes in ℋ.

Step 4. Reconstruct the projected component

Since ℋ is in general position in ℝd, all d of the hyperplanes intersect at a unique point, and each hyperplane contains exactly (r−1d−1) points from 𝒮. Moreover, any hyperplane in ℝd\ℋ contains fewer points from 𝒮 (see Lemma 27) so we can reconstruct ℋ by iterating over subsets of d points, generating a list of candidate hyperplanes, and checking which ones contain the characterizing number of points (see Algorithm HP-reconstruction for a randomized version of this algorithm).

Recovering ℋ is equivalent to the recovery of Q=A⁢V(i) (up to scaling and permutation of the columns), because by definition, the columns of Q are the normals to the hyperplanes.

Step 5. Restore the remaining components 𝑽(𝒋) and 𝑽(ℓ)

With the projected component Q:=A⁢V(i) in hand we recover the components V(j) and V(ℓ) by essentially reducing the problem to the “Jennrich setting,” and applying the Jennrich decomposition algorithm (see Lemma 9).

Specifically, we find a matrix P∈ℝ2×k perpendicular to the first k−2 columns of Q, and effectively “zero out” those summands in the decomposition. The resulting tensor now has a reduced number of components, r′=r−k+2. Furthermore, this reduction transforms the tensor from an initial (k1,k2,k3) configuration to a “Jennrich state” where the Kruskal ranks of the j and ℓ modes are r′, while the projected i-mode’s Kruskal rank becomes 2. The fact that the truncated components have full rank r′ follows from the Kruskal condition (∑ki≥2⁢r+2 implies that kj,kℓ≥r′=r−k+2). This ensures that the truncated tensor satisfies the requirement for Jennrich’s algorithm. The remaining component can be recovered directly by applying pj with ⟨pj,qi⟩=δi,j (see Algorithm 5). Thus, we can uniquely recover V(j) and V(ℓ) (in polynomial time).

Step 6. Recover 𝑽(𝒊)

We recover V(i) by solving the linear system obtained from the entries of T (which become linear in the entries of V(i) once we plug in the values of V(j) and V(ℓ)).

Generalization to 𝒎-way tensors

The transition from 3-way to general m-way tensors is handled via a 3-reshaping process. Any m-dimensional tensor satisfying the Kruskal condition can be partitioned into a 3-way tensor through Khatri-Rao products of its components.

We prove that if the original tensor satisfies the Kruskal condition, there exists a partition I⊔J⊔L=[m] such that the resulting 3-way tensor also satisfies the Kruskal condition. After decomposing this reshaped 3-way tensor to find the partitioned components (e.g., VI,VJ,VL), we recover the original components {V(i)}i∈[m] by exploiting the “injectivity” (up to scaling) of the Khatri-Rao product. This extends our “worst-case” guarantee to any m-dimensional tensor satisfying the Kruskal bound. In terms of time complexity, the search for a 3-partition satisfying the Kruskal condition contributes a multiplicative factor of O⁢(3m) which is polynomial in the input size.

2.1 Outline of the paper

Section 3 collects definitions and notation (including some given previously). In Section 4 we formally present the decomposition algorithm, analyze it, and provide auxiliary lemmas and subroutines. Section 5 contains the proof of the main theorem. Finally, the Appendix contains more technical details, some of which may be well-known.

3 Preliminaries

For matrix M∈ℝm×n, row indices R⊆[m], and column indices C⊆[n], we denote M⁢(R,C) the submatrix of M induced by R and C. We also denote M⁢(C):=M⁢([m],C). We denote the rows of M by mi∗ for i∈[m] and the columns by mj for j∈[n] (generally, rows/columns are denoted by the lower case letter of the matrix notation).

Definition 3 (Kruskal rank).

The Kruskal rank (or k-rank) of a matrix M∈ℝm×n is the largest integer k such that every subset of k columns is linearly independent.

Definition 4 (mode-i slices).

We define the mode-i slices of a m-way tensor T as the (m−1)-way tensors obtained from T by fixing the i-th index.

Definition 5 (mode-i product).

The mode-i product of a tensor T∈ℝn1×⋯×nm with a matrix A∈ℝJ×ni, denoted by T×iA, is a tensor of dimensions n1×⋯×ni−1×J×ni+1×⋯×nm. For all j∈J, the j-th slice of T×iA is a linear combination of the slices of T along the same mode:

(T×iA):,…,:,j,:,…,:=∑k=1niAj⁢k⁢T:,…,:,k,:,…,:,

Note that if T has decomposition

T=∑j=1rvj(1)⊗vj(2)⊗⋯⊗vj(m)

then

T×iA=∑j=1rvj(1)⊗vj(2)⊗⋯⊗A⁢vj(i)⊗⋯⊗vj(m).
Definition 6 (mode-V contractions).

Let T∈ℝnu×nv×nw be a 3-way tensor with decomposition T=⟦U,V,W⟧. Define for all z∈ℝnv−1,

M⁢(z)=U⁢D⁢([z,1]⁢V)⁢WT=T×v[z,1].

In other words, M⁢(z)∈ℝnu×nw is the mode-v product of T with [z,1]∈ℝnv.

Definition 7 (Khatri-Rao product).

Let matrices A=[a1,a2,…,ar]∈ℝna×r and B=[b1,b2,…,br]∈ℝnb×r. The Khatri-Rao product A⊙B is the column-wise Kronecker product:

A⊙B=[a1⊗b1,a2⊗b2,…,ar⊗br]∈ℝ(na⁢nb)×r.

Note that the rows of A⊙B are the Hadamard products ai∗⊙bj∗∈ℝr, for all i∈[na] and j∈[nb] (the Hadamard product is the element-wise product of two vectors of the same length).

Definition 8 (d-reshaping).

Let T∈ℝn1×⋯×nm be an m-way tensor and let 𝒮={S1,…,Sd} be a partition of the index set [m]. The d-reshaping of T with respect to 𝒮 is a d-way tensor T𝒮∈ℝN1×⋯×Nd, where the dimension of the i-th mode is given by Ni=∏j∈Sinj.

If T has a CP decomposition T=∑j=1rvj(1)⊗⋯⊗vj(m), then T𝒮 admits a d-way decomposition:

T𝒮=∑j=1ruj(1)⊗uj(2)⊗⋯⊗uj(d), (4)

where each factor vector uj(i) is the vectorization of the Kronecker product of the original factor vectors in the i-th partition set:

uj(i)=vec⁢(⨂k∈Sivj(k)). (5)

Equivalently, the component U(i)=[u1(i),…,ur(i)] is the Khatri-Rao product of the matrices {V(k)}k∈Si following the order of indices in Si.

Lemma 9 (Jennrich).

Let T=⟦U,V,W⟧ be a 3-way tensor of rank r, with components U∈ℝnu×r, V∈ℝnv×r, and W∈ℝnw×r. Suppose that:

  1. 1.

    U and W have full column rank r.

  2. 2.

    Every pair of columns of V is linearly independent (i.e., kv≥2).

Then the rank decomposition of T is unique. Furthermore, there exists an algorithm to recover the components in time O⁢(nu⁢nv⁢nw+r3).

Lemma 10 (Making Tensor Concise).

[2, Lemma 5.1] Let T∈ℝn1×n2×⋯×nm be an m-th order tensor of rank r. There exists a randomized algorithm with runtime poly⁢(∏ni,r) that outputs a compressed tensor T′∈ℝn1′×n2′×⋯×nm′ and transformation matrices Ai∈ℝni′×ni of full row rank such that, with probability 1:

  1. 1.

    T′=(A1⊗A2⊗⋯⊗Am)⋅T.

  2. 2.

    The dimensions of T′ are reduced such that ni′=ranki⁢(T)≤r for each mode i∈{1,…,m}, where r⁢a⁢n⁢ki⁢(T) refers to the rank of the 2-reshaping associated with the partition [m]={i}⊔[m]\{i}.

  3. 3.

    The rank r and the Kruskal ranks {ki}i=1m are preserved; consequently, T satisfies the Kruskal condition if and only if T′ does.

The factor matrices {Ui}i=1m of the original tensor T are recovered from the factor matrices {Ui′}i=1m of T′ via an efficiently computable linear mapping, requiring at most poly⁢(∏ni,r) time.

Theorem 11 (Lakshman–Lazard [20]).

Let I⊂K⁢[x1,…,xn] be a zero-dimensional ideal over a field K of characteristic zero, with generators of degree at most d. Under any admissible monomial ordering, the Gröbner basis of radical of I can be computed in dO⁢(n) time, with probability at least 1−1/2dn. Moreover, zero-dimensionality of I can be decided at the same cost and with the same probability of success.

Lemma 12 (Shape Lemma [9, Ch. 2, Exercise 16]).

Let I⊂K⁢[x1,…,xn] be a zero-dimensional radical ideal such that the xn-coordinates of the points in V⁢(I) are distinct. Let G be the reduced Gröbner basis of I with respect to a lex monomial order with xn as the last variable. If |V⁢(I)|=m, then G consists of n polynomials

g1=x1+h1(xn),…,gn−1=xn−1+hn−1(xn),gn=xnm+hn(xn),

where h1,…,hn are polynomials in xn of degree at most m−1. In particular, all m points of V⁢(I) are recovered by finding the roots of gn and back-substituting into g1,…,gn−1.

4 The decomposition algorithm

In this section, we formally present the decomposition algorithm for 3-way tensors. For tensors with m>3 modes, we first reshape the tensor to a 3-way tensor, apply our decomposition algorithm, and finally recover the components of the original m-way tensor. We defer the discussion of this process to Section 5.

Let T be a 3-way tensor with an r-decomposition ⟦U,V,W⟧, and corresponding Kruskal ranks ku,kv,kw. Recall that

k:=2⁢r+2−max⁡{ku+kv,ku+kw,kv+kw}. (6)

Under the Kruskal condition, k is independent of the specific decomposition (see Lemma 19).

Because k and its associated mode are unknown a priori, the full decomposition procedure (Algorithm 1) iteratively tests candidate values k′≥2 and modes i∈{u,v,w} until it reaches k (and its corresponding mode).777To optimize performance, the search over k′ is the outermost loop, as time complexity is exponential in k. We restate a single round of this process as Algorithm DECOMPOSE, which takes parameters (T,k′,i). We show that Algorithm DECOMPOSE(T,k,i) returns the unique decomposition with high probability, provided that T satisfies the Kruskal condition and ki=min⁡{ku,kv,kw}.

Here, we will assume that our tensor is concise: A concise tensor is an m-th order tensor where each mode’s dimension ni is equal to the rank of its mode-i flattening. Consequently, ni≤r for all i∈{1,…,m}. Due to Lemma 10, this “compression” and the consequent “decompression” can be performed efficiently, so we may assume WLOG that the input tensor is given in concise form.

Algorithm DECOMPOSE single iteration of Algorithm 1.

Input: (T,k′,i), where T is concise. Denote d′=k′−1.

  1. 1.

    Apply a random matrix A∈ℝk′×ni on the i-th mode of T to obtain T′=T×iA.

  2. 2.

    For variables z=(z1,…,zd′) form the matrix

    M⁢(z):=T′×i[z,1].
  3. 3.

    Run MinRank(M,k′).

    1. (a)

      if the algorithm fails, continue to the next iteration of Algorithm 1. Otherwise,

    2. (b)

      denote m:=MinRank⁢(M,k′) and

      𝒮={z∈ℝd′|rank⁡M⁢(z)=m}
    3. (c)

      Set r′=m+d′.

  4. 4.

    Apply HP-reconstruction(𝒮,r′,k′) to recover the i-th component of T′.

    1. (a)

      If the algorithm fails, continue to the next iteration of Algorithm 1. Otherwise,

    2. (b)

      denote the output Q.

  5. 5.

    To recover the two other components U′ and W′, apply Algorithm 5 on input (T′,Q).

  6. 6.

    Solve the system obtained by the v-mode reshaping

    T{v},{u,w}=V′⁢(U′⊙W′).
  7. 7.

    Verify that T=⟦U′,V′,W′⟧ (if not, continue to the next iteration).

Output: (U′,V′,W′).

The algorithm follows the exact procedure described in the proof overview (Section 2). We begin by elaborating on each step, followed by the proof of correctness of DECOMPOSE in Theorem 21.

As discussed in the proof overview (Section 2), a key initial step to recover the factor matrices is to ensure that the hyperplanes perpendicular to the columns of the projected matrix V lie in general position. In Step 1 of the algorithm, a projection matrix A∈ℝk×nv is drawn from an absolutely continuous distribution (e.g., Gaussian) and applied to the v-mode. The resulting tensor 𝒯′∈ℝnu×k×nw possesses an r-rank decomposition ⟦U,Q,W⟧, where Q:=A⁢V∈ℝk×r represents the projected factor matrix.

This ensures that Q inherits the necessary structural properties for recovery, as formalized in the following lemma:

Lemma 13 (General Position of Projected Hyperplanes).

Let V∈ℝnv×r be a matrix with Kruskal rank kv≥k. Let A∈ℝk×nv be drawn from an absolutely continuous distribution (e.g., Gaussian), and let Q=A⁢V. Then, with probability 1:

  1. 1.

    Q has full Kruskal rank k.

  2. 2.

    The first d=k−1 rows of Q have full Kruskal rank d.

  3. 3.

    The set of hyperplanes ℋ={Hj}j=1r defined by Hj={z∈ℝd∣⟨[z,1],qj⟩=0} lies in general position in ℝd.

We defer the proof of the above to the full version [3]. If A is a “good” projection, then the set ℋ lies in general position in ℝd, meaning that the intersection of any t≤d hyperplanes forms a (d−t)-dimensional flat in ℝd (and the intersection of t>d hyperplanes is empty). In particular, d-wise intersections are points, and every point is unique. There are exactly (rd) such points (see Lemma 27).

4.1 Step 2: finding the 𝒅-wise intersection points

In Step 2 of the decomposition algorithm, we list all the points z, for which

rank⁡M⁢(z)≤r−d.

We need to show that:

  1. 1.

    These points are precisely the d-wise intersection points of hyperplanes in ℋ. Namely,

    {z:rank⁡M⁢(z)≤r−d}=⋃I:|I|=d⋂i∈IHi,

    and

  2. 2.

    there are exactly (rd) points in this set.

The first fact is proved in the following lemma, while the second fact follows from a standard property of general position; we defer its statement and proof to Appendix Lemma 25.

Theorem 14.

Let T=⟦U,Q,W⟧ be an r-decomposition satisfying the Kruskal condition. Suppose that Q∈ℝk×r has full Kruskal rank k, and let M⁢(z)=T×q[z,1] denote the Q-mode contraction for z∈ℝd where d=:k−1.

Then,

  1. 1.

    for all z∈ℝd,

    rank⁡M⁢(z)≥r−d,
  2. 2.

    rank⁡M⁢(z)=r−d if and only if Wt⁡([z,1]⁢Q)=r−d, and

  3. 3.

    The minimizers of rank⁡M are precisely the d-wise intersection points:

    {z:rank⁡M⁢(z)=r−d}=⋃I:|I|=d⋂i∈IHi.
Proof.

For all z∈ℝd, denote Iz⊆[r] the indices where [z,1]⁢Q is zero. Thus, [z,1]⁢Q⁢(Iz)=0 for all z. Since Q has Kruskal rank k, the subsystem Q⁢(I) has full rank k whenever |I|≥k, so |Iz|≤d=k−1 for all z. That is, [z,1]⁢Q has at most d zeros,

⟹Wt⁡([z,1]⁢Q)≥r−d (7)

for all z.

For Items 1 and 2 we prove the following:

  1. 1.

    Wt⁡([z,1]⁢Q)=r−d⟹rank⁡M⁢(z)=r−d,

  2. 2.

    Wt⁡([z,1]⁢Q)>r−d⟹rank⁡M⁢(z)>r−d.

Fix z∈ℝd, and denote t:=|Iz|. Let A,D and B be the subsystems of U,D⁢([z,1]⁢Q) and WT (respectively) induced by [r]\Iz. We have A∈ℝnu×(r−t),D∈ℝ(r−t)×(r−t), B∈ℝ(r−t)×nw, and

rank⁡(M⁢(z))=rank⁡(A⁢D⁢B).

By Sylvester’s rank inequality, for pair of matrices M1∈ℝn1×r and M2∈ℝr×n2,

rank⁡(M1⁢M2)≥rank⁡M1+rank⁡M2−r.

Therefore,

rank⁡(A⁢D⁢B) ≥rank⁡(A)+rank⁡(D⁢B)−(r−t)
=rank⁡(A)+rank⁡(B)−(r−t)
≥min⁡{r−t,ku}+min⁡{r−t,kw}−(r−t)
⟹rank⁡(M⁢(z))≥min⁡{r−t,ku}+min⁡{r−t,kw}−(r−t). (8)
  1. 1.

    Suppose that t=d.888Recall that t:=|Iz| is the number of zeros in [z,1]⁢V. The Kruskal condition implies that r−d<ku,kw.999If ku≤r−d=r−k+1 then ku+k+kw≤2⁢r+1 which violates the Kruskal condition.

    Equation 8 then implies that

    rank⁡(M⁢(z))≥r−d.

    On the other hand,

    rank⁡(M)=rank⁡(A⁢D⁢B)≤min⁡{nu,r−d,nw}=r−d.

    We have shown that

    Wt⁡([z,1]⁢Q)=r−d⟹rank⁡(M⁢(z))=r−d.
  2. 2.

    Let t<d. By symmetry, we may assume without loss of generality that ku≥kw. The proof proceeds in two cases:

    Case 1: 𝒓−𝒕≤𝒌𝒖.

    Equation 8 implies that

    rank⁡(M)≥min⁡{r−t,kw}>r−d.
    Case 2: 𝒓−𝒕>𝒌𝒖.
    rank⁡(M)≥ku+kw−(r−t)≥(2⁢r+2−k)−(r−t)=(r−d)+(t+1)>r−d.

    In both cases, rank⁡(M)>r−d.

This concludes the proof of Items 1 and 2.

Finally, we prove Item 3:

{z:rank⁡M⁢(z)≤r−d}=⋃I:|I|=d⋂i∈IHi.

Recall the definition of the hyperplanes,

Hi={z∈ℝd|⟨[z,1],qi⟩=0},∀i∈[r].

For every z it holds that z∈⋂i∈IzHi.

  • ■

    If rank⁡M⁢(z)≤r−d, then by Items 1 and 2, Wt⁡([z,1]⁢Q)=r−d

    ⟹|Iz|=d
    ⟹z∈⋃I:|I|=d⋂i∈IHi.
  • ■

    Conversely, if z∈⋃I:|I|=d⋂i∈IHi then there exists I of size d s.t. z∈⋂i∈IHi, so |Iz|≥d. By Equation 7, |Iz|=d. We apply again Item 2 and obtain

    rank⁡M⁢(z)=r−d.

◀

Algorithm MinRank Find minimizers of rank⁡M⁢(z).

Input: (M⁢(z),k′) where

  • ■

    M⁢(z)=U⁢D⁢([z,1]⁢V)⁢WT (this is an affine linear combination of matrices with coefficients {z1,…,zk′−1,1}),

  • ■

    U∈ℝnu×r,V∈ℝk′×r, W∈ℝnw×r, and

  • ■

    nu,nw≤r.

Denote n=min⁡{nu,nw}.
Termination Condition: If at any point the total number of operations performed in the steps below exceeds nc⁢k′, for a large fixed determinable constant c, the algorithm aborts and outputs fail.

  1. 1.

    For t=n,n−1,…,1:

    1. (a)

      Generate the system ℰt by setting all t-minors of M⁢(z) to zero.

    2. (b)

      If ℰt is zero-dimensional, compute the lexicographical Gröbner basis of the radical of ℰt using Theorem 11. Then, determine whether a solution exists by applying Lemma 12.

    3. (c)

      If the solution set is empty:

      1. i.

        Find 𝒮=𝕍⁢(ℰt+1), the solution set of the previous system.

      2. ii.

        Output (t,𝒮).

  2. 2.

    Output fail.

Before proving the correctness of Algorithm MinRank, we emphasize that Theorem 11 determines dimensionality over ℂ. It is possible in principal, that ℰt has only solutions in ℂ\ℝ. In that case, MinRank may output a complex set of minimizers. However, by Theorem 14

MinRank⁢(M⁢(z),k)=r−d,

and all the minimizers are solutions of a real linear system, and are therefore real. Thus, when Algorithm MinRank is called with k′=k, the solution set is real.

Lemma 15.

Let d:=k−1 as in Theorem 14, with k as defined by Equation 6, and suppose that r>9. If V has full Kruskal rank k′, then with probability at least 1−2−r,

  1. 1.

    if k′<k, then MinRank outputs MinRank(M⁢(z),k′) and a list of all the minimizers, or fail, and

  2. 2.

    if k′=k, then MinRank outputs MinRank(M⁢(z),k) and a list of all the minimizers.

Either way, the time complexity is at most nc⁢k′=rO⁢(d).

We defer the proof of the above to the full version [3].

4.2 Step 3: Reconstructing 𝓗

Suppose that Step 2 of the decomposition algorithm was successfully executed. Then we have obtained the set of points

𝒮={z∈ℝd|rank⁡M⁢(z)≤r−d}.

By Theorem 14 part 3, this is also the set of d-wise intersection points of hyperplanes in ℋ, i.e.,

𝒮=⋃I:|I|=d⋂i∈IHi.

Our goal at Step 4 is to reconstruct ℋ. We use the following criterion for identifying hyperplanes in ℋ, the proof of which is deferred to the Appendix: Lemma 27.

Lemma 16 (Hyperplane Identification Criterion).

Let ℋ be a set of r hyperplanes in ℝd in general position, and let 𝒮 be the set of all d-wise intersections of hyperplanes in ℋ. If H is a (d−1)-dimensional hyperplane such that

|H∩𝒮|>rd⁢(r−2d−2),

then H∈ℋ.

Thus, to reconstruct ℋ we can iterate through all the d-size subsets of 𝒮. For each subset, we calculate its affine rank. If it is d−1, we find the normal to the hyperplane through the points. Then we can count the number of intersection points on the hyperplane to determine whether or not the hyperplane belongs to ℋ. However, there are (|𝒮|d)=rO⁢(d2) subsets to explore, which dominate any other step of the decomposition algorithm, so we again apply randomness to reduce time complexity to rO⁢(d). We describe our algorithm for faster sampling of the points to recover these hyperplanes, followed by an analysis of its time complexity and success probability.

Algorithm HP-reconstruction Reconstruct ℋ given 𝒮, the set of d-wise intersection points.

Input: (𝒮,r,t)
Initialize ℋ=∅.
Repeat r2⁢d+t times:

  1. 1.

    Choose d points D⊂𝒮 uniformly at random.

  2. 2.

    If D has affine rank d−1, find the unique hyperplane H through D (up to scaling).

  3. 3.

    If |H∩𝒮|>rd⁢(r−2d−2), and H∉ℋ, add H to ℋ.

If |ℋ|≠r, output fail. Otherwise, output ℋ.

Lemma 17.

Let ℋ={H1,…,Hr} be a set of hyperplanes in general position in ℝd,101010Meaning that the affine dimension of ⋂i∈IHi is max⁡{d−|I|,0} for all I⊆[r]. and let 𝒮 be the set of d-wise intersection points of hyperplanes in ℋ, then when invoked with parameters (𝒮,r,t), Algorithm HP-reconstruction outputs ℋ with probability at least 1−r⋅2−rt and runtime rO⁢(d)+t. In particular, when t=1 Algorithm HP-reconstruction runs in time rO⁢(d) and has success probability at least 1−2−r.

We defer the proof of the above to the full version [3].

4.3 Step 4: restoring 𝑼 and 𝑾 given 𝑸

Once Step 3 outputs Q, which is a random projection of V, Step 4 aims to learn U and W using Q. That is, we are given a tensor T=∑i=1rui⊗vi⊗qi∈ℝn1×n2×n3, one of its components Q, and we want to recover U and V.

We assume that Q has full Kruskal rank kq=nq≤min⁡{ku,kv}, and that the factor matrices satisfy the Kruskal condition. Our approach is to use Q to eliminate kq−2 terms from the decomposition while maintaining a Kruskal rank of 2 in the third mode. This reduction transforms the problem into a Jennrich instance, which can then be solved directly. We formally present the algorithm for this recovery process, followed by its proof of correctness.

Algorithm 5 Decompose T given a component Q.

Input: (T,Q) where

  • ■

    T=∑i=1rui⊗vi⊗qi.

  • ■

    kq=nq.

  • ■

    Kruskal’s Condition: ku+kv+kq≥2⁢r+2

For any matrix M with r columns, denote M1 the sub-matrix induced by columns {1,…,kq−2}, and M2 the sub-matrix induced by columns {kq−1,…,r}.

  1. 1.

    If kq=2 apply Jennrich’s Algorithm and output the decomposition.

  2. 2.

    Else, find full rank P∈ℝ2×kq perpendicular to Q1∈ℝkq×(kq−2) (i.e. the rows of P form a basis for the left kernel of Q1).

  3. 3.

    Apply Jennrich’s algorithm to T×qP=∑i=1rui⊗vi⊗P⁢qi=∑i=kq−1rui⊗vi⊗P⁢qi to obtain a decomposition

    T×qP=⟦U~2,V~2,P⁢Q2⟧.

    Matching U,V,Q: Resolve permutation and scaling ambiguities by aligning the output with the ground-truth matrix P⁢Q2. Apply the resulting permutation and inverse scaling to U~2 and V~2 to maintain consistency with Q.

  4. 4.

    Compute T2:=∑i=kq−1ru~i⊗v~i⊗qi=∑i=kq−1rui⊗vi⊗qi and

    T1:=T−T2=∑i=1kq−2ui⊗vi⊗qi.
  5. 5.

    For each ℓ∈[kq−2], find pℓ∈ℝkq such that ⟨pℓ,qj⟩=𝟙ℓ=j for all j∈[kq−2] (for example, take pℓ to be the ℓ-th row of Q1†). Decompose the rank-1 matrices

    T1×qpℓ=uℓ⊗vℓ

    via SVD to find U~1 and V~1.

Output: U~ and V~.

Lemma 18.

Let T=∑i=1rui⊗vi⊗qi be a rank-r tensor. Given the tensor T and the factor matrix Q=[q1,…,qr] as input, if kq=nq and ku+kv+kq≥2⁢r+2, then Algorithm 5 outputs a decomposition of T in time poly⁢(nu⁢nv⁢nq). Furthermore, under these conditions, the decomposition {ui,vi,qi}i=1r is the only possible rank-r representation of T up to permutation and scaling.

We defer the proof of the above to the full version [3].

4.4 Analysis of DECOMPOSE

As mentioned in Section 2 (Proof Overview), the correctness of DECOMPOSE follows from these 3 facts:

  1. 1.

    Lemma 19: If an r-decomposition satisfies the Kruskal condition, then the tensor has rank r and every rank-decomposition has the same k.

  2. 2.

    Lemma 20: At iterations k′≤k, the decomposition algorithm does not output a “redundant” decomposition. i.e. an r′-decompositions with r′>r.

  3. 3.

    Theorem 21: At iteration k, the algorithm outputs a (unique) r-decomposition (with high probability).

We now prove these 3 statements by order.

Lemma 19 (Uniqueness of “good” parameters).

Let T be a tensor with r-decomposition ⟦V(1),V(2),V(3)⟧ and Kruskal ranks (k1,k2,k3). If (k1,k2,k3) satisfies the Kruskal condition then:

  1. 1.

    rank⁡T=r, and

  2. 2.

    all rank decompositions have Kruskal ranks (k1,k2,k3) (and therefore satisfy the Kruskal condition).

We defer the proof of the above to the full version [3].

Lemma 20.

If k′≤k, then the candidate rank r′ computed by the algorithm at Step 3c satisfies r′≤r, where r is the true rank of the tensor.

We defer the proof of the above to the full version [3]. With the ingredients from the previous sections established, we are now ready to analyze the decomposition algorithm. The full decomposition algorithm includes loops iterating over guesses of k and corresponding mode.

In order to simplify the probability estimates, we assume that r>9. Indeed, for cases where r≤9, one can simply apply a worst-case tensor decomposition algorithm [2, 4], which remains functional even without assuming the Kruskal condition.

Theorem 21 (Decomposition algorithm).

Let T∈ℝnu×nv×nw be a 3-way tensor with an r-decomposition ⟦U,V,W⟧ with U,V,W over ℝ of appropriate dimensions. Suppose that the decomposition satisfies the Kruskal condition and r>9. Let k be as in Eq.6. Then, we have the following.

  1. 1.

    Completeness: By iteration k of the outer-loop, the algorithm outputs the components U,V,W (up to scaling and permutation) with probability at least 1−2−0.1⁢r.

  2. 2.

    Soundness: With probability at least 1−2−0.1⁢r, the algorithm does not output a redundant decomposition ( r′-decomposition with r′>r) at any iteration k′∈{2,…,k} of the outer-loop.

  3. 3.

    Time Complexity: The overall time complexity is

    max⁡(poly⁢(N),(r)O⁢(k))

    where N=nu⁢nv⁢nw.

We defer the proof of the above to the full version [3].

5 Decomposing high dimensional tensors

In the previous section we presented a decomposition algorithm for 3-way tensors. To deal with m-way tensors with m>3, we reshape the tensor and apply the decomposition algorithm on the reshaped tensor. We now introduce some notation to formalize this process. Let T=⟦V(1),V(2),…,V(m)⟧. For all J⊆[m] with |J|≥2 we denote VJ the Khatri-Rao product of {V(j)}j∈J (see Definition 7), and kJ the Kruskal rank of VJ. For any partition I⊔J⊔L=[m], the 3-reshaping of T (see Definition 8) has components VI,VJ and VL. If these matrices satisfy the Kruskal condition, then the decomposition algorithm succeeds (w.h.p) and outputs VI,VJ,VL up to permutation and scaling of the columns. The columns of these matrices are vectorized rank one tensors, so we can directly recover the components of the m-way tensor (via SVD).

Theorem 22 (Reshaped Kruskal Theorem).

Let r≥2, and m≥3. Consider an m-way tensor T∈ℝn1×n2×⋯×nm with decomposition:

T=∑j=1rvj(1)⊗vj(2)⊗⋯⊗vj(m) (9)

For each subset ∅≠I⊆[m], define VI as the Khatri-Rao product of matrices {V(i)}i∈I. If there exists a partition I⊔J⊔L=[m] and the Kruskal ranks of these partitioned matrices satisfy the Kruskal condition:

kI+kJ+kL≥2⁢r+2 (10)

then,

  1. 1.

    rank⁡T=r, and

  2. 2.

    the decomposition is unique up to scaling and permutation.

Furthermore, for r>9, there exists a randomized algorithm that recovers the unique components with probability at least 1−2−0.1⁢r and runtime NO⁢(k), where k:=2⁢r+2−max⁡{kI+kJ,kI+kL,kJ+kL} and N:=∏i∈[m]ni is the input size.

We defer the proof of the above to the full version [3].

If the input m-dimensional tensor satisfies the Kruskal condition, then there is an efficient way to perform the 3-reshaping. As a direct consequence, we obtain the following corollary. We defer the reshaping procedure and its correctness to the Appendix, see Lemma 24.

Corollary 23 (Kruskal’s Uniqueness Theorem).

Consider an m-way tensor T∈ℝn1×n2×…×nm with decomposition:

T=∑j=1rVj(1)⊗Vj(2)⊗⋯⊗Vj(m) (11)

If the Kruskal condition holds:

∑t=1mkt≥2⁢r+m−1 (12)

then,

  1. 1.

    rank⁡T=r, and

  2. 2.

    the decomposition is unique up to scaling and permutation.

Furthermore, for r>9, there exists a randomized algorithm that outputs the unique components with probability at least 1−2−0.1⁢r and runtime NO⁢(k), where k:=2⁢r+2−max⁡{kI+kJ,kI+kL,kJ+kL} and N:=∏i∈[m]ni.

Again, for cases where r≤9, one can simply apply a worst-case tensor decomposition algorithm [2, 4], which will decompose these m-tensors even without the Kruskal condition in poly⁢(m,max⁡ni) time.

References

  • [1] M. Bardet, M. Bros, D. Cabarcas, P. Gaborit, R. A. Perlner, D. Smith-Tone, J.-P. Tillich, and J. A. Verbel. Improvements of algebraic attacks for solving the rank decoding and minrank problems. In ASIACRYPT, 2020.
  • [2] Vishwas Bhargava, Shubhangi Saraf, and Ilya Volkovich. Reconstruction algorithms for low-rank tensors and depth-3 multilinear circuits. Technical Report 045, Electronic Colloquium on Computational Complexity (ECCC), 2021. Version with full proofs; subsequently appeared in STOC 2021. URL: https://eccc.weizmann.ac.il/report/2021/045/.
  • [3] Vishwas Bhargava, Leonard Schulman, and Shiri Sivan. An algorithmic proof of kruskal’s tensor decomposition theorem. Technical Report TR26-073, Electronic Colloquium on Computational Complexity (ECCC), 2026. URL: https://eccc.weizmann.ac.il/report/2026/073/.
  • [4] Vishwas Bhargava and Devansh Shringi. Faster & Deterministic FPT Algorithm for Worst-Case Tensor Decomposition. In 52nd International Colloquium on Automata, Languages, and Programming (ICALP 2025), volume 334 of Leibniz International Proceedings in Informatics (LIPIcs), pages 28:1–28:20. Schloss Dagstuhl – Leibniz-Zentrum für Informatik, 2025. doi:10.4230/LIPIcs.ICALP.2025.28.
  • [5] Aditya Bhaskara, Moses Charikar, Ankur Moitra, and Aravindan Vijayaraghavan. Open problem: Tensor decompositions: Algorithms up to the uniqueness threshold? In Conference on Learning Theory, pages 1280–1282. PMLR, 2014. URL: http://proceedings.mlr.press/v35/bhaskara14b.html.
  • [6] J.D. Carroll and J.J. Chang. Analysis of individual differences in multidimensional scaling via an n-way generalization of “Eckart-Young” decomposition. Psychometrika, 35:283–319, 1970. doi:10.1007/BF02310791.
  • [7] Luca Chiantini and Giorgio Ottaviani. On generic identifiability of 3-tensors of small rank. SIAM Journal on Matrix Analysis and Applications, 33(3):1018–1037, 2012. doi:10.1137/110829180.
  • [8] Luca Chiantini, Giorgio Ottaviani, and Nick Vannieuwenhoven. Effective criteria for specific identifiability of tensors and forms. SIAM Journal on Matrix Analysis and Applications, 38(2):656–681, 2017. doi:10.1137/16M1090132.
  • [9] David A. Cox, John Little, and Donal O’Shea. Using Algebraic Geometry, volume 185 of Graduate Texts in Mathematics. Springer, New York, 1998.
  • [10] Harm Derksen. Kruskal’s uniqueness inequality is sharp. Linear Algebra and its Applications, 438(2):708–712, 2013. Tensors and Multilinear Algebra. doi:10.1016/j.laa.2011.05.041.
  • [11] J. Faugère, M. Safey El Din, and P. Spaenlehauer. Computing loci of rank defects of linear matrices using Gröbner bases and applications to cryptology. In International Symposium on Symbolic and Algebraic Computation (ISSAC), pages 257–264, 2010. doi:10.1145/1837934.1837984.
  • [12] Richard A. Harshman. Foundations of the PARAFAC procedure: Models and conditions for an “explanatory” multi-model factor analysis. In UCLA working papers in phonetics, 16(1):84, 1970. URL: https://api.semanticscholar.org/CorpusID:6816804.
  • [13] Johan Håstad. Tensor rank is NP-complete. Journal of Algorithms, 11(4):644–654, 1990. doi:10.1016/0196-6774(90)90014-6.
  • [14] Christopher J Hillar and Lek-Heng Lim. Most tensor problems are NP-hard. Journal of the ACM (JACM), 60(6):1–39, 2013. doi:10.1145/2512329.
  • [15] Doug Ierardi. Quantifier elimination in the theory of an algebraically-closed field. In Proceedings of the twenty-first annual ACM symposium on Theory of computing, pages 138–147, 1989. doi:10.1145/73007.73020.
  • [16] Pravesh K Kothari, Ankur Moitra, and Alexander S Wein. Overcomplete tensor decomposition via Koszul-Young flattenings. arXiv preprint arXiv:2411.14344, 2024. doi:10.48550/arXiv.2411.14344.
  • [17] J. B. Kruskal. More factors than subjects, tests and treatments: an indeterminacy theorem for canonical decomposition and individual differences scaling. Psychometrica, 41(3):281–293, 1976.
  • [18] J. B. Kruskal. Rank, decomposition and uniqueness for 3-way and n-way arrays. In Multiway data analysis, pages 7–18. Elsevier Science Publishers B.B. (North-Holland, 1989.
  • [19] Joseph B. Kruskal. Three-way arrays: Rank and uniqueness of trilinear decompositions, with application to arithmetic complexity and statistics. Linear Algebra and its Applications, 18(2):95–138, 1977. doi:10.1016/0024-3795(77)90069-6.
  • [20] Yagati N Lakshman and Daniel Lazard. On the complexity of zero-dimensional algebraic systems. In Effective methods in algebraic geometry, pages 217–225. Springer, 1991.
  • [21] Joseph M Landsberg. Tensors: geometry and applications: geometry and applications, volume 128. American Mathematical Soc., 2011.
  • [22] Benjamin Lovitz and Fedor Petrov. A generalization of Kruskal’s theorem on tensor decomposition. In Forum of Mathematics, Sigma, volume 11, page e27. Cambridge University Press, 2023.
  • [23] Tengyu Ma, Jonathan Shi, and David Steurer. Polynomial-time tensor decompositions with sum-of-squares. In 2016 IEEE 57th Annual Symposium on Foundations of Computer Science (FOCS), pages 438–446. IEEE, 2016. doi:10.1109/FOCS.2016.54.
  • [24] S. Rabanser, O. Shchur, and S. Günnemann. Introduction to tensor decompositions and their applications in machine learning, 2017. arXiv:1711.10781.
  • [25] John A. Rhodes. A concise proof of Kruskal’s theorem on tensor decomposition. Linear Algebra and its Applications, 432(7):1818–1824, 2010. doi:10.1016/j.laa.2009.11.033.
  • [26] Marcus Schaefer and Daniel Štefankovič. The complexity of tensor rank. Theory of Computing Systems, 62(5):1161–1174, 2018. doi:10.1007/S00224-017-9800-Y.
  • [27] Nicholas D Sidiropoulos and Rasmus Bro. On the uniqueness of multilinear decomposition of n-way arrays. Journal of Chemometrics: A Journal of the Chemometrics Society, 14(3):229–239, 2000.
  • [28] Charles Spearman. “General intelligence” objectively determined and measured. In J. J. Jenkins & D. G. Paterson, editor, Studies in individual differences: The search for intelligence, pages 59–73. Appleton-Century-Crofts, 1961. Reprint from American Journal of Psychology, 1904, Vol 15[2], 201-293. doi:10.1037/11491-006.
  • [29] Alwin Stegeman and Nicholas D Sidiropoulos. On Kruskal’s uniqueness condition for the Candecomp/Parafac decomposition. Linear Algebra and its applications, 420(2-3):540–552, 2007.

Appendix A Reshaping into 3-tensors

Lemma 24.

If T is an m-way tensor satisfying the Kruskal condition, then there exists a partition of [m] into disjoint subsets J⊔K⊔L=[m] s.t. ⟦VJ,VK,VL⟧ satisfies the Kruskal condition, namely,

kJ+kK+kL≥2⁢r+2.

We defer the proof of the above to the full version [3].

A.1 General position statements

In the first step of the decomposition algorithm (Algorithm DECOMPOSE) we apply a random matrix on one of the components of the tensor we are trying to decompose. We obtain a component Q∈ℝk×r and define, for every column, the orthogonal hyperplane: Hi:={z∈ℝd|⟨[z,1],qi|⟩=0}. We denote the full set of hyperplanes by

ℋ={H1,…,Hr}.
Lemma 25.

If Q∈ℝk×r has Kruskal rank k, and Q⁢([d],[r]) has Kruskal rank d=k−1, then ℋ is in general position in ℝd. That is, for all I⊆[r], the affine dimension of ⋂i∈IHi is max⁡{d−|I|,0}.

We defer the proof of the above to the full version [3].

When ℋ is in general position, we can count the number of d-wise intersection points on the flats defined by ℋ:

Corollary 26.

If ℋ is in general position in ℝd, then

𝒮:={⋂i∈IHi|I⊆[r],|I|=d}

has size (rd). More generally, for any I⊆[r], |I|=t≤d,

|⋂i∈IHi∩𝒮|=(r−td−t).

We defer the proof of the above to the full version [3].

Lemma 27.

Let r>d, d≥2. We have

  • ■

    ℋ, a set of r hyperplanes of dimension d−1 in general position in ℝd.

  • ■

    𝒮, the set of d-wise intersections of hyperplanes in ℋ.

If H is a (d−1)-dimensional hyperplane, and

|H∩𝒮|>rd⁢(r−2d−2),

then H∈ℋ.

We defer the proof of the above to the full version [3].

Corollary 28.

For any hyperplane H,

|H∩𝒮|≤(r−1d−1)

We defer the proof of the above to the full version [3].

Lemma 29.

Let 𝒫={s1,…,sd} be a set of points chosen uniformly at random from 𝒮. Then, for any Hℓ∈ℋ, 𝒫 spans Hℓ with probability at least 1r2⁢d.

We defer the proof of the above to the full version [3].