Canonical Polyadic Decomposition of Third-Order Tensors of Which a Mode Rank is Equal to the Rank

(1)

Canonical Polyadic Decomposition of

Third-Order Tensors of Which a Mode Rank is

Equal to the Rank

∗

Lieven De Lathauwer

December 1, 2011

1 Introduction

We consider the Canonical Polyadic Decomposition (CPD) of a rank-R tensor T ∈ KI1×I2×I3, with R = rank 3(T ) 6 I3: T = R X r=1 ar ⊗br ⊗cr. (1) Contributions:

• A new, relaxed, algebraic condition on ar, br, 1 6 r 6 R that is

sufficient for essential uniqueness of the CPD.

• Under this condition the CP factors can be obtained by means of stan-dard linear algebra (set of linear equations and matrix Eigenvalue De-composition (EVD)).

∗_{Available in SISTA publication engine as: L. De Lathauwer, CPD Study, Part III,}

(2)

• Algorithm for the estimation of the CP factors in the case of an ap-proximate CPD.

• Generic uniqueness for a given set (I1, I2, R) is established by checking

the algebraic condition for one random sample. • Verification for (3 × 7 × 12) tensors of rank 12.

• Similar approach for ar = br, 1 6 r 6 R. (Generic bound will be

different.)

Preceding work. In matrix format we have:

T[1,2;3] = (A ⊙ B) · CT. (2)

A necessary condition for uniqueness of the CPD is that A ⊙ B has full column rank. Our mode-3 rank assumption implies that rank(C) = R. In [1] we showed that (1) can be obtained by means of standard linear algebra (set of linear equations and matrix Eigenvalue Decomposition (EVD)) if the tensors in the set U = {Φ(as1b

T s1, as2b

T

s2)|1 6 s1 < s2 6 R} are linearly

independent. In the case of an approximate decomposition, after solving a set of linear equations in least-squares sense, the computation of (1) can be reduced to the computation of a CPD of the form

M =

R

X

r=1

wr ⊗wr ⊗lr, (3)

in which M ∈ KR×R×R _{and in which W ∈ K}R×R _{and L ∈ K}R×R _are

nonsin-gular. Taking into account symmetries [2], the tensors in the set U can be represented by the columns of a matrix U ∈ KCI12 C

2 I2×C

2 R

, which needs to have full rank. Consequently, in this approach R is bounded by C2

R 6CI21C

2 I2. In

the next section we propose a technique that is less restrictive.

2 Symmetrized mapping

Let F ∈ KI1I2×R

and G ∈ KI3×R _{be such that}

(3)

is a rank-revealing decomposition of T[1,2;3]. Then we have:

F = (A ⊙ B) · W−₁

, (5)

GT = W · CT, (6)

in which W ∈ KR×R _{is nonsingular.}

Theorem 1 Consider the mapping Φ(3) _{: (X}(r1)_{, X}(r2)_{, X}(r3)

) ∈ KI×J _×

KI×J _{× K}I×J _{→ Φ}(3)_(X(r1)_{, X}(r2)_{, X}(r3)) ∈ KI×I×I×J ×J ×J _{defined by}

