A New Convex Relaxation for Tensor Completion

Introduction

During the recent years, there has been a growing interest on the problem of learning a tensor from a set of linear measurements, such as a subset of its entries, see and references therein. This methodology, which is also referred to as tensor completion, has been applied to various fields, ranging from collaborative filtering , to computer vision , to medical imaging , among others. In this paper, we propose a new method to tensor completion, which is based on a convex regularizer which encourages low rank tensors and develop an algorithm for solving the associated regularization problem.

Arguably the most widely used convex approach to tensor completion is based upon the extension of trace norm regularization to that context. This involves computing the average of the trace norm of each matricization of the tensor . A key insight behind using trace norm regularization for matrix completion is that this norm provides a tight convex relaxation of the rank of a matrix defined on the spectral unit ball . Unfortunately, the extension of this methodology to the more general tensor setting presents some difficulties. In particular, we shall prove in this paper that the tensor trace norm is not a tight convex relaxation of the tensor rank.

The above negative result stems from the fact that the spectral norm, used to compute the convex relaxation for the trace norm, is not an invariant property of the matricization of a tensor. This observation leads us to take a different route and study afresh the convex relaxation of tensor rank on the Euclidean ball. We show that this relaxation is tighter than the tensor trace norm, and we describe a technique to solve the associated regularization problem. This method builds upon the alternating direction method of multipliers and a subgradient method to compute the proximity operator of the proposed regularizer. Furthermore, we present numerical experiments on one synthetic dataset and two real-life datasets, which indicate that the proposed method improves significantly over tensor trace norm regularization in terms of estimation error, while remaining computationally tractable.

The paper is organized in the following manner. In Section 2, we describe the tensor completion framework. In Section 3, we highlight some limitations of the tensor trace norm regularizer and present an alternative convex relaxation for the tensor rank. In Section 4, we describe a method to solve the associated regularization problem. In Section 5, we report on our numerical experience with the proposed method. Finally, in Section 6, we summarize the main contributions of this paper and discuss future directions of research.

Preliminaries

where $\gamma$ is a positive parameter which may be chosen by cross validation. The role of the regularizer $R$ is to encourage tensors ${\boldsymbol{\cal W}}$ which have a simple structure in the sense that they involve a small number of “degrees of freedom”. A natural choice is to consider the average of the rank of the tensor’s matricizations. Specifically, we consider the combinatorial regularizer

Finding a convex relaxation of this regularizer has been the subject of recent works . They all agree to use the trace norm for tensors as a convex proxy of $R$. This is defined as the average of the trace norm of each matricization of ${\boldsymbol{\cal W}}$, that is,

A rational behind the regularizer (3) is that the trace norm is the tightest convex lower bound to the rank of a matrix on the spectral unit ball, see [8, Thm. 1]. This lower bound is given by the convex envelope of the function

where $\|\cdot\|_{\infty}$ is the spectral norm, namely the largest singular value of $W$. The convex envelope can be derived by computing the double conjugate of $\Psi$. This is defined as

We refer to for a detailed discussion of these ideas. We will use this equivalence between spectral and gauge functions repeatedly in the paper.

Alternative Convex Relaxation

In this section, we show that the tensor trace norm is not a tight convex relaxation of the tensor rank $R$ in equation (2). We then propose an alternative convex relaxation for this function.

Note that due to the composite nature of the function $R$, computing its convex envelope is a challenging task and one needs to resort to approximations. In , the authors note that the tensor trace norm $\|\cdot\|_{\rm tr}$ in equation (3) is a convex lower bound to $R$ on the set

The key insight behind this observation is summarized in Lemma 4, which we report in Appendix A. However, the authors of leave open the question of whether the tensor trace norm is the convex envelope of $R$ on the set $\mathcal{G}_{{\rm\infty}}$. In the following, we will prove that this question has a negative answer by showing that there exists a convex function $\Omega\neq\|\cdot\|_{\rm tr}$ which underestimates the function $R$ on $\mathcal{G}_{{\rm\infty}}$ and such that for some tensor ${\boldsymbol{\cal W}}\in\mathcal{G}_{{\rm\infty}}$ it holds that $\Omega({\boldsymbol{\cal W}})>\|{\boldsymbol{\cal W}}\|_{\rm tr}$.

To describe our observation we introduce the set

