Beyond Unfolding: Exact Recovery of Latent Convex Tensor Decomposition
under Reshuffling
Chao Li,1 Mohammad Emtiyaz Khan,1 Zhun Sun,1
Gang Niu,1 Bo Han,1 Shengli Xie,2,∗ Qibin Zhao1,∗
2

1
RIKEN Center for Advanced Intelligence Project (AIP), Japan
School of Automation, Guangdong University of Technology, China

Abstract
Exact recovery of tensor decomposition (TD) methods is a
desirable property in both unsupervised learning and scientific data analysis. The numerical defects of TD methods,
however, limit their practical applications on real-world data.
As an alternative, convex tensor decomposition (CTD) was
proposed to alleviate these problems, but its exact-recovery
property is not properly addressed so far. To this end, we focus on latent convex tensor decomposition (LCTD), a practically widely-used CTD model, and rigorously prove a sufficient condition for its exact-recovery property. Furthermore,
we show that such property can be also achieved by a moregeneral model than LCTD. In the new model, we generalize the classic tensor (un-)folding into reshuffling operation,a more flexible mapping to relocate the entries of the matrix into a tensor. Armed with the reshuffling operations andexact-recovery property, we explore a totally novel application for (generalized) LCTD, i.e., image steganography. Experimental results on synthetic data validate our theory, andresults on image steganography show that our method outperforms the state-of-the-art methods.

Introduction
Tensor decomposition (TD), a multi-linear extension of matrix factorization, has been successfully employed on various applications (Rabusseau and Kadri 2016; Sharan andValiant 2017; He et al. 2017). More importantly, TD methods are also crucial tools for unsupervised discovery of
structures hidden behind the data, e.g., localizing the regionsof brain from EEG waveforms (Becker et al. 2014), user detection from mobile-communication data, and understanding the kinetic-theory description of materials (González etal. 2010). One of the reasons, behind the success of TD
methods in these tasks, is due to the exact-recovery of theirsolutions, i.e., TD methods are able to achieve a demixing ofdata in which individual component can tightly correspondwith physical interpretation (Williams et al. 2018).
The numerical properties of TD methods, however, are
not as promising as their exact-recovery property. The mostpopular canonical polyadic decomposition (CPD) provides∗

Corresponding authors. Email:{chao.li, qibin.zhao}@riken.jp.
Copyright c 2020, Association for the Advancement of Artificial
Intelligence (www.aaai.org). All rights reserved.

a uniqueness solution up to permutation, but the opennessof its solution space and the weak express power are too
restrictive for higher-dimensional problems (Comon, Luciani, and De Almeida 2009). Other alternatives, such as(heretical) Tucker decomposition (Cichocki et al. 2009;
Song et al. 2013), tensor-train decomposition (Oseledets
2011) and tensor-ring decomposition (Zhao et al. 2016), unfortunately do not possess the exact-recovery property, i.e.,their components can be arbitrarily rotated without changing
the resultant. Block term decomposition (BTD) (De Lathauwer 2008), a constrained version of Tucker decomposition, inherits the uniqueness property from CPD but has
more flexible decomposition form. However, the model selection for BTD, such as the determination of the multilinear rank for each block, would be highly challenging.
Furthermore, the non-convexity of the aforementioned TD
methods generally leads to unstable convergence to global
minimum (Kolda and Bader 2009).
Convex-optimization-based approaches were proposed to alleviate these unsatisfactory numericalproblems (Imaizumi, Maehara, and Hayashi 2017;
Tomioka, Hayashi, and Kashima 2010; Yu et al. 2019). In
these methods, it is not necessary to specify the rank of thedecomposition beforehand, and the convexity of the models
guarantees both the convergence to the global minimum
and their statistical performance. However, there is an
important issue that is not properly addressed so far: Are
the convex approaches able to exactly recover the low-rankcomponents like their non-convex counterparts?
To answer this question, we theoretically prove a sufficient condition for exact-recovery of latent convex tensor decomposition (LCTD), a practically widely-used convex TD
method (Yamada et al. 2017; Nimishakavi, Jawanpuria, andMishra 2018) that decomposes a tensor into a mixture of
low-rank components. Armed with the notion of incoherence among the components, we rigorously prove that the
low-rank components can be exactly recovered when a typeof incoherence measure is sufficiently small. Moreover, weshow that the exact-recovery property can be owned by amore general class of models than LCTD. In the new model,
we introduce the reshuffling operation to replace the conventional tensor (un-)folding used in LCTD, and the new reshuf-

fling operations can give the model the capacity to explore
richer low-rank structures. Last, by leveraging the reshuffling operations and the exact-recovery property, we explore
a totally novel application of (generalized) LCTD, i.e., image steganography, a classic task in both computer visionand information security. Experimental results not only validate the correctness of our theory, but also demonstrate themodel’s effectiveness in real-world application. Supplementary materials are available at: http://qibinzhao.github.io.