(Φ(3)_(X(r1)_{, X}(r2)_{, X}(r3)₎₎ i1i2i3j1j2j3 = P σ∈P(r1,r2,r3)(x (σ(r1)) i1j1 x (σ(r2)) i2j2 − x (σ(r1)) i1j2 x (σ(r2)) i2j1 )x (σ(r3)) i3j3 ,

in which P(r1, r2, r3) is the set of permutations of (r1, r2, r3). Then Φ(3)(X, X, X) =

O if and only if X has at most rank 1.

Note that Φ(3)_{(·, ·, ·) is symmetric in its arguments. Define}

Pr1r2r3 = Φ (3)_((f r1)[1;2],(fr2)[1;2],(fr3)[1;2]), 1 6 r1 6r2 6r3 6R. (7) We have that Pr1r2r3 = R X s1,s2,s3=1 (W−₁ )s1r1(W −₁ )s2r2(W −₁ )s3r3Φ (3)_(a s1b T s1, as2b T s2, as3b T s3). (8) We look for symmetric tensors X ∈ KR×R×R _{that satisfy}

R X r1,r2,r3=1 xr1r2r3Pr1r2r3 = O. (9) We stack Pr1r2r3, 1 6 r1 6 r2 6 r3 6 R, in a matrix ¯P ∈ K I3 1I 3 2×C 3 R+2.

We determine the vectors x in the null space of ¯P, which can then be mapped to tensors X . We make the crucial assumption that the set of tensors {Φ(3)_(a s1b T s1, as2b T s2, as3b T s3)|1 6 s1 < s2 6 s3 6 R} is linearly independent.

Combination of (7)–(9) shows that, under this assumption, the null space of ¯P has dimension R. Moreover, the vectors in the null space represent

(4)

tensors that are linear combinations of wr ⊗wr ⊗wr, 1 6 r 6 R. If we stack

the tensors obtained from a basis of the null space of ¯P in a super-tensor Y ∈ KR×R×R×R_{, then this super-tensor admits the CPD}

Y =

R

X

r=1

wr ⊗wr ⊗wr ⊗vr. (10)

The vectors vr, 1 6 r 6 R, are linearly independent since they generate a

basis of the null space of ¯P. Since additionally W is nonsingular, CPD (10) is unique. In the case of an exact decomposition, (10) follows from a matrix EVD [3, 4]. In the case of an approximate decomposition, the EVD may be used to initialize an algorithm for fitting (10).

Remark 1. Instead of working with ¯P_{∈ K}I13I 3 2×C

3

R+2, one may also estimate

the null space of ¯PH _{· ¯}_P_{∈ K}C3 R+2×C

3

R+2. The latter matrix can be computed

efficiently by computing (fr1⊗fr2⊗fr3) H_(f s1⊗fs2⊗fs3) as (f H r1fs1)(f H r2fs2)(f H r3fs3).

The disavantage is that the condition number is squared.

Example 1. Consider tensors of dimension (3 × 7 × R) and rank R. The number of free parameters in the CPD equals (3+7+R−2)R = (8+R)R. The number of tensor entries is 21R. Consequently, for R = 13 there generically exist a finite number of decompositions. Moreover, we expect that for R 6 12 the CPD is generically unique. For R = 12 uniqueness is not guaranteed by the result in [1]. However, it can be verified that the 352 tensors in the set {Φ(3)_(a s1b T s1, as2b T s2, as3b T

s3)|1 6 s1 < s2 6s3 612} are linearly independent

for a random example. Hence, for R = 12 the decomposition is generically unique and it can be computed as explained above.

Remark 2. It is possible to generalize the trilinear mapping Φ(3) _{to a}

K-linear mapping Φ(K) _{that is supposed to yield a matrix ¯}_P _{∈ K}IK 1 I

K 2 ×C

K R+K−1

that has an R-dimensional null space.

(The trick here is that by increasing K we work with a matrix of which the number of rows increases faster than the number of columns. Since the columns have a multiplicative rather than a linear structure, this will hopefully lead to more relaxed conditions and algorithms.)

(5)

3 Alternative mapping

Theorem 2 Consider the mappingΨ(3) _{: (X, Y, Z) ∈ K}I×J

×KI×J ×KI×J _→ Ψ(3)_{(X, Y, Z) ∈ K}I×I×J ×J ×6 _{defined by} (Ψ(3)_{(X, Y, Z))} i1i2j1j21 = xi1j1yi2j2 − xi1j2yi2j1, (Ψ(3)_{(X, Y, Z))} i1i2j1j22 = yi1j1xi2j2 − yi1j2xi2j1, (Ψ(3)_{(X, Y, Z))} i1i2j1j23 = xi1j1zi2j2 − xi1j2zi2j1, (Ψ(3)_{(X, Y, Z))} i1i2j1j24 = zi1j1xi2j2 − zi1j2xi2j1, (Ψ(3)_{(X, Y, Z))} i1i2j1j25 = yi1j1zi2j2 − yi1j2zi2j1, (Ψ(3)_{(X, Y, Z))} i1i2j1j26 = zi1j1yi2j2 − zi1j2yi2j1.

Then Ψ(3)_{(X, X, X) = O if and only if X has at most rank 1.}

Define Pr1r2r3 = Ψ (3)_((f r1)[1;2],(fr2)[1;2],(fr3)[1;2]), 1 6 r1 6r2 6r3 6R. (11) We have that Pr1r2r3 = R X s1,s2,s3=1 (W−1₎ s1r1(W −1₎ s2r2(W −1₎ s3r3Ψ (3)_(a s1b T s1, as2b T s2, as3b T s3). (12) We look for tensors X ∈ KR×R×R _{that satisfy}

R X r1,r2,r3=1 xr1r2r3Pr1r2r3 = O. (13) We stack Pr1r2r3, 1 6 r1, r2, r3 6 R, in a matrix ¯P ∈ K 6I2 1I 2 2×R 3 . We de-termine the vectors x in the null space of ¯P, which can then be mapped to tensors X . We make the crucial assumption that the set of tensors {Ψ(3)_(a s1b T s1, as2b T s2, as3b T s3)|1 6 s1, s2, s3 6 R,¬(s1 = s2 = s3)} is linearly

independent. Combination of (11)–(13) shows that, under this assumption, the null space of ¯Phas dimension R. Moreover, the vectors in the null space represent tensors that are linear combinations of wr ⊗wr ⊗wr, 1 6 r 6 R.

(6)

If we stack the tensors obtained from a basis of the null space of ¯P in a super-tensor Y ∈ KR×R×R×R_{, then this super-tensor admits the CPD}

Y =

R

X

r=1

wr ⊗wr ⊗wr ⊗vr. (14)

The vectors vr, 1 6 r 6 R, are linearly independent since they generate a

basis of the null space of ¯P. Since additionally W is nonsingular, CPD (14) is unique. In the case of an exact decomposition, (14) follows from a matrix EVD [3, 4]. In the case of an approximate decomposition, the EVD may be used to initialize an algorithm for fitting (14).

Remark 3. Since all tensors X are symmetric, we might as well have worked with a matrix ¯¯P _{∈ K}6I12I22×CR+23 in which the columns corresponding to the

same entry are added together. In other words, we might have symmetrized Ψ(3)_.

Remark 4. It seems that the mapping Ψ(3) _{does not yield a uniqueness}

condition that is more relaxed than the one in [1]. However, Ψ(3) _might

perhaps be used as a starting point for further research.

References

[1] L. De Lathauwer, “A Link between the Canonical Decomposition in Mul-tilinear Algebra and Simultaneous Matrix Diagonalization”, SIAM J. Ma-trix Anal. Appl., Vol. 28, No. 3, 2006, pp. 642–666.

[2] A link between the decomposition of a third-order tensor in rank-(L,L,1) terms and simultaneous matrix diagonalization, Tech. report.

[3] R.A. Harshman, Foundations of the PARAFAC procedure: Model and conditions for an“explanatory” multi-mode factor analysis, UCLA Work-ing Papers in Phonetics, 16 (1970), pp. 1–84.

[4] S.E. Leurgans, R.T. Ross, and R.B. Abel, A decomposition for three-way arrays, SIAM J. Matrix Anal. Appl., 14 (1993), pp. 1064–1083.