where $\|\cdot\|_{2}$ is the Euclidean norm for tensors, that is,

In particular, if $s=kx$ for some $k>0$ this inequality becomes

If $k$ is large enough, the maximum is attained at $r={\rm card}(x)$. Consequently,

By the definition of the convex envelope, it also holds that $\omega_{\alpha}^{**}(x)\leq{\rm card}(x)$. The result follows. ∎

The next lemma provides, together with Lemma 1, a sufficient condition for the existence of a tensor ${\boldsymbol{\cal W}}\in{\cal G}_{\rm\infty}$ at which the proposed regularizer is strictly larger than the tensor trace norm.

Without loss of generality we assume that $p_{1}\leq\cdots\leq p_{N}$. By hypothesis $p_{1}<p_{N}$. First we consider the special case

We define a class of tensors ${\boldsymbol{\cal W}}$ by choosing a singular value decomposition for their mode-$N$ matricization,

The case $n=N$ is immediate. If $n=1$ we have

We are now ready to present the main result of this section.

By construction $\Omega_{\alpha}({\boldsymbol{\cal W}})\leq R({\boldsymbol{\cal W}})$ for every ${\boldsymbol{\cal W}}\in\alpha{\cal G}_{2}$. Since ${\cal G}_{{\rm\infty}}\subset\alpha\mathcal{G}_{{\rm 2}}$ then $\Omega_{\alpha}$ is a convex lower bound for the tensor rank $R$ on the set $\mathcal{G}_{{\rm\infty}}$ as well. The first claim now follows by Lemmas 1 and 2. Indeed, all tensors obtained following the process described in Lemma 2 have the property that

Furthermore there are infinitely many such tensors which satisfy this claim since the left singular vectors can be arbitrarily chosen in equation (9). To prove the second claim, we note that since $\omega_{1}^{**}$ is the convex envelope of the cardinality ${\rm card}$ on the Euclidean unit ball, then it holds that $\omega_{1}^{**}(\sigma)\geq\|\sigma\|_{1}$ for every vector $\sigma$ such that $\|\sigma\|_{2}\leq 1$. Consequently,

The above result stems from the fact that the spectral norm is not an invariant property of the matricization of a tensor, whereas the Euclidean (Frobenius) norm is. This observation leads us to further study the function $\Omega_{\alpha}$.

Optimization Method

In this section, we explain how to solve the regularization problem associated with the proposed regularizer (6). For this purpose, we first recall the alternating direction method of multipliers (ADMM) , which was conveniently applied to tensor trace norm regularization in .

To explain ADMM we consider a more general problem comprising both tensor trace norm regularization and the regularizer we propose,

where $E({\boldsymbol{\cal W}})$ is an error term such as $\|y-{\cal I}({\boldsymbol{\cal W}})\|^{2}_{2}$ and $\Psi$ is a convex spectral function. It is defined, for every matrix $A$, as

The corresponding augmented Lagrangian (see e.g. ) is given by

ADMM is based on the following iterative scheme

Step (16) is straightforward, whereas step (14) is described in . Here we focus on the step (15) since this is the only problem which involves function $\Psi$. We restate it with more explanatory notations as

By completing the square in the right hand side, the solution of this problem is given by

where $X=W_{\left(n\right)}-\frac{1}{\beta}A_{n\left(n\right)}$. By using properties of proximity operators (see e.g. [2, Prop. 3.1]) we know that if $\psi$ is a gauge function then

where $U_{X}$ and $V_{X}$ are the orthogonal matrices formed by the left and right singular vectors of $X$, respectively.

If we choose $\psi=\left\|\cdot\right\|_{1}$ the associated proximity operator is the well-known soft thresholding operator, that is, ${\rm prox}_{\frac{1}{\beta}\left\|\cdot\right\|_{1}}\left(\sigma\right)=v$, where the vector $v$ has components

On the other hand, if we choose $\psi=\omega^{**}_{\alpha}$, we need to compute ${\rm prox}_{\frac{1}{\beta}\omega^{**}_{\alpha}}$. In the next section, we describe a method to accomplish this task.

2 Computation of the Proximity Operator