Related Works
The notion of convex tensor decomposition (CTD) was
firstly introduced in (Tomioka, Hayashi, and Kashima2010), where the decomposition was implemented by minimizing a type of tensor nuclear norm. Subsequently, thestudies on CTD, especially on latent convex tensor decomposition (LCTD), were continually concerned on both theoretical and practical sides. LCTD was theoretically provedto achieve tighter upper-bound on the reconstruction performance than its overlapped counterpart (Tomioka andSuzuki 2013), and had achieved the state-of-the-art results in many tasks (Guo, Yao, and Kwok 2017; Wanget al. 2019). Works on variants of LCTD are richly proposed recently (Nimishakavi, Jawanpuria, and Mishra 2018;
Wimalawarne, Yamada, and Mamitsuka 2017), but surprisingly the exact-recovery property of LCTD is not paid muchattention so far. More interestingly, in the first paper proposing LCTD (Tomioka, Hayashi, and Kashima 2010), the authors stated that LCTD might not be able to exactly recover
the components due to its relation with Tucker model.
For the exact-recovery property of tensor decomposition,
the solution of CPD is unique up to permutation under
mild condition, while Tucker and TT decompositions arenot. However, as aforementioned, the numerical problemsof these methods, such as the rank determination, somehow
lead to difficult implementation in practice. On the otherside, many approaches have focused on restricting the modelsuch that the exact-recovery property is ensured. For example, one way to eliminate the ambiguity in Tucker decomposition is to incorporate additional constraints on the latent
factors, e.g., by forcing them to be independent, sparse, orsmooth (Cichocki et al. 2009). This could work in practice,
but in some cases, these constraints might be too strong forthe data in hand. In contrast, we focus on exactly recovering
the components by only exploiting the low-rank structuresof the tensor. Additional assumptions such as sparsity, independence are temporarily out of the scope in our work.
The (un-)folding operation used in TD methods builds a
connection of the low-rank structures between a matrix andits higher-order form. In the existing TD methods, this operation is defined by various manipulations (Mu et al. 2014; Yu
et al. 2019), but one basic principal behind them is to keep
more low-rankness of the tensor along the modes. However,
the recent studies on exploiting the low-rank tensor decomposition under general linear transformations (Li et al. 2019;
Lu, Peng, and Wei 2019) inspire us that the low-rank structures of a tensor can be explored by more flexible operations
than tensor (un-)folding. In this paper, we therefore generalize the conventional tensor (un-)folding into reshuffling,which gets rid of the stereotype in the existing definitionsof tensor (un-)folding. In contrast to using arbitrary linear
transformation, the proposed reshuffling operations only relocate the entries (without up-sampling) of the data without addition and multiplication, which would result in lower
FLOPS in practice.

Latent Convex Tensor Decomposition
Throughout the paper, we will denote matrices by boldface
capital letters, e.g., X ∈ Rm×n is a matrix of size m × n.
We will denote tensors by X ∈ RI1 ×I2 ×···×IK , where K is
the order of the tensor. Given data in the form of a tensor X ,in this paper we consider low-rank tensor decomposition as
the sum of multiple components:
X = A1 + A2 + · · · + AN ,(1)

where Ai are low-rank components of the same size as
X . One may argue that the form (1) is not low-rank tensor decomposition as it is not in a multiplication of latent
factors (Kolda and Bader 2009). Analogous to the singularvalue decomposition (SVD) of a matrix, note that both themost popular CPD and Tucker decomposition can be trivially rewritten as (1), and their ranks determine the numberof components (Cichocki et al. 2009).
Latent convex tensor decomposition (LCTD) can be also
formulated as similar as (1). LCTD decomposes a tensor into
the sum of components but incorporating the (un-)folding
operations on each component, i.e.,
X = Φ1 (A1 ) + Φ2 (A2 ) + · · · + ΦK (AK ),(2)

where the matrices Ai are assumed to be low-rank, andthe operations Φi (·) are also called folding or tensorization in literature. In contrast to CPD and Tucker decomposition, LCTD explore the low-rank structures of the matricized form of each component, and its algorithm is based onconvex optimization.
The exact-recovery property of LCTD, however, is still akind of open problem. The consistency of Ai in (2) was discussed in (Tomioka and Suzuki 2013), but the reconstruction bound does not tend to be zero when decreasing thestrength of the noise. It implies that there is still a theoreticalgap between consistency and identiability, a.k.a., the exactrecovery of components.
In this paper, we will try to fill this gap by imposing
the notion of incoherence as the condition for the exactrecovery. We find that the different (un-)folding operations
on components could bring us incoherent low-rank structures, which has been proved as an important characteristic in compressed sensing and matrix completion (Candès
and Recht 2009). Meanwhile, we also find by the numerical experiments that the unbalance between the number ofrows and columns of the unfolded tensor is probably one
key reason that restricts the capacity for exactly recovering
the components. To this end, we generalize the conventional
tensor (un-)folding into reshuffling, which is able to result in
more balanced and incoherent model.

Reshuffled Tensor Decomposition
In this section, we introduce a generalization of LCTD,
where low-rank matrices are mapped to tensors by the
reshuffling operations. To be self-contained of the paper, we
also briefly derive an algorithm with stable convergence.

Formulation
We consider the following model similar to (2):
X = R1 (A1 ) + R2 (A2 ) + · · · + RN (AN ),(3)