In order to compute the proximity operator of the function ${\frac{1}{\beta}\omega_{\alpha}^{**}}$ we will use several properties of proximity calculus. First, we use the formula (see e.g. ) ${\rm prox}_{g^{*}}\left(x\right)=x-{\rm prox}_{g}\left(x\right)$ for $g^{*}=\frac{1}{\beta}\omega_{\alpha}^{**}$. Next we use a property of conjugate functions from , which states that $g(\cdot)=\frac{1}{\beta}\omega_{\alpha}^{*}(\beta\cdot)$. Finally, by the scaling property of proximity operators , we have that ${\rm prox}_{g}\left(x\right)=\frac{1}{\beta}{\rm prox}_{\beta\omega_{\alpha}^{*}}\left(\beta x\right)$.

It remains to compute the proximity operator of a multiple of the function $\omega_{\alpha}^{*}$ in equation (7), that is, for any $\beta>0$, $y\in\mathcal{S}$, we wish to compute

In order to solve this problem we employ the projected subgradient method, see e.g. . It consists in applying two steps at each iteration. First, it advances along a negative subgradient of the current solution; second, it projects the resultant point onto the feasible set $\mathcal{S}$. In fact, according to , it is sufficient to compute an approximate projection, a step which we describe in Appendix B. To compute a subgradient of $h$ at $w$, we first find any integer ${k}$ such that ${{k}}\in\underset{r=0}{\overset{d}{{\rm argmax}}}\left\{\alpha\left\|w_{1:r}\right\|_{2}-r\right\}$. Then, we calculate a subgradient ${g}$ of the function $h$ at $w$ by the formula

Now we have all the ingredients to apply the projected subgradient method, which is summarized in Algorithm 1. In our implementation we stop the algorithm when an update of ${\hat{w}}$ is not made for more than $10^{3}$ iterations.

Experiments

We have conducted a set of experiments to assess whether there is any advantage of using the proposed regularizer over the tensor trace norm for tensor completion. First, we have designed a synthetic experiment to evaluate the performance of both approaches under controlled conditions. Then, we have tried both methods on two tensor completion real data problems. In all cases, we have used a validation procedure to tune the hyper-parameter $\gamma$, present in both approaches, among the values $\left\{10^{j}\,:\,j=-7,-6,\ldots,0\right\}$. In our proposed approach there is one further hyper-parameter, $\alpha$, to be specified. It should take the value of the Frobenius norm of any matricization of the underlying tensor. Since this is unknown, we propose to use the estimate

where each entry of the Tucker decomposition components is sampled from the standard Gaussian distribution $\mathcal{N}(0,1)$. We then created the ground truth tensor ${\boldsymbol{\cal W}}^{0}$ by the equation

where ${\rm mean}({\boldsymbol{\cal W}})$ and ${\rm std}({\boldsymbol{\cal W}})$ are the mean and standard deviation of the elements of ${{\boldsymbol{\cal W}}}$ and the $\xi_{i_{1},i_{2},i_{3}}$ are i.i.d. Gaussian random variables with zero mean and variance $\sigma^{2}$. We have randomly sampled $10\%$ of the elements of the tensor to compose the training set, $45\%$ for the validation set, and the remaining $45\%$ for the test set. After repeating this process $20$ times, we report the average results in Figure 1 (Left). Having conducted a paired $t$-test for each value of $\sigma^{2}$, we conclude that the visible differences in the performances are highly significant, obtaining always $p$-values less than $0.01$ for $\sigma^{2}\leq 10^{-2}$.

2 School Dataset

We have selected randomly $5\%$ of the instances to make the test set and another $5\%$ of the instances for the validation set. From the remaining instances, we have randomly chosen $m$ of them for several values of $m$. This procedure has been repeated $20$ times and the average performance is presented in Figure 2 (Left). There is a distinguishable improvement of our approach with respect to tensor trace norm regularization. To check whether this gap is significant, we have conducted a set of paired $t$-tests for each value of $m$. In all cases we obtained a $p$-value below $0.01$.

3 Video Completion

Conclusion

In this paper, we proposed a convex relaxation for the average of the rank of the matricizations of a tensor. We compared this relaxation to a commonly used convex relaxation used in the context of tensor completion, which is based on the trace norm. We proved that this second relaxation is not tight and argued that the proposed convex regularizer may be advantageous. Empirical comparisons indicate that our method consistently improves in terms of estimation error over tensor trace norm regularization, while being computationally comparable on the range of problems we considered. In the future it would be interesting to study methods to speed up the computation of the proximity operator of our regularizer and investigate its utility in tensor learning problem beyond tensor completion such as multilinear multitask learning .