where Ai ∈ Rmi ×ni , i ∈ [N ] denotes the latent components. However, in contrast to (2), we arbitrarily choose thenumber of components, and replace the folding operations
Φi (·) with a more general and flexible type for operations
Ri : Rmi ×ni → RI1 ×I2 ×···×IK , which are a group of linear
operators called reshuffling and defined as follows:
Definition 1 (Reshuffling) The reshuffling operation, denoted by R, is defined as a mapping that maps a matrix A ∈ Rm×n into a real tensor A := R(A) of size
I1 × I2 × · · · × IK , such that the number of elements in
A is equal to the number of elements in A, i.e., mn =I1 I2 . . . IK , and every entry of A correspond to one andonly one entry in A.
The main idea behind our model is to employ a variety ofreshuffling operations for every component, i.e., for eachcomponent Ai , we reshuffle it using a distinct operator Ri to
get a tensor Ai := Ri (Ai ). Meanwhile, we assume that therank of Ai to be small, and would like to be able to capture
a variety of low-rank structures within the tensor.
We build a convex optimization problem for recovering
the components Ai from Eq. (3). Assuming that the reshuffling operations Ri are known for each component, Ai , i ∈[N ] are recovered by minimizing the following optimization
problem:min

A1 ,...,AN

N
X
i=1

kAi k∗ , s.t., X =N
Xto exploit some specific characteristics of the data. For example, the unfolding operation exploits the physical meaning associated with the modes to convert the tensor into amatrix. Reshuffling operations could be designed in a similar fashion to exploit other types of structural information
about the tensor. In this paper, we focus on discussing theexact-recovery (conditions) of the method and assume thereshuffling operations to be known in advance.

Algorithm
Below, we derive an algorithm called reshuffled tensor decomposition (Reshuffled-TD) to solve (4).
Due to the existence of the equality constraints, we apply
the augmented Lagrangian method for solving (4), of which
the Lagrangian function is given by
L (A1 , . . . , AN , Y, κ) =N
X

kAi k∗

i=1*+

Y, X −NX
i=1+
Ri (Ai )

N
X
κ
+
X−
Ri (Ai )
2
i=1

2

,
F

(5)
where the tensor Y denotes the Lagrangian dual and κ > 0
is a positive scalar. In the algorithm, we sequentially update
Ai , i ∈ [N ], Y and κ in each iteration. As the key step of
the algorithm, we update Ai by minimizing the following
sub-problem:
1A+
i = arg min kAi k∗
Ai κ

 2
X
1
Ai − Ri? X −Rj (Aj ) + κ−1 Y  ,
+
2
j6=i

F

(6)
where Ri? denotes the adjoint operator of Ri . In contrast toRi (Ai ),

(4)

i=1

where k · k∗ denotes the matrix nuclear norm, which equals
the sum of the singular values of matrix. It has been provedthat the matrix nuclear norm is the convex envelope of matrix rank (Fazel, Hindi, and Boyd 2001). Hence, we canroughly say that solving the problem (4) is equivalent tolooking for the most low-rank latent components from theobserved tensor X .
Remark 1: Similarly to the conventional folding operation, reshuffling maps a matrix into a tensor, and is a linear and reversible operator. However, reshuffling can operate
more flexibly and handle the matrices of arbitrary size. The
conventional tensor (un-)folding can be therefore obtained
as a special case of reshuffling. The flexibility of reshuffling
could enable recovery of low-rank structures that were previously unrecoverable by using folding operations.
Remark 2: The reshuffling operations for particular data
could be difficult to find. In some applications, such as
stenography (shown in Experiment Section), the operations
are known beforehand. In general, one could design them

the renowned Alternating Direction Method of Multipliers(ADMM) methods, we update the scalar κ in each iteration
by multiplying a constant ρ larger than 1 (for example ρ =1.01), i.e., κ+ = ρκ, and the work in (Lin, Chen, and Ma
2010) shows that such trick could efficiently accelerate theconvergence compared to ADMM. The complete procedureof Reshuffled-TD is given in Alg. 1, where Dα ( · ) denotes
the soft-thresholding operation on the singular values (Cai,
Candès, and Shen 2010).
Next, we show that the algorithm results in convergenceto the optimal point of (4). For the brevity of the proof, weonly consider the case when N = 2 in the theorem, and fromthe experimental results we also find the stable convergenceeven in the case with more components (N ≥ 3).
Theorem 1 Using Reshuffled-TD and assuming the numberof the components N = 2, if the sequence {κ(k) } is nonP∞(k)(k)
decreasing, and k=1 κ(k),−1 = +∞, then (A1 , A2 )converges to an optimal solution of (4).
The proof is given in the supplemental material.
It is worth noting that the convergence of Reshuffled-TD
does not ensure whether the solution is equal to the “true”

Algorithm 1 Reshuffled-TD
Initialize: The observation X , Y = sgn(X ), Ai =Ri? (X )/N, ∀i, and ρ > 1, κ0 > 0
Iteration until convergence://STEP 1: update every latent component Ai , i ∈ [N ]for i = 1, . . . , N do
Update Ai by using X
Rj (Aj ) + κ−1 Y 
Ai ← Dκ−1 Ri? X −j6=i