References

Appendix

In this appendix, we describe an auxiliary result and present the main steps for the computation of an approximate projection.

Appendix A A Useful Lemma

Since the restriction of $g$ on ${\cal D}^{N}\subseteq\prod_{n=1}^{N}{\cal C}_{n}$ equals to $h$, the convex envelope of $g$ when evaluated on the smaller set ${\cal D}^{N}$ cannot be larger than the convex envelope of $h$ on ${\cal D}$. ∎

Using this result it is immediately possible to derive a convex lower bound for the function $R$ in equation (2). Since the convex envelope of the rank function on the unit ball of the spectral norm is the trace norm, using Lemma 4 with ${\cal C}_{n}=\{{\boldsymbol{\cal W}}:\|W_{(n)}\|_{\infty}\leq 1\}$ and

we conclude that the convex envelope of the function $R$ on the set ${\cal G}_{{\rm\infty}}$ is bounded from below by $\frac{1}{N}\sum_{n=1}^{N}\|W_{(n)}\|_{\rm tr}$. Likewise the convex envelope of $R$ on the set $\alpha{\cal G}_{{\rm 2}}$ is lower bounded by the function $\Omega_{\alpha}$ in equation (6).

Appendix B Computation of an Approximated Projection

Here, we address the issue of computing an approximate Euclidean projection onto the set

$\mathcal{S}=\mathcal{S}_{d}\subseteq\mathcal{S}_{d-1}\subseteq\ldots\subseteq\mathcal{S}_{1}$,

$\left\|P_{\mathcal{S}_{i}}\left(\cdots P_{\mathcal{S}_{1}}\left(v\right)\right)-z\right\|_{2}\leq\left\|v-z\right\|_{2}$

$\left\|P_{\mathcal{S}_{i+1}}\left(P_{\mathcal{S}_{i}}\left(\cdots P_{\mathcal{S}_{1}}\left(v\right)\right)\right)-z\right\|_{2}\leq\left\|P_{\mathcal{S}_{i}}\left(\cdots P_{\mathcal{S}_{1}}\left(v\right)\right)-z\right\|_{2}\leq\left\|v-z\right\|_{2}$,

since $z$ is also contained in $\mathcal{S}_{i+1}$.

Note that to evaluate the right hand side of equation (18) we do not require full knowledge of $P_{{\cal S}_{i}}$, we only need to compute $P_{\mathcal{S}_{i+1}}(v)$ for $v\in\mathcal{S}_{i}$. The next proposition describes a recursive formula to achieve this step.

For any $v\in\mathcal{S}_{i}$, we express its first $i$ elements as $v_{1:i}=\left[v_{1:i-j},\,v_{i}\boldsymbol{1}^{j}\right]$, where the last $j\in[i]$ is the largest integer such that $v_{i-j+1}=v_{i-j+2}=\cdots=v_{i}$. It holds that

The first case is straightforward. In the following we prove the remaining two. In both cases it will be useful to recall that the projection operator $P_{\mathcal{C}}$ on any convex set $\mathcal{C}$ is characterized as

To prove the second case, we use property (19) and apply simple algebraic transformations to obtain, for all $z\in\mathcal{S}_{i+1}$, that

Finally we prove the third case. We want to show that if $x=P_{\mathcal{S}_{i+1}}(v)$ then

By using property (19), the last equation is equivalent to the statement that if

A way to show that it holds true is to prove that the term in the left hand side of (21) is upper bounded by the corresponding term in (20). That is, for every $z\in\mathcal{S}_{i+1}$, we want to show that

A direct computation yields the equivalent inequality

Since $x=P_{\mathcal{S}_{i+1}}\left(v\right)$, $v_{i-j+1}=v_{i-j+2}=\cdots=v_{i}$ and $v_{i+1}>v_{i}$, then $x_{i-j+1}=x_{i-j+2}=\cdots=x_{i+1}$. Consequently, the left hand side of inequality (22) is equivalent to

Note that the first factor is negative and the second is positive because $z$ and $v$ are in $\mathcal{S}_{i+1}$. The result follows. ∎

Algorithm 2 summarizes our method to compute the approximated projection operator onto the set $\mathcal{S}$, based on Proposition 5.