end for//STEP 2: updatedual the LagrangianYPNY ← Y + κ X − i=1 Ri (Ai )//STEP 3: update the scalar κκ ← ρκ
Output: (Ai , . . . , AN )components that give rise to the observed tensor. Therefore, for the guarantee of the exact-recovery, we need to takemore structural assumptions of the “true” components into
account.

Exact Recovery with Reshuffled-TD
In this section, we derive and prove the exact-recovery conditions when using the Reshuffled-TD method. We start witha formal statement of the problem.
Problem 1 (Conditions for Exact Recovery) Given a tensor X , suppose there exist low-rank matrices A∗i with rankPN∗
ki such that X =i=1 Ri (Ai ). Under what conditions
on A∗i , Ri and ki , the estimated Âi , obtained by using
Reshuffled-TD, will be equal to A∗i for all i?
Our solution for the problem is stated in Theorem 2. The
main result relies on an incoherence (defined in Definition 3)
which measures the change in the rank of a component A∗iwhen the operation Ri is replaced by any other operator Rj ,i.e. from Ri? (A∗i ) to Rj? (A∗i ) where A∗i := Ri (A∗i ) and R?
denotes the adjoint of R. To be able to measure this change,
we first need to define a low-rank manifold over tensor fora given Ri (see Definition 2), and a neighborhood in thismanifold. For the latter, we will show a type of the tangent
space in this manifold (see Proposition 1). We start with theformal definition of the manifold.
Definition 2 (Tensor manifold under reshuffling) Givena reshuffling operation Ri , the following set of tensors Y
such that the rank of the matrix Y = Ri? (Y) is equal to kiPi := Y ∈ RI1 ×···×IK |rank(Y) = ki .
(7)
defines a smooth manifold (Hosseini, Luke, and Uschmajew2019).
We now define a neighborhood in Pi using a type of tangent
space. In the derivation, the tangent space around a tensor
A ∈ Pi is obtained by the truncated singular-value decomposition of A = Ri? (A) where A ∈ Rmi ×ni is a rank-ki

matrix. Truncated SVD of A with first ki leading singularvalues is given by UΛV> , where U and V are matrices ofsize mi × ki and ni × ki respectively, and Λ is a diagonal
matrix that contains the ki singular values as its diagonal.
By considering all possible real matrices of size mi × ki andni × ki , the tangent space of Pi at the point A = Ri (A) is
given by the following proposition.
Proposition 1 (Tangent Space) Given a rank-ki matrix
A ∈ Rmi ×ni , which generates a tensor A by Ri , i.e.,
A = Ri (A), the tangent space of the manifold Pi at thegiven tensor A is formalized as the following,
Ti (A) := {Y|(8)Ri? (Y) = UV̄> + ŪV> , Ū ∈ Rmi ×ki , V̄ ∈ Rni ×ki.

The proof of the proposition can be trivially achieved
from Eq. (3.2) in (Chandrasekaran et al. 2011). The tangent
space Ti (A) gives us an approximation of the manifold ina neighborhood of A. Due to the relationship A = Ri (A),Ti (A) can be used to analyze how perturbation influences
the rank of A. This is captured in the following incoherence
measure, which we define next.
Definition 3 (Reshuffled-low-rank incoherence)
Consider the tangent space Ti (A) in the manifold Piof (true) rank ki . Given a different operation Rj , we look at
all the tensors Y ∈ Ti (A), and find the maximum spectral
norm kRj? (Y)k2 while kRi? (Y)k2 < 1 for the i’th operator.
The incoherence of a tensor A = Ri (A) is then defined tobe the maximum spectral norm obtained for all operations
Rj 6= Ri . Formally,µi (A) := max
j6=i

max
Y ∈ Ti (A),kRi? (Y)k2 ≤ 1

Rj? (Y) 2 .

(9)

The above incoherence measure captures the change in therank when the operation is changed from Ri to any other Rj .
This is due to a relationship between the spectral norm andthe rank. The spectral norm is the dual of the nuclear norm
which is a convex surrogate of the matrix rank. Roughly
speaking, when the spectral norm under Ri is constrained, asmall spectral norm obtained under Rj would imply a largechange in the rank of the reshuffled matrices. Therefore, asmall value of the incoherence measurement would imply
an increase in the rank when the true operator is replaced bya different one.
Our main result is to show that bounding the incoherence
measurement ensures exact recovery.
Theorem 2 (Exact-Recovery Condition) The estimatedÂi , obtained by Reshuffled-TD, are equal to the true A∗i for
all i, when1
max µi (A∗i ) <,(10)
i=1,...,N
3N − 2where N denotes the number of the components.
The above condition states that if incoherence measurements
are small enough, then exact-recovery is possible. Roughly,this implies that, for exact recovery, the rank of components must increase drastically whenever we switch its corresponding reshuffling operation to other else. From the geometric view, it implies that the tangent spaces Ti (A∗i ), ∀ineed to be well “separated” to each other.
Using Theorem 2 and imposing more constraints on thereshuffling operation Ri , we can get a more intuitive condition for the exact recovery:
PN∗
Corollary 1 Assume that X =i=1 Ri (Ai ) is a Kthorder tensor with the size I × . . . × I, and the the reshuffling
operations Ri : Rn×n → RI×...×I for each component. In
addition, suppose that (a) the rank of A∗i equals r; (b) it isfull-rank for all matrices Rj? (Ri (A∗i )) , ∀j 6= i; (c) For eachmatrix Rj? (Ri (A∗i )) , ∀i, j, its non-zero singular values areequal to each other. Then (A1 , . . . , AN ) = (A∗1 , . . . , A∗N )2is the unique solution of (4) if n > (3N − 2) r.
It implies from Corollary 1 that, the lower bound of the size
n will linearly changed with the rank r, but quadraticallychanged with the number of components N for the exact recovery. Although assumptions in Corollary 1 is quite strict,the result still reveals an intuitive fact that the true components can be more likely exactly recovered by our method if
the data size n is large enough.
Our proof builds upon some of the techniques usedin (Chandrasekaran et al. 2011) to prove similar results for a
type of matrix decomposition. Our proof extends these techniques to the tensor with multiple reshuffled-low-rank structures. Compared to the theoretical studies in (Tomioka andSuzuki 2013), we focus on the conditions for the exact recovery while they mainly analyze the statistical performance
influenced by the perturbation like Gaussian noise. Although
the Theorem 2 and 3 in (Tomioka and Suzuki 2013) showsthat the upper bound for the sum of the reconstruction error of components tends to be tighter with decreasing thestrength of the perturbation, the upper bound is not guaranteed to go to zero even though the strength of the perturbation goes to zero. However, we rigorously prove that the decomposition can exactly recover the latent components, andexplicitly give the incoherence condition on exact recovery
for the first time.

Experimental Results
In the experiments, we specify the reshuffling operations by
uniformly random permutation for simplicity. In practice,
the reshuffling operations can be determined by more practical rules or prior knowledge.

Validation of Exact-Recovery Conditions
We firstly perform an experiment using synthetic data tovalidate the theoretical results in the paper. We generate
data by using N square matrices A∗i ∈ Rn×n , i ∈ [N ].
Each A∗i is generated by multiplying two random semiorthonormal matrices with rank r, i.e., A∗i = Ui Vi> inwhich Ui , Vi ∈ Rn×r denote the random semi-orthonormal
matrices.

We measure the performance using the total signal-tointerference ratio (tSIR) defined as follows:
N
N
X
X
2
2
tSIR :=(11)kA∗i kF /A∗i − Âi
i=1

i=1

F

Fig. 1 shows the phase transition of Reshuffled-TD with different parameters, such as the rank r and size n of the matrices A∗i and the number of components N . In each plot,
the white blocks indicate tSIR ≥ 25dB which implies very
good recovery, and the black blocks indicate tSIR ≤ 15dB
which implies no recovery. The gray area corresponds theresults in between and indicates the phase transition fromexact recovery to partial or no recovery. This can be compared with the theoretical bound given in Corollary 1 whichis shown with the red line. From Corollary 1 we can find, fora fix N , the relationship between n and r is linear, and, whenn is fixed, the relationship between N and r is quadratic.
This matches the relationship shown from experimental results. Our bound is a bit conservative, but correctly captures
a major chunk of the area where exact recovery is possible.
Noise Robustness of Reshuffled-TD. We impose theGaussian noise to evaluate the impact on the performance ofReshuffled-TD. Specifically, we fix the size of the components n = 100, the number of the components N = 10 and
set the rank of each component by r = 1, . . . , 4. Then, we
add the zero-mean i.i.d. Gaussian noise to the data, and thevariance of the noise is controlled by the signal to noise ratio(SNR). Fig. 2 (a) illustrate the performance of ReshuffledTD when 5dB ≤ SN R ≤ 35dB.
As shown in 2 (a), four performance curves are split into
2 groups. we know from Fig. 1 (d) that group 1 corresponds
the rank which satisfies the exact-recovery condition (in thewhite area), while group 2 corresponds the rank whose values do not satisfy the conditions (in the back area). Hence,
the two groups have different trend with the variety of SNR.
In addition, we can find that, in group 1, tSIR is larger
than 20dB when SN R ≥ 20dB. It implies that our method
works smoothly with high SNR. Because we do not consider the noise in our model (it can be seen from the equality
constraint), the performance of our method becomes worse
when SNR is low.
Robustness on the Number of the Components. We
consider the case that the proposed method does not exactly
know how many components are contained in the observation. To simulate this situation, we randomly remove somecomponents from original 10 components, and the removal
probability for each component satisfies the Bernoulli distribution (the mean value equals 0.5). As to the proposed
method, we still assume that all 10 components are contained in the data. To estimate the number of components by
Reshuffled-TD, we compare the norm of the recovered components with a threshold (we choose η = 0.1 for numericalconsideration). Fig. 2 (b) illustrates the estimation accuracy.
Besides the accuracy, the corresponding tSNR performanceis also shown above the bar-plot in the figure.
As shown in Fig. 2 (b), the proposed method is able to reconstruct the contained components with high performanceeven if the the true number of the components is less than expectation. With high SNR value (SNR≥ 20), the estimation

Theoretical bound

3
5

5

5

11

11

13

13

15

20 40 60 80 100 120 140

Theoretical bound
20 40 60 80 100 120 140

Size (n)Size (n)

3

7
9
11

9

9

1

Rank (r)

7

7

15

1

3

Rank (r)Rank (r)31(a) N = 2(b) N = 4

Rank (r)

1

13
15
2

Theoretical bound
4

6

8

10

12

14

No. of the components (N)

5
7
9
11
13
15
2

(c) n = 50

Theoretical bound
4

6

8

10

12

14

No. of the components (N)(d) n = 100

45
40tSIR (dB)35

Rank r=1
Rank r=2
Rank r=3
Rank r=4group 130
25

group 2

20
15
105

10

15

20

25

SNR (dB)

30(a) Gaussian Noise

35

Accuracy percentage (%)

Figure 1: Results on synthetic data to validate our exact-recovery results. We vary different experiment parameters, such as
rank r, size n and number of the components N . In each plot, the darker areas denote the worse reconstruction (tSIR ≤ 15dB)
while the white areas denote a good recovery (tSIR ≥ 25dB). The gray boundary shows the phase transition, while the redline denotes the phase transition predicted by our theoretical bound derived in Corollary 1.
14.92dB 13.87dB 12.01dB 10.89dB
29.83dB 26.48dB 23.68dB 20.78dB

100 43.71dB 40.51dB 35.33dB 25.51dB
80
60
40
SNR=5dB
SNR=20dB
SNR=35dB

20
0

1

2

3

Rank of components

Table 1: Performance comparison of image steganography.
In the experiment, we use SIR (dB) to quantify the distortionof both cover(C) and secret(S) images, where larger valueof SIR indicates better performance. For each row, the SIR
value will be highlighted if it outperforms other methods under a given parameter.

4(b) Components’ number determination

Figure 2: Illustration on robustness of Reshuffled-TD on
non-ideal experimental conditions. Sub-figure (a) shows thetSIR performance of the method which is influenced by
Gaussian noise. In the sub-figure, group 1 corresponds therank which satisfies the exact-recovery condition, while therank in group 2 is not. Sub-figure (b) shows the estimation accuracy of the number of components, and the values shown above the bars are the corresponding tSIR performance of the proposed method.

accuracy of the number of the components achieves 100%,and the accuracy decreases when choosing a large rank. The
high accuracy of our method is due to the fact that the exactrecovery conditions can be still theoretically satisfied as longas we assign an incoherent reshuffling operation Ri , even if
the norms of some components equal zero.

Image Steganography using Reshuffled-TD
Steganography is about concealing a secret message within
an ordinary message and then extracting it at its destination.
In this experiment, we will use Reshuffled-TD for image
steganography, i.e., to conceal a “secret” image in an ordinary “cover” image.
Image steganography is a classic problem for both computer vision and information security. In the existing methods, the most popular one is the least-significant-bits (LSB)
method, which uses the least significant bits of the cover
to hide the most significant bits of the image. In addition,
the similar idea is also extended to transform domains like
Fourier and wavelet (Holub, Fridrich, and Denemark 2014).

LSB

DWT

Datasets

1 bit/
chn

2 bits/chn

2 bits/chn

3 bits/
chn

DPS

Ours

DTD(C)
CART.(S)

26.70
6.92

9.66
14.42

25.17
12.90

23.45
17.81

—
14.04

20.40
21.64

DTD(C)DTD(S)

23.77
3.38

7.53
7.84

22.81
5.27

19.40
9.05

—
3.43

23.69
11.36

DTD(C)FIVEK(S)

24.05
1.12

7.76
6.00

24.70
4.69

22.27
8.48

—
8.97

23.36
10.87

FIVEK(C)FIVEK(S)

23.02
3.37

6.56
7.52

21.54
5.48

18.57
8.74

—
8.96

21.86
6.67

FVC(C)FIVEK(S)

18.19
3.32

3.27
6.42

24.47
4.84

19.95
8.30

—
8.90

20.25
12.80

LIVE(C)FIVEK(S)

24.50
4.08

7.66
7.58

24.49
5.32

20.93
9.46

—
9.50

24.71
11.49

Some recent approaches have used deep neural networks
to hide and recover images (Baluja 2017), but these methods
require lots of training data, and they are generally sensitive
to the images not present in the training data. The computational requirement is also heavy. In contrast, our method ismuch simpler. It does not require any training, and therefore
does not have any such sensitivity issues.
We tried various ways to make the problem challenging
for the methods. We try to conceal a full-size RGB image(8 × 3 bits per pixel) into a grayscale image (8 bits per
pixel). Meanwhile, we choose different types of images forsteganography, e.g., natural, cartoon and fingerprint.
The datasets we used in the experiment include texture (DTD), natural (LIVE and FIVEK (Bychkovsky et al.
2011)), cartoon (CART. (Royer et al. 2017)) and fingerprint(FVC (Maltoni et al. 2009)) datasets. For different datasets,
we unify the shape of all images to 2000×2000, and convert
the image to grayscale when the cover image is colored.

Figure 3: An example to illustrate the performance of image steganography by different methods. In the figure, the first row
shows the original cover images (the first column) and the container images generated by different methods; The second row
shows the original secret images (the first column) and its recovery by different methods.
A sketch of our Reshuffled-TD method is shown in supplemental materials. During the concealing phase, we consider each channel of the secret image as one component, and they are randomly reshuffled. Then, we added thereshuffled “components” to the cover image to obtain the agrayscale “container” image (a.k.a. observation). The difference between the container and cover images will tend tozero as we decrease the strength of the secret components
by multiplying a scalar. Therefore, we expect that the secretimage can be hidden well if we choose a appropriate valueof this “strength” scalar. In the recovery phase, we use thereshuffling operations as key, and recover the RGB components of the secret image by Reshuffled-TD.
Experimental results are shown in Table 1 as measured
by the signal to interference ratio (SIR). A higher valueof SIR indicate better performance. The experiment is conducted on 10 randomly chosen image pairs. We compare to
three existing methods: (a) the LSB method (b) the discrete
wavelet transform based method (DWT), and deep stego(DPS) (Baluja 2017). Because DPS converts the grayscale
cover image into a RGB image as the container, we just show
the SIR for secret image in the table.
As shown in Table 1, Reshuffled-TD significantly outperforms all the state-of-the-art methods in the experiment. For
example, in FVC+FIVEK dataset, Reshuffled-TD achieves
20.25dB on the cover images and 12.80dB on the secret images. With the similar SIR on the cover image, LSB only
achieves 3.32dB. Meanwhile, DWT and DPS achieve 8.30
and 8.90dB, respectively. As the worst performance, the
FIVEK+FIVEK dataset shows an exception in the experiment. This is because FIVEK is a dataset of natural images,
which contains full of detail information. Hence, there is lessroom in the cover images to hide additional information. Fig.
3 shows two examples of reconstructed images obtained in
the experiment. More examples for visual comparison are
shown in the supplementary material.

Conclusion
By leveraging the (generalized) latent convex tensor decomposition, a.k.a. Reshuffled-TD, we rigorously proved thatthe low-rank components an be exactly recovered when theincoherence is sufficiently upper-bounded. In addition, weapplied the generalized model to a totally novel task, i.e.

image steganography. Experimental results on various realworld datasets demonstrate that Reshuffled-TD outperforms
both the classic state-of-the-arts but also deep-learningbased methods.
As potential works in the future, first we consider to takethe low-tensor-rank assumption into the model instead of thecurrent multi-linear rank. Second, it would be an interesting topic to seek for the “optimal” reshuffling operations if
there exist training data. If such operations are learnable, it
implies hat we might find lower-rank representation for thedata, such that the well-developed low-rank-based methods
can be employed on the transformed data.

Acknowledgments
This work was partially supported by JSPS KAKENHI(Grant No. 17K00326) and NSFC (Grant No. 61673124).

References
Baluja, S. 2017. Hiding images in plain sight: Deepsteganography. In Advances in Neural Information Processing Systems, 2069–2079.
Becker, H.; Albera, L.; Comon, P.; Haardt, M.; Birot, G.;
Wendling, F.; Gavaret, M.; Bénar, C.-G.; and Merlet, I.
2014. EEG extended source localization: tensor-based vs.
conventional methods. NeuroImage 96:143–157.
Bychkovsky, V.; Paris, S.; Chan, E.; and Durand, F. 2011.
Learning photographic global tonal adjustment with a
database of input / output image pairs. In The Twenty-FourthIEEE Conference on Computer Vision and Pattern Recognition.
Cai, J.-F.; Candès, E. J.; and Shen, Z. 2010. A singularvalue thresholding algorithm for matrix completion. SIAM
Journal on Optimization 20(4):1956–1982.
Candès, E. J., and Recht, B. 2009. Exact matrix completion via convex optimization. Foundations of Computational
mathematics 9(6):717.
Chandrasekaran, V.; Sanghavi, S.; Parrilo, P. A.; and Willsky, A. S. 2011. Rank-sparsity incoherence for matrix
decomposition. SIAM Journal on Optimization 21(2):572–
596.

Cichocki, A.; Zdunek, R.; Phan, A. H.; and Amari, S.-i.
2009. Nonnegative matrix and tensor factorizations: applications to exploratory multi-way data analysis and blindsource separation. John Wiley & Sons.
Comon, P.; Luciani, X.; and De Almeida, A. L. 2009. Tensor decompositions, alternating least squares and other tales.
Journal of Chemometrics: A Journal of the Chemometrics
Society 23(7-8):393–405.
De Lathauwer, L. 2008. Decompositions of a higher - order tensor in block terms, part II: Definitions and uniqueness. SIAM Journal on Matrix Analysis and Applications
30(3):1033–1066.
Fazel, M.; Hindi, H.; and Boyd, S. P. 2001. A rank minimization heuristic with application to minimum order system approximation. In American Control Conference, 2001.
Proceedings of the 2001, volume 6, 4734–4739. IEEE.
González, D.; Ammar, A.; Chinesta, F.; and Cueto, E. 2010.
Recent advances on the use of separated representations. International Journal for Numerical Methods in Engineering
81(5):637–659.
Guo, X.; Yao, Q.; and Kwok, J. T.-Y. 2017. Efficientsparse low-rank tensor completion using the Frank-Wolfe algorithm. In AAAI, 1948–1954.
He, L.; Lu, C.-T.; Ma, G.; Wang, S.; Shen, L.; Yu, P. S.; andRagin, A. B. 2017. Kernelized support tensor machines.
In Proceedings of the 34th International Conference on Machine Learning-Volume 70, 1442–1451. JMLR. org.
Holub, V.; Fridrich, J.; and Denemark, T. 2014. Universal
distortion function for steganography in an arbitrary domain.
EURASIP Journal on Information Security 2014(1):1.
Hosseini, S.; Luke, D. R.; and Uschmajew, A. 2019. Tangent and normal cones for low-rank matrices. In Nonsmooth
Optimization and Its Applications. Springer. 45–53.
Imaizumi, M.; Maehara, T.; and Hayashi, K. 2017. On tensor train rank minimization: Statistical efficiency and scalable algorithm. In Advances in Neural Information Processing Systems, 3933–3942.
Kolda, T. G., and Bader, B. W. 2009. Tensor decompositions
and applications. SIAM review 51(3):455–500.
Li, C.; He, W.; Yuan, L.; Sun, Z.; and Zhao, Q. 2019. Guaranteed matrix completion under multiple linear transformations. In Proceedings of the IEEE Conference on Computer
Vision and Pattern Recognition, 11136–11145.
Lin, Z.; Chen, M.; and Ma, Y. 2010. The augmented lagrange multiplier method for exact recovery of corrupted
low-rank matrices. arXiv preprint arXiv:1009.5055.
Lu, C.; Peng, X.; and Wei, Y. 2019. Low-rank tensor completion with a new tensor nuclear norm induced by invertible
linear transforms. In Proceedings of the IEEE Conference on
Computer Vision and Pattern Recognition, 5996–6004.
Maltoni, D.; Maio, D.; Jain, A. K.; and Prabhakar, S. 2009.
Handbook of fingerprint recognition. Springer Science &Business Media.
Mu, C.; Huang, B.; Wright, J.; and Goldfarb, D. 2014.
Square deal: Lower bounds and improved relaxations fortensor recovery. In International Conference on Machine
Learning, 73–81.
Nimishakavi, M.; Jawanpuria, P. K.; and Mishra, B. 2018.
A dual framework for low-rank tensor completion. In Advances in Neural Information Processing Systems, 5484–
5495.
Oseledets, I. V. 2011. Tensor-train decomposition. SIAM
Journal on Scientific Computing 33(5):2295–2317.
Rabusseau, G., and Kadri, H. 2016. Low-rank regression
with tensor responses. In Advances in Neural Information
Processing Systems, 1867–1875.
Royer, A.; Bousmalis, K.; Gouws, S.; Bertsch, F.; Mosseri,
I.; Cole, F.; and Murphy, K. 2017. Xgan: Unsupervised image-to-image translation for many-to-many mappings. arXiv preprint arXiv:1711.05139.
Sharan, V., and Valiant, G. 2017. Orthogonalized ALS: A
theoretically principled tensor decomposition algorithm forpractical use. arXiv preprint arXiv:1703.01804.
Song, L.; Ishteva, M.; Parikh, A.; Xing, E.; and Park, H.
2013. Hierarchical tensor decomposition of latent tree
graphical models. In International Conference on Machine
Learning, 334–342.
Tomioka, R., and Suzuki, T. 2013. Convex tensor decomposition via structured Schatten norm regularization. In
Advances in neural information processing systems, 1331–
1339.
Tomioka, R.; Hayashi, K.; and Kashima, H. 2010. Estimation of low-rank tensors via convex optimization. NIPS2010
workshop “Tensors, Kernels, and Machine Learning”.
Wang, A.; Song, X.; Wu, X.; Lai, Z.; and Jin, Z. 2019. Latent
schatten tt norm for tensor completion. In ICASSP 20192019 IEEE International Conference on Acoustics, Speechand Signal Processing (ICASSP), 2922–2926. IEEE.
Williams, A. H.; Kim, T. H.; Wang, F.; Vyas, S.; Ryu, S. I.;
Shenoy, K. V.; Schnitzer, M.; Kolda, T. G.; and Ganguli, S.
2018. Unsupervised discovery of demixed, low-dimensional
neural dynamics across multiple timescales through tensor
component analysis. Neuron 98(6):1099–1115.
Wimalawarne, K.; Yamada, M.; and Mamitsuka, H. 2017.
Convex coupled matrix and tensor completion. arXiv
preprint arXiv:1705.05197.
Yamada, M.; Lian, W.; Goyal, A.; Chen, J.; Wimalawarne,
K.; Khan, S. A.; Kaski, S.; Mamitsuka, H.; and Chang, Y.
2017. Convex factorization machine for toxicogenomics
prediction. In Proceedings of the 23rd ACM SIGKDD International Conference on Knowledge Discovery and Data
Mining, 1215–1224. ACM.
Yu, J.; Li, C.; Zhao, Q.; and Zhao, G. 2019. Tensorring nuclear norm minimization and application for visual:
Data completion. In ICASSP 2019-2019 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 3142–3146. IEEE.
Zhao, Q.; Zhou, G.; Xie, S.; Zhang, L.; and Cichocki,
A. 2016. Tensor ring decomposition. arXiv preprint
arXiv:1606.05535.