Tensor p-shrinkage nuclear norm for low-rank tensor completion

TL;DR

This paper introduces a tensor p-shrinkage nuclear norm (p-TNN) based on tensor SVD, which better approximates tensor rank for low-rank tensor completion, along with an efficient algorithm and theoretical guarantees.

Contribution

It proposes a novel tensor p-shrinkage nuclear norm and a low-rank tensor completion model with convergence guarantees, outperforming existing methods.

Findings

p-TNN provides a better approximation of tensor rank than traditional nuclear norm.

The proposed algorithm converges globally under certain conditions.

Numerical results show superior performance on synthetic and real data.

Abstract

In this paper, a new definition of tensor p-shrinkage nuclear norm (p-TNN) is proposed based on tensor singular value decomposition (t-SVD). In particular, it can be proved that p-TNN is a better approximation of the tensor average rank than the tensor nuclear norm when p < 1. Therefore, by employing the p-shrinkage nuclear norm, a novel low-rank tensor completion (LRTC) model is proposed to estimate a tensor from its partial observations. Statistically, the upper bound of recovery error is provided for the LRTC model. Furthermore, an efficient algorithm, accelerated by the adaptive momentum scheme, is developed to solve the resulting nonconvex optimization problem. It can be further guaranteed that the algorithm enjoys a global convergence rate under the smoothness assumption. Numerical experiments conducted on both synthetic and real-world data sets verify our results and demonstrate the superiority of our p-TNN in LRTC problems over several state-of-the-art methods.

Tables

Table 1. Table 1: Comparison of the recovered results on color image data sets with different sampling rates. The optimal values of PSNR, RSE, and SSIM for each test are highlighted in bold.

sr	method	lena			facade			baboon			house
		PSNR	RSE	SSIM	PSNR	RSE	SSIM	PSNR	RSE	SSIM	PSNR	RSE	SSIM
0.1	SPC	25.5846	0.0948	0.7742	26.4615	0.0921	0.8189	21.0733	0.1638	0.5082	25.8000	0.0842	0.7979
	STDC	19.7222	0.1862	0.4401	20.0552	0.1925	0.6117	16.9850	0.2622	0.3892	19.7105	0.1697	0.4099
	Tmac-TT	23.2412	0.1242	0.6972	19.7929	0.1985	0.6098	18.8118	0.2096	0.4046	23.1784	0.1145	0.6474
	TT-SGD	15.3518	0.2282	0.4245	20.8763	0.1731	0.6134	16.0909	0.2946	0.3520	16.8428	0.2448	0.3781
	TNN	24.1932	0.1080	0.7105	22.1981	0.1505	0.6571	20.6019	0.1729	0.4156	22.3312	0.1255	0.6060
	proposed	25.9028	0.0916	0.8015	24.3383	0.1180	0.7328	21.3046	0.1595	0.5145	25.9565	0.0828	0.8152
0.2	SPC	27.4107	0.0768	0.8474	26.7079	0.0891	0.8740	22.2023	0.1438	0.6126	27.5046	0.0692	0.8480
	STDC	25.8928	0.0915	0.8066	24.4484	0.1161	0.7727	19.6754	0.1924	0.4094	28.1505	0.0642	0.8673
	Tmac-TT	25.8833	0.0916	0.8015	22.0045	0.1538	0.6506	21.5952	0.1542	0.5884	25.9585	0.0826	0.8230
	TT-SGD	20.8364	0.1677	0.6224	25.2864	0.1067	0.7953	18.7745	0.2125	0.3959	21.4953	0.1363	0.4762
	TNN	26.0373	0.0900	0.8151	26.4776	0.0919	0.8213	21.9099	0.1487	0.6032	25.5434	0.0867	0.7918
	proposed	28.4658	0.0681	0.8820	27.9356	0.0783	0.8790	22.7415	0.1352	0.6565	29.3077	0.0566	0.8894
0.3	SPC	28.7736	0.0657	0.8871	28.6894	0.0713	0.8927	23.2435	0.1276	0.7150	29.0217	0.0581	0.8828
	STDC	27.9464	0.0722	0.8553	25.6342	0.1013	0.8092	21.7230	0.1520	0.5933	30.4959	0.0490	0.9129
	Tmac-TT	26.8824	0.0816	0.8341	23.6913	0.1267	0.7285	22.8057	0.1342	0.6657	26.9618	0.0736	0.8375
	TT-SGD	24.0597	0.1113	0.7066	28.1859	0.0750	0.8900	20.7091	0.1709	0.4982	24.1978	0.1010	0.6512
	TNN	28.0295	0.0715	0.8649	28.8579	0.0699	0.9100	22.9498	0.1320	0.6867	27.8168	0.0667	0.8652
	proposed	30.4150	0.0545	0.9324	30.1530	0.0602	0.9270	24.0340	0.1165	0.7626	32.1177	0.0407	0.9225
0.4	SPC	30.0387	0.0568	0.9147	29.0100	0.0687	0.9142	24.2049	0.1142	0.7897	29.8044	0.0531	0.9070
	STDC	29.1445	0.0629	0.8960	26.6845	0.0898	0.8721	23.4226	0.1250	0.7475	31.5988	0.0432	0.9180
	Tmac-TT	27.8189	0.0733	0.8511	25.5229	0.1026	0.8081	23.7994	0.1197	0.7497	27.9107	0.0660	0.8663
	TT-SGD	26.6475	0.0856	0.8272	29.6810	0.0619	0.9248	22.3022	0.1412	0.6204	25.8701	0.0831	0.8015
	TNN	29.8159	0.0582	0.9005	30.4659	0.0581	0.9337	23.9500	0.1176	0.7580	29.4037	0.0556	0.8934
	proposed	32.5232	0.0427	0.9582	31.3102	0.0527	0.9466	25.3032	0.1006	0.8308	33.6028	0.0343	0.9417

Table 2. Table 2: Comparison of the recovered results on Hyperspectral image data sets with different sampling rates. The optimal values of PSNR and RSE for each test are highlighted in bold.

sr	method	WDC Mall		PaviaU
		PSNR	RSE	PSNR	RSE
0.1	SPC	22.9355	0.1926	26.1087	0.2270
	STDC	16.1068	0.4221	21.4573	0.3877
	Tmac-TT	16.5548	0.4009	21.6522	0.3791
	TT-SGD	17.4878	0.3606	21.0885	0.4045
	TNN	19.6377	0.2813	23.2825	0.3142
	proposed	22.9442	0.1924	27.3111	0.1976
0.2	SPC	26.8823	0.1222	31.2002	0.1263
	STDC	20.8547	0.2447	24.2979	0.2796
	Tmac-TT	18.7733	0.3105	23.4255	0.3091
	TT-SGD	20.0943	0.2668	24.7403	0.2657
	TNN	21.1869	0.2353	24.9034	0.2608
	proposed	27.6564	0.1118	31.4136	0.1234
0.3	SPC	29.8279	0.0871	34.2700	0.0887
	STDC	22.8627	0.1941	26.4978	0.2171
	Tmac-TT	20.5384	0.2534	25.0036	0.2577
	TT-SGD	21.7456	0.2206	26.2153	0.2242
	TNN	24.5025	0.1608	27.7402	0.1881
	proposed	30.0019	0.0854	34.4901	0.0865
0.4	SPC	31.1864	0.0745	35.6849	0.0754
	STDC	24.2289	0.1658	28.0926	0.1807
	Tmac-TT	22.0513	0.2129	26.6126	0.2142
	TT-SGD	23.2092	0.1863	27.5421	0.1924
	TNN	26.6422	0.1257	29.6721	0.1507
	proposed	32.1133	0.0669	36.9243	0.0653

Equations

$$\underset{\mathcal{X}}{\min }rank\left(\mathcal{X}\right)s.t.{\mathcal{P}}_{\Omega }\left(\mathcal{X}\right)={\mathcal{P}}_{\Omega }\left(\mathcal{T}\right).$$

$$\underset{\mathcal{X}}{\min }\ell \left(\mathcal{X}\right)+\lambda f\left(\mathcal{X}\right).$$

$${\Vert \mathcal{X}\Vert }_{CNN}=\inf \left\{\sum _{i=1}^r\left|{\lambda }_i\right||\mathcal{X}=\sum _{i=1}^r{\lambda }_i{\mathbf{u}}_{1,i}\circ {\mathbf{u}}_{2,i}\circ \cdots \circ {\mathbf{u}}_{n,i}\right\}.$$

$${\Vert \mathcal{X}\Vert }_{SNN}=\sum _{i=1}^n{\Vert {\mathbf{X}}_{\left(i\right)}\Vert }_{\ast }.$$

$${\Vert \mathcal{X}\Vert }_{TTNN}=\sum _{i=1}^{n-1}{\Vert {\mathbf{X}}_{\left[i\right]}\Vert }_{\ast }.$$

$$unfold\left({\cal X}\right)=\left[{\begin{array}[]{*{20}{c}}{{{\mathbf{X}}^{\left(1\right)}}}\\ {{{\mathbf{X}}^{\left(2\right)}}}\\ \vdots\\ {{{\mathbf{X}}^{\left({{I_{3}}}\right)}}}\end{array}}\right],fold\left({unfold\left({\cal X}\right)}\right)={\cal X}.$$

$$bcirc\left({\cal X}\right)=\left[{\begin{array}[]{*{20}{c}}{{\mathbf{X}^{\left(1\right)}}}&{{\mathbf{X}^{\left(I_{3}\right)}}}&\cdots&{{\mathbf{X}^{\left(2\right)}}}\\ {{\mathbf{X}^{\left(2\right)}}}&{{\mathbf{X}^{\left(1\right)}}}&\cdots&{{\mathbf{X}^{\left(3\right)}}}\\ \vdots&\vdots&\ddots&\vdots\\ {{\mathbf{X}^{\left(I_{3}\right)}}}&{{\mathbf{X}^{\left(I_{3}-1\right)}}}&\cdots&{{\mathbf{X}^{\left(1\right)}}}\end{array}}\right].$$

$$\mathcal{C}=\mathcal{A}\ast \mathcal{B}=fold\left(bcirc\left(\mathcal{A}\right)\cdot unfold\left(\mathcal{B}\right)\right).$$

$$\mathcal{X}=\mathcal{U}\ast \mathcal{S}\ast {\mathcal{V}}^{\top }.$$

$$ran{k}_t\left(\mathcal{X}\right)={\sum }_{k=1}^{\min \left\{I_1,I_2\right\}}{\mathbf{I}}_{\mathcal{S}\left(k,k,:\right)\ne 0}.$$

$${\Vert \mathcal{X}\Vert }_{\ast }=\frac{1}{I_3}\sum _{k=1}^{I_3}{\Vert {\bar{X}}^{\left(k\right)}\Vert }_{\ast }.$$

$$s_p^{\mu }\left(x\right)=sign\left(x\right)max\left\{\left|x\right|-\mu {\left|x\right|}^{p-1},0\right\}.$$

$${\Vert \mathcal{X}\Vert }_p:=\frac{1}{I_3}\sum _{k=1}^{I_3}{\Vert {\bar{\mathbf{X}}}^{\left(k\right)}\Vert }_p:=\frac{1}{I_3}\sum _{k=1}^{I_3}\sum _{i=1}^{\min \left\{I_1,I_2\right\}}{S}_p^{\mu }\left(\bar{\mathcal{S}}\right).$$

$${\left(S_p^{\mu }\left(\bar{\mathcal{S}}\right)\right)}_{iik}=s_p^{\mu }\left({\bar{\mathcal{S}}}_{iik}\right)=max\left\{\left|{\bar{\mathcal{S}}}_{iik}\right|-\mu {\left|{\bar{\mathcal{S}}}_{iik}\right|}^{p-1},0\right\}.$$

$${\Vert \rho {\mathcal{X}}_1+\left(1-\rho \right){\mathcal{X}}_2\Vert }_p\le \rho {\Vert {\mathcal{X}}_1\Vert }_p+\left(1-\rho \right){\Vert {\mathcal{X}}_2\Vert }_p.$$

$${\Vert \mathcal{X}\Vert }_p={\Vert \mathcal{U}\ast \mathcal{X}\ast {\mathcal{V}}^{\top }\Vert }_p.$$

$$pro{x}_{\tau {\Vert \mathcal{X}\Vert }_p}\left(\mathcal{Z}\right)=\arg \underset{\mathcal{X}}{\min }\frac{1}{2}{\Vert \mathcal{X}-\mathcal{Z}\Vert }_F^2+\tau {\Vert \mathcal{X}\Vert }_p.$$

$$pro{x}_{\tau {\Vert \mathcal{X}\Vert }_p}\left(\mathcal{Z}\right)=\mathcal{U}\ast \mathcal{D}\ast {\mathcal{V}}^{\top }.$$

$$\underset{\mathcal{X}}{\min }\frac{1}{2}{\Vert {\mathcal{P}}_{\Omega }\left(\mathcal{X}-\mathcal{T}\right)\Vert }_F^2+\lambda {\Vert \mathcal{X}\Vert }_p.$$

$$\mathrm{R}\left(\mathcal{X}\right)=\frac{1}{I_1{I}_2{I}_3}{\Vert \mathcal{X}-{\mathcal{X}}^{\ast }\Vert }_F^2.$$

$$\frac{1}{I_1{I}_2{I}_3}{\Vert \mathcal{X}-{\mathcal{X}}^{\ast }\Vert }_F^2\le Cmax\left\{\frac{r{I}_1{I}_3\log \left(3/\alpha \right)}{\left|\Omega \right|},\sqrt{\frac{\log \left(3/\alpha \right)}{\left|\Omega \right|}}\right\}.$$

$$\mathcal{L}\left(\mathcal{X},\mathcal{Y},\mathcal{Z},\beta \right)=\lambda {\Vert \mathcal{Y}\Vert }_p+\frac{1}{2}{\Vert {\mathcal{P}}_{\Omega }\left(\mathcal{X}-\mathcal{T}\right)\Vert }_F^2+⟨\mathcal{Z},\mathcal{X}-\mathcal{Y}⟩+\frac{\beta }{2}{\Vert \mathcal{X}-\mathcal{Y}\Vert }_F^2.$$

$$\{\begin{array}{rl}{\mathcal{Y}}^{t+1}& =\arg \underset{\mathcal{Y}}{\min }\mathcal{L}\left({\mathcal{X}}^t,\mathcal{Y},{\mathcal{Z}}^t,{\beta }^t\right)\\ {\mathcal{X}}^{t+1}& =\arg \underset{\mathcal{X}}{\min }\mathcal{L}\left(\mathcal{X},{\mathcal{Y}}^{t+1},{\mathcal{Z}}^t,{\beta }^t\right)\\ {\mathcal{Z}}^{t+1}& ={\mathcal{Z}}^t+{\beta }^t\left({\mathcal{Y}}^{t+1}-{\mathcal{X}}^{t+1}\right)\\ {\beta }^{t+1}& =\min \left(\eta {\beta }^t,{\beta }^{\max }\right)\end{array}$$

$$\begin{array}{rl}{\mathcal{Y}}^{t+1}& =\arg \underset{\mathcal{Y}}{\min }\lambda {\Vert \mathcal{Y}\Vert }_p+⟨{\mathcal{Z}}^t,{\mathcal{X}}^t-\mathcal{Y}⟩+\frac{\beta }{2}{\Vert {\mathcal{X}}^t-\mathcal{Y}\Vert }_F^2\\ & =\arg \underset{\mathcal{Y}}{\min }\lambda {\Vert \mathcal{Y}\Vert }_p+\frac{\beta }{2}{\Vert {\mathcal{X}}^t-\mathcal{Y}+\frac{1}{\beta }{\mathcal{Z}}^t\Vert }_F^2\end{array}.$$

$${\mathcal{Y}}^{t+1}=pro{x}_{\lambda {\Vert \mathcal{Y}\Vert }_p}\left({\mathcal{X}}^t-\frac{1}{\beta }{\mathcal{Z}}^t\right).$$

$$\begin{array}{rl}{\mathcal{X}}^{t+1}& =\arg \underset{\mathcal{X}}{\min }\frac{1}{2}{\Vert {\mathcal{P}}_{\Omega }\left(\mathcal{X}-\mathcal{O}\right)\Vert }_F^2+⟨{\mathcal{Z}}^t,\mathcal{X}-{\mathcal{Y}}^{t+1}⟩+\frac{\beta }{2}{\Vert \mathcal{X}-{\mathcal{Y}}^{t+1}\Vert }_F^2\\ & =\arg \underset{\mathcal{X}}{\min }\frac{1}{2}{\Vert {\mathcal{P}}_{\Omega }\left(\mathcal{X}-\mathcal{O}\right)\Vert }_F^2+\frac{\beta }{2}{\Vert \mathcal{X}-{\mathcal{Y}}^{t+1}+\frac{1}{\beta }{\mathcal{Z}}^t\Vert }_F^2\end{array}.$$

$${\mathcal{X}}^{t+1}={\mathcal{P}}_{\bar{\Omega }}\left({\mathcal{Y}}^{t+1}+\frac{1}{\beta }{\mathcal{Z}}^t\right)+\mathcal{T}.$$

$${\mathcal{Y}}^{t+1}=pro{x}_{\lambda {\Vert \mathcal{Y}\Vert }_p}\left({\mathcal{W}}^t-\frac{1}{{\beta }^t}{\mathcal{Z}}^t\right);$$

$${\mathcal{X}}^{t+1}={\mathcal{P}}_{\bar{\Omega }}\left({\mathcal{Y}}^{t+1}+\frac{1}{{\beta }^t}{\mathcal{Z}}^t\right)+\mathcal{T};$$

$${\mathcal{Z}}^{t+1}={\mathcal{Z}}^t+{\beta }^t\left({\mathcal{Y}}^{t+1}-{\mathcal{X}}^{t+1}\right);$$

Full text

1 Introduction

In the fields of computer vision and signal processing, there are massive multi-dimensional data that needs to be analyzed and processed. In particular, multi-dimensional arrays (i.e., tensors) provide a natural representation form for these data. Tensor, which is regarded as a multi-linear generalization of matrix/vector, can mathematically model the multi-dimensional data structures, making tensor learning so attractive that there are increasing applications in computer vision [1, 2, 3, 4], machine learning [5, 6], signal processing [7, 8], and pattern recognition [9] in recent years. Unfortunately, the challenge of missing elements in the actually observed tensors limits its applications. Inspired by matrix completion [10], tensor completion attempts to recover the underlying tensor, that on the low-rank assumption, from its incomplete observations. Mathematically, the low-rank tensor completion (LRTC) problem can be formulated as the following model:

[TABLE]

where ${{\cal X}},{{\cal T}}\in{{}^{{I_{1}}\times\cdots\times{I_{n}}}}$ denotes the underlying tensor and its observation tensor, respectively, $rank\left({{\cal X}}\right)$ denotes the rank function of ${\cal X}$, $\Omega$ denotes the index set of observed entries, and ${{{\cal P}}_{\Omega}}\left({{\cal X}}\right)$ denotes the projection operator that ${\left[{{{{\cal P}}_{\Omega}}\left({{\cal X}}\right)}\right]_{{i_{1}}\cdots{i_{n}}}}={{{\cal X}}_{{i_{1}}\cdots{i_{n}}}}$ if $\left({{i_{1}},\cdots,{i_{n}}}\right)\in\Omega$ and 0 otherwise. Different from matrix case, however, the definition of the tensor rank is not unique, and we can see that each definition corresponds to the particular tensor decomposition [11]. For example, tensor n-rank [12] is related to the Tucker decomposition [13], CP-rank [14] is related to the CANDECOMP/PARAFAC (CP) decomposition [15], tensor train (TT) rank is related to the TT decomposition [16], and tubal rank [2] and tensor multi-rank are related to tensor singular value decomposition (t-SVD) [17].

Due to the non-convex and non-smoothness of the rank function $rank\left(\cdot\right)$, the calculation of it is usually NP-hard problem that cannot be solved within polynomial time. Therefore, the function $rank\left({{\cal X}}\right)$ is usually relaxed as its convex/non-convex surrogate, and (1) can be rewritten as follows:

[TABLE]

where $f$ denotes the surrogate function, $\lambda$ denotes the regularization parameter, and $\ell$ denotes the smooth loss function. Obviously, the difference among the existing LRTC models mainly lies in the choice of $f$. Although various tensor rank surrogates [3, 6, 18] are proposed to approximate the tensor rank, they all face some challenges in practical applications.

According to the CP decomposition [15], Friedland et al. [18] claimed that the CP-rank could be relaxed it with CP tensor nuclear norm (CNN), which is a convex surrogate that can be defined as:

[TABLE]

where ${{\cal X}}\in{{}^{{I_{1}}\times\cdots\times{I_{n}}}}$, $\left\|{{\mathbf{u}_{j,i}}}\right\|=1,j=1,\cdots,n$, and $\circ$ denotes the outer product. Moreover, some mathematical properties of CNN are also presented. Under the framework of Frank-Wolfe method, Yang et al. [19] directly relaxed the CP-rank with CNN and pointed out that the proposed model can be solved by converting it to calculate the corresponding spectral norm and rank-one tensors. In addition, Yuan et al. [20] developed the sub-differential of CNN and attempted to recover the underlying tensor with a dual certificate.

Although the CNN minimization defined above is convex, any efficient implementation has not been developed so far. Indeed, it is generally NP-complete to compute the CP-rank of tensor, as well as the tensor spectral norm [21]. Furthermore, the tightness of CNN relative to the CP-rank is hard to measure. Therefore, the application of CNN in LRTC problem is limited.

To avoid the NP-complete CP-rank calculations, tensor n-rank [12] is another commonly used tensor rank. For a given tensor ${{\cal X}}\in{{}^{{I_{1}}\times\cdots\times{I_{n}}}}$, the tensor n-rank of $\cal X$ is defined as ${{rank}_{tc}}\left({{\cal X}}\right)=\left({rank\left({{\mathbf{X}_{\left(1\right)}}}\right),\cdots,rank\left({{\mathbf{X}_{\left(n\right)}}}\right)}\right)$, where ${\mathbf{X}_{\left(i\right)}}$ denotes the mode-$i$ unfolding matrix of ${\cal X}$. Generalized from the nuclear norm in matrix case, Liu et al. [1] first defined the sum of nuclear norm (SNN) as a convex relaxation of the tensor n-rank:

[TABLE]

where ${\left\|\cdot\right\|_{*}}$ denotes the nuclear norm of a matrix. In addition to the basic definition of SNN, an average version was proposed by Gandy et al. [22] for tensor completion. In particular, Signoretto et al. [23] further generalized SNN to the Shatten-$p,q$ norm. Considering the subspace structure in each mode, Kasai et al. [24] developed a LRTC method based on the Riemannian manifold. Moreover, many other norms based on SNN such as latent trace norm and scaled latent trace norm were proposed to be an approximation of the tensor n-rank. More recently, combining the tensor total variation and SNN, Yokota et al. [25] proposed a new LRTC model that can be solved based on primal-dual splitting framework.

However, several limitations of SNN occur with the increasing tensor dimension and scale, including: 1) the operation that simply unfold tensor into matrices along each mode ignores the tensor’s intrinsic structure. 2) The SNN model is neither the tightest convex lower bound of tensor n-rank nor the optimal solution with the dimension increasing. 3) The unfolding and folding are expensive.

Unlike the unfolding method of tensor n-rank, the tensor train (TT) rank [26] employs a well-balanced matricization scheme to capture the global correlation of tensor entries. For a given tensor ${{\cal X}}\in{{}^{{I_{1}}\times\cdots\times{I_{n}}}}$, Imaizumi et al. [27] and Bengua et al. [26] matricize the tensor along permutations of modes and define the tensor train nuclear norm (TTNN) based on TT rank as:

[TABLE]

where ${X_{\left[i\right]}}\in{{}^{{m_{1}}\times{m_{2}}}}$, ${m_{1}}=\prod\nolimits_{k=1}^{i}{{I_{k}}}$ and ${m_{2}}=\prod\nolimits_{k=i+1}^{n}{{I_{k}}}$ denotes the mode-$\left({1,\cdots,i}\right)$ unfolding matrix of ${\cal X}$. Imaizumi et al. [27] established the statistical theory and developed a scalable algorithm for the TTNN model. They provided a statistical error bound for TTNN which achieves the same efficiency as SNN. Moreover, Bengua et al. [26] stated that tensor train nuclear norm is more tractable than SNN for LRTC, and they introduced ket augmentation scheme to obtain a higher order tensor from a given tensor. Furthermore, algorithmic development for the solution of TTNN model was also devoted in the research.

Although the matricization scheme ${\mathbf{X}_{\left[i\right]}}$ is more balance than the unfolding ${\mathbf{X}_{\left(i\right)}}$ , the matricize operation may also destroy the original tensor’s internal structure. Moreover, the TTNN model is efficient for higher order tensors.

To avoid the tensor matricization and maintain the intrinsic data structure, Kilmer et al. [17] proposed tensor singular value decomposition (t-SVD) and this motivates the new tensor multi-rank and tubal rank. Another advantage of such method is that the resultant algebra and analysis are very close to those of matrix case. Zhang et al. [28] give the definition of a new tensor nuclear norm (TNN) corresponding to tubal rank and t-SVD. Leveraging the conclusion of matrix case, they state that TNN is the tightest convex approximation of tensor average rank. Furthermore, they derived the exact recovery conditions for LRTC problems in [29].

Due to the advantages of TNN, it has received more and more attention in recent years. Nevertheless, the computational complexity of TNN increases dramatically for the large scale tensor. Additionally, TNN penalizes all the singular values with the same weight. In fact, the larger the singular value, the more information it contains, and the less penalized it should be. To conquer the challenge of high complexity, Zhou et al. [30] focus on utilizing the technique of tensor factorization for LRTC problems. Only two smaller tensors are maintained in the optimization process, which can be used to preserve the low-rank structure of the underlying tensor. In addition, extending Schatten-$p$ norm to tensor space, Kong et al. [31] state that their propose tensor Schatten-$p$ norm can better approximate the tensor average rank.

Recently, Voronin et al. [32] generalize the iterative soft thresholding method to $p$-shrinkage thresholding for solving sparse signal recovery (SSR) problems. Both theoretical and empirical results in SSR problem prove that the $p$-shrinkage thresholding function is a good alternative to the Schatten-$p$ norm, which can achieve better recovery performance than the existing surrogates. Motivated by this success application, we attempt to introduce the $p$-shrinkage scheme to LRTC problem and develop an efficient algorithm to solve the resulting model.

In this paper, we propose a new tensor $p$-shrinkage nuclear norm ($p$-TNN), which is defined in (13) based on the t-SVD and $p$-shrinkage scheme. When $-\infty<p<1$, the proposed norm can achieve a better approximation of the tensor average rank than TNN, i.e., $p$-TNN is a tighter surrogate. Extending the matrix norm surrogate to the tensor case, we establish the $p$-TNN based LRTC model. Additionally, the recovery guarantee of the proposed model is provided. To cope with the challenges of the resultant non-convex optimization problem, we develop an efficient algorithm under the alternating direction method of multipliers (ADMM) framework. Furthermore, we also incorporate the adaptive momentum scheme to accelerate the empirical convergence for the proposed algorithm. Subsequently, the resulting algorithm is analyzed in detail from the aspects of time complexity and convergence, respectively.

In summary, the primary contributions of our work include:

• •

We propose a new definition of tensor $p$-shrinkage nuclear norm with some desirable properties, for example, positivity and unitary invariance. The proposed $p$-TNN ($-\infty<p<1$) is a tighter envelope of the tensor average rank than TNN within the unit ball of the spectral norm, which is beneficial to improve the recovery performance of LRTC problem.

• •

By employing the tensor $p$-shrinkage nuclear norm, we propose a novel LRTC model and provide a strong guarantee for tensor recovery, i.e., the error in recovering a ${I_{1}}\times{I_{2}}\times{I_{3}}$ tensor is ${\cal O}\left({{{r{I_{1}}{I_{3}}\log\left({{3\mathord{\left/{\vphantom{3\alpha}}\right.\kern-1.2pt}\alpha}}\right)}}/{{\left|\Omega\right|}}}\right)$, where $\alpha={\textstyle{3\over{\left({{I_{1}}+{I_{2}}}\right){I_{3}}}}}$, $r$ denotes the tensor tubal rank, and $\left|\Omega\right|$ denotes the cardinal number of the index set.

• •

Incorporating the adaptive momentum scheme, we develop an efficient algorithm which establishes a unified framework for algorithms that applying the soft and hard thresholding shrinkage. Under the smoothness assumption, the Convergence guarantee to critical points is provided.

2 Notations and preliminaries

In this section, we first introduce some main basic notations and then briefly give some necessary definitions which will be used later.

Tensors are denoted by uppercase calligraphy letters, e.g., ${\cal{A}}$. Matrices are denoted by uppercase boldface letters, e.g., $\mathbf{A}$. Vectors are denoted by lowercase boldface letters, e.g., $\mathbf{a}$, while scalars are denoted by lowercase letters, e.g., $a$.

The fields of complex numbers and real numbers are denoted by ${\mathbb{C}}$ and ${\mathbb{R}}$, respectively. We simply represent ${\left\{1,2,\cdots n\right\}}$ by ${\left[n\right]}$. For a n-dimensional tensor ${{\cal{A}}\in\mathbb{C}{{}^{{I_{1}}\times{I_{2}}\times\cdots\times{I_{n}}}}}$, we use ${\cal A}_{i_{1}i_{2}\cdots i_{n}}$ to represent its $\left({i_{1}i_{2}\cdots i_{n}}\right)$-th entry, where ${i_{k}\in\left[I_{k}\right]}$ and ${k\in{\left[n\right]}}$. Let ${\mathbf{A}^{\left(k\right)}}$ denote the $k$-th frontal slice of the tensor ${\cal{A}}$. We use ${\cal A}^{\top}$ to denote its transpose tensor. The inner product of ${\cal{A}}$ and ${\cal{B}}$ is defined as $\left\langle{{\cal A},{\cal B}}\right\rangle=\sum\nolimits_{{i_{1}},{i_{2}},\cdots,{i_{n}}}{{{\cal A}_{{i_{1}},{i_{2}},\cdots,{i_{n}}}}}{{\cal B}_{{i_{1}},{i_{2}},\cdots,{i_{n}}}}$, and the Frobenius norm of ${\cal A}$ is defined as ${{\left\|{\cal A}\right\|_{F}}=\sqrt{\left\langle{{\cal A},{\cal A}}\right\rangle}=\sqrt{\sum\nolimits_{{i_{1}},{i_{2}},\cdots,{i_{n}}}{\left(A_{{i_{1}},{i_{2}},\cdots,{i_{n}}}\right)^{2}}}}$.

For a given 3-dimensional tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, we use ${\bar{\cal X}}\in{\mathbb{C}^{I_{1}\times I_{2}\times I_{3}}}$ to denote the Discrete Fourier Transformation (DFT) of ${\cal X}$ along the third dimension, i.e., ${\bar{\cal X}}=fft{\left({\cal X},[{\kern 1.0pt}{\kern 1.0pt}],3\right)}$. Correspondingly, the $k$-th frontal slice of ${\bar{\cal X}}$ is denoted by ${{\bar{\mathbf{X}}}^{\left(k\right)}}$. The unfold and its inverse operator of $\cal X$ [2] are defined as

[TABLE]

And the block circulant matrix of ${\cal X}$ is further defined as ${bcirc\left({\cal X}\right)}\in\mathbb{R}^{I_{1}I_{3}\times I_{2}I_{3}}$ [2]:

[TABLE]

Then we summarize some necessary definitions and results.

Definition 1. (t-product [17]) For two given tensors ${\cal A}\in{\mathbb{R}^{I_{1}\times d\times I_{3}}}$ and ${\cal B}\in{\mathbb{R}^{d\times I_{2}\times I_{3}}}$, the t-product of $\cal A$ and $\cal B$ is defined as:

[TABLE]

where ${\cal C}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$. Therefore, the t-product is homologous in form to the matrix multiplication except that the operation of circulant convolution and unfolding. Additionally, we can also note that if $I_{3}=1$, the t-product reduces to matrix multiplication.

In addition, the definitions of orthogonal tensor, identity tensor, frontal-slice-diagonal tensor (f-diagonal tensor)and tensor transpose can be found in the Appendix. Subsequently, the Tensor Singular Value Decomposition (t-SVD) can be defined as follows by using these above definitions.

Definition 2. (t-SVD [17]) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, there exist ${\cal U}\in{\mathbb{C}^{I_{1}\times I_{1}\times I_{3}}}$, ${\cal V}\in{\mathbb{C}^{I_{2}\times I_{2}\times I_{3}}}$ and ${\cal S}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ such that:

[TABLE]

where ${\cal U}$ and ${\cal V}$ are orthogonal tensor, i.e., ${\cal U}*{\cal U}^{\top}={\cal V}*{\cal V}^{\top}={\cal I}$, and ${\cal S}$ is a f-diagonal tensor.

Definition 3. (Tensor tubal rank [28] and multi-rank [29]) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, then the tensor tubal rank of ${\cal X}$, herein denoted by $rank_{t}\left({\cal X}\right)$, is defined to be the number of nonzero singular tubes of ${\cal S}$, i.e.,

[TABLE]

where ${\cal S}$ is defined in Definition 2, ${\cal S}\left({k,k,:}\right)$ denotes the $k$-th diagonal tube of ${\cal S}$, and $\mathbf{I}_{{\cal S}\left({k,k,:}\right)\neq 0}$ is an indicator function, i.e., $\mathbf{I}_{{\cal S}\left({k,k,:}\right)\neq\mathbf{0}}=1$ if ${{\cal S}\left({k,k,:}\right)\neq\mathbf{0}}$ is true and $\mathbf{I}_{{\cal S}\left({k,k,:}\right)\neq\mathbf{0}}=0$ otherwise. The tensor multi-rank of ${\cal X}$ is a vector $\mathbf{r}\in{\mathbb{R}^{I_{3}}}$ consisting of the rank of the frontal slice of $\bar{\cal X}$, i.e., ${\mathbf{r}}_{k}=rank\left({{{\bar{\mathbf{X}}}^{\left(k\right)}}}\right)$.

Definition 4. (Tensor nuclear norm and spectral norm [33]) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, the tensor nuclear norm of ${\cal X}$, denoted by ${\left\|{\cal X}\right\|_{*}}$, is defined as the average of the nuclear norm of all the frontal slices of ${\bar{\cal X}}$:

[TABLE]

Furthermore, the tensor spectral norm of $\cal X$, denoted by ${\left\|{{\bar{\cal X}}}\right\|}$, is defined as ${\left\|{{\bar{\cal X}}}\right\|}=\mathop{\max}\limits_{i}\left\|{{{\bar{\mathbf{X}}}^{\left(k\right)}}}\right\|$.

According to the Von Neumann’s inequality, it can be seen that the tensor spectral norm is the dual norm of the tensor nuclear norm and vice versa.

Definition 5. ($p$-shrinkage thresholding operator [32]) For a given scalar $x\in{\mathbb{R}}$, $\forall\mu>0$ and $p\leq 1$, the $p$-shrinkage thresholding operator, denoted by $s_{p}^{\mu}$, is defined as:

[TABLE]

where $sign\left(x\right)$ denotes the sign function.

Several $p$-shrinkage functions with different $p$ values are shown in Fig.1. Note that if $p=1$, the $p$-shrinkage thresholding transforms to the soft thresholding, which imposes no different penalty for all the inputs and is commonly used in matrix case, and if $p\to-\infty$, it turns into the hard thresholding with no penalty for large inputs. However, hard thresholding is discontinuous. Moreover, when $-\infty<p<1$, the $p$-shrinkage function satisfies the property that the larger the input, the less penalized it should be.

3 Model

In this section, we propose a new definition of tensor $p$-shrinkage nuclear norm ($p$-TNN) based on $p$-shrinkage scheme and t-SVD. Subsequently, we propose a novel LRTC model by employing our $p$-TNN.

Definition 5. (Tensor $p$-shrinkage nuclear norm, $p$-TNN) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, let ${\cal X}={\cal U}*{\cal S}*{\cal V}^{\top}$ be the tensor singular value decomposition of ${\cal X}$. Then the tensor $p$-shrinkage nuclear norm is defined as:

[TABLE]

where $-\infty<p<1$, and ${S_{p}^{\mu}\left({{\bar{\cal S}}}\right)}:{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}\to{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ is defined as follows:

[TABLE]

where ${{\bar{\cal S}_{kki}}}$ denotes the $\left(k,k,i\right)$-th tensor singular value of ${\bar{\cal X}}$.

The following are some of the properties of our proposed $p$-TNN that we use in this paper. For the proofs, please refer to the Appendix.

Proposition 1. (Positivity) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, the $p$-TNN of ${\cal X}$, denoted by ${\left\|{\cal{X}}\right\|_{p}}$. Obviously, ${\left\|{\cal{X}}\right\|_{p}}\geq 0$ with equality holding if and only if ${\cal X}$ is zero. Moreover, ${\left\|{\cal{X}}\right\|_{p}}$ is non-decreasing.

Proposition 2. (Non-convexness) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, if $-\infty<p<1$, ${\left\|{\cal{X}}\right\|_{p}}$ is non-convex w.r.t. ${\cal X}$. i.e., ${\left\|{\cal{X}}\right\|_{p}}$ cannot satisfy the inequality for any $\rho\in\left(0,1\right)$:

[TABLE]

where ${\cal X}_{1},{\cal X}_{2}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ and ${\cal X}_{1}\neq{\cal X}_{2}$.

Proposition 3. (Unitary invariance) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, there exist orthogonal tensors ${\cal U}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ and ${\cal V}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, such that:

[TABLE]

In addition, ${\left\|{{\cal X}}\right\|_{p}}={\left\|{{{\cal U}}*{{\cal X}}}\right\|_{p}}={\left\|{{{\cal X}}*{{{\cal V}}^{\top}}}\right\|_{p}}$ also holds. Therefore, the tensor $p$-shrinkage unclear norm is unitary invariance.

Extending the proximal operator to the tensor case, then for any $\tau>0$, we can define the proximal operator of the tensor $p$-shrinkage nuclear norm ${\left\|{\cdot}\right\|_{p}}$ as:

[TABLE]

The following definition extends the generalized singular value thresholding (GSVT) to the tensor space. Similar to the matrix case, the following shows that the proximal operator of the tensor $p$-shrinkage nuclear norm has a closed-form solution.

Definition 6. (Tensor generalized singular value thresholding, t-GSVT) Let ${\cal X},{\cal Z}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, for any $\tau>0$, the solution of the proximal operator of $p$-TNN can be calculated as:

[TABLE]

where ${\cal U}*{\cal S}*{\cal V}^{\top}$ is the t-SVD of ${\cal Z}$, and the f-diagonal tensor ${\cal D}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ that satisfies ${\bar{\cal D}}=S_{p}^{\mu}\left({\bar{\cal S}}\right)$.

For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, Lu et al. [33] defined the tensor average rank as $rank_{a}\left({\cal X}\right)=\frac{1}{{I_{3}}}\sum\nolimits_{k=1}^{I_{3}}{rank\left({{{\bar{\mathbf{X}}}^{\left(i\right)}}}\right)}$, and they pointed out that TNN is the convex envelope of the tensor average rank within the unit ball of the tensor spectral norm. Recently, Kong et al. [31] proposed a tensor Schatten-$p$ norm and claimed that it can be a tighter non-convex approximation of the tensor average rank than TNN. Here, we prove that the proposed $p$-TNN ${\left\|{\cal X}\right\|_{p}}$ is even a tighter envelope within the same unit ball.

Proposition 4. (Tightness) For a given tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, when $-\infty<p<1$, ${\left\|{\cal X}\right\|_{p}}$ is an non-convex envelope of the tensor average rank within the unit ball of the spectral norm, which is tighter than TNN in the sense of $\left\|{{\cal X}}\right\|_{*}\leq{\left\|{{\cal X}}\right\|_{p}}\leq rank_{a}\left({\cal X}\right)$.

Considering the LRTC model defined in (2), the non-smoothness rank function is usually relaxed as its convex/non-convex surrogate. By using Definition 4 and its Proposition 4, we can rewrite the tensor completion problem into the following form:

[TABLE]

where ${\cal T}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ denotes the observed tensor, $\Omega$ denotes the index set of the observed entries, and ${\cal P}_{\Omega}$ denotes the projection operator, i.e., ${\left[{{{\cal P}_{\Omega}}\left({\cal A}\right)}\right]_{i_{1}i_{2}i_{3}}}={{\cal A}_{i_{1}i_{2}i_{3}}}$ if $\left(i_{1},i_{2},i_{3}\right)\in\Omega$ and 0 otherwise. We aim to recover the completed tensor (underlying tensor) ${\cal X}$ based on the observed missing tensor ${\cal T}$.

Recovery guarantee. Following the definition of the recovery error in [34], here the the recovery error can be formulated as:

[TABLE]

where ${\cal X},{\cal X}^{*}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ are the underlying tensor and its recovery tensor, respectively. The following theorem establishes the upper bound of the recover error based on our proposed $p$-TNN.

Theorem 1. (Main Result 1) Let ${\cal X}^{*}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ be the solution to (19), then with a probability of at least$1-\alpha$, we have

[TABLE]

where $C$ is an absolute constant, $r$ denotes the tensor tubal rank of the underlying tensor, ${\left|\Omega\right|}$ denotes the cardinal number of the index set $\Omega$, and $\alpha={\textstyle{3\over{\left({I_{1}+I_{2}}\right)I_{3}}}}$. The proof can be found in Appendix.

According to (21), it can be seen that the upper bound of the recovery error is related to the tubal of rank the underlying tensor, as well as the sampling rate $sr$, which can be calculated by $sr={\textstyle{{\left|\Omega\right|}\over{I_{1}I_{2}I_{3}}}}$.

When using a tighter surrogate, one can get a better solution. From the above narrative, we can see that the main advantage of $p$-TNN is the outstanding tightness (Proposition 4), i.e., our $p$-TNN is a better approximation of the tensor average rank, which can lead to a better solution. However, $\left\|{\cal X}\right\|_{p}$ is non-convex when $-\infty<p<1$. The resulting non-convex optimization problem is much more challenging. A strong performance guarantee would be hard to get as in the convex case.

4 Algorithm

In this section, we will show that the proposed non-convex model (19) can be solved effectively based on the ADMM framework. Moreover, we will show that the proposed algorithm can be further accelerated with adaptive momentum.

Let $F\left({{\cal X}}\right)\equiv\frac{1}{2}\left\|{{{{\cal P}}_{\Omega}}\left({{{\cal X}}-{{\cal T}}}\right)}\right\|_{F}^{2}+\lambda{\left\|{{\cal X}}\right\|_{p}}$. By introducing an auxiliary variable ${\cal Y}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, such that ${\cal Y}={\cal X}$, the augmented Lagrangian function of (19) is given as follows:

[TABLE]

where ${\cal Z}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$ denotes the Lagrangian multiplier, and $\beta$ denotes the penalty parameter. In the following, (22) can be solved by applying the classical ADMM framework. Therefore, ${\cal X}$, ${\cal Y}$, ${\cal Z}$ and $\beta$ are iteratively update as:

[TABLE]

1) For ${\cal Y}$-subproblem.

By fixing ${\cal X}^{t}$, ${\cal Z}^{t}$ and $\beta^{t}$, we can obtain ${\cal Y}^{t+1}$ by:

[TABLE]

It can be seen that (24) satisfies the form of the proximal operator of $p$-TNN. Therefore, according to (17), ${\cal Y}^{t+1}$ can be rewritten as:

[TABLE]

2) For ${\cal X}$-subproblem.

By fixing ${\cal Y}^{t+1}$, ${\cal Z}^{t}$ and $\beta^{t}$, we can obtain ${\cal X}^{t+1}$ by:

[TABLE]

The above equation has a closed-form solution as following:

[TABLE]

where ${\bar{\Omega}}$ denotes the complementary set of the index set $\Omega$.

Next, we will show the above update step can be accelerated by the adaptive momentum, which has been popularly used for stochastic gradient descent and proximal algorithms. The idea of the adaptive momentum is to apply historical iterations to accelerate convergence. The whole procedure is shown in Algorithm 1. Let the underlying tensor ${\cal X}\in{\mathbb{R}^{I_{1}\times I_{2}\times I_{3}}}$, we suppose that the tubal rank of $\cal X$ is $rank_{t}\left({\cal X}\right)=r$. Obviously, the main per-iteration cost lies in the update of ${\cal Y}^{t+1}$ (step 9 in Algorithm 1), which takes ${\cal O}\left(r\left(I_{1}+I_{2}\right)I_{3}logI_{3}+r\left|\Omega\right|\right)$ time. As for TNN in [33], the computational complexity at each iteration is ${\cal O}\left(I_{1}I_{2}I_{3}\left(logI_{3}+min\left\{I_{1},I_{2}\right\}\right)\right)$. With $r<min\left\{I_{1},I_{2}\right\}$ and $\left|\Omega\right|<I_{1}I_{2}I_{3}$, it can be seen that our algorithm is much more efficient than TNN based methods in each iteration.

Convergence analysis. In this section, we investigate the convergence of the Algorithm 1. $\left({\cal X}^{*},{\cal Y}^{*},{\cal Z}^{*}\right)$ is the KKT point of the problem (14), if it satisfies the following system:

[TABLE]

where ${\cal U}^{*}*{\cal S}^{*}*{{\cal V}^{*}}^{\top}$ denotes the t-SVD of $\left({\cal X}^{*}-\frac{1}{\beta}{{{\cal Z}}^{*}}\right)$, and $\bar{\cal D}^{*}=S_{p}^{\mu}\left({\bar{\cal S}^{*}}\right)$. Let $F\left({{\cal X}}\right)\equiv\frac{1}{2}\left\|{{{{\cal P}}_{\Omega}}\left({{{\cal X}}-{{\cal T}}}\right)}\right\|_{F}^{2}+\lambda{\left\|{{\cal X}}\right\|_{p}}$, the following lemma shows that $F\left({\cal X}\right)$ is always non-increasing ($\beta>0$) as the iterations proceed.

Lemma 1. According to the properties of the $p$-TNN, and when $-\infty<p<1$, let $\left\{{\cal X}^{t}\right\}$ be the sequence generated by Algorithm 1. Then, we have the following inequality:

[TABLE]

Theorem 2. Let $\left\{{\cal X}^{t}\right\}$ be the sequence produced by algorithm 1, we say $\left\{{\cal X}^{t}\right\}$ is a bounded iterative sequence, i.e., $\sum\nolimits_{t=1}^{\infty}{\left\|{{{\cal X}^{t+1}}-{{\cal X}^{t}}}\right\|}_{F}^{2}<\infty$.

According to Theorem 2, it can be seen that the proposed algorithm generates a bounded iterative sequence. Moreover, from the above theorem, we must have $\mathop{\lim}\limits_{t\to\infty}\left\|{{{\cal{X}}^{t+1}}-{{\cal{X}}^{t}}}\right\|_{F}^{2}=0$. Hence, the iteration sequence has limit point.

In (19), we choose $\frac{1}{2}\left\|{{{\cal P}_{\Omega}}\left({{\cal X}-{\cal T}}\right)}\right\|_{F}^{2}$, which is commonly used in LRTC problem and low-rank matrix completion problem, as the loss function. Particularly, the loss function we choose is $\rho$-Lipschitz smooth, i.e., it satisfies $\left\|{\nabla\ell\left({{{\cal{X}}_{1}}}\right)-\nabla\ell\left({{{\cal{X}}_{2}}}\right)}\right\|_{F}^{2}\leq\rho\left\|{{{\cal{X}}_{1}}-{{\cal{X}}_{2}}}\right\|_{F}^{2}$. Similar with [37], we use ${\left\|{{{\cal X}^{t+1}}-{{\cal X}^{t}}}\right\|}_{F}^{2}$ to perform the convergence analysis of the proposed algorithm. The convergence of Algorithm 1 is shown in the following theorem, and the proof can be founded in Appendix.

Theorem 3. Let $\left\{{\cal X}^{t}\right\}$ be the sequence produced by Algorithm 1. For the consecutive elements ${\cal X}^{t}$ and ${\cal X}^{t+1}$, we have

[TABLE]

where ${\rm inf}{\kern 1.0pt}F$ denotes the minimizer of the objective function $F$. The above theorem shows that Algorithm 1 converge to a critical point at the rate of ${\cal O}\left(1/T\right)$.

5 Experiments

In this section, we perform numerical experiments on synthetic and real-world data sets to demonstrate the effectiveness of our proposed algorithm. Furthermore, the experimental results show the superiority of our method. Each experiment is repeated ten times, and the average results are reported. All experiments are implemented in Matlab on Windows 10 with Intel Xeon 2.8GHz CPU and 128GB memory.

5.1 Evaluation metrics

Let ${\cal X},{{\cal{X}}^{*}}\in{{\mathbb{R}}^{{I_{1}}\times{I_{2}}\times{I_{3}}}}$ denote the underlying tensor and its recovery tensor (i.e., the output of the LRTC methods), respectively. The following metrics are chosen to evaluate the recovery performance of the LRTC algorithms.

1. Relative Square Error, denoted by RSE, is defined as:

[TABLE]

2. Peak Signal-to-Noise Ratio, denoted by PSNR, is defined as:

[TABLE]

3. Additionally, for the experiments on real-world data sets, another metric, SSIM, is defined as:

[TABLE]

where ${\mathbf{X}}$ and ${\mathbf{X}}^{*}$ denote the greyscale images for the original image and its recovery image, $c_{1}$ and $c_{2}$ are constants, $\mu_{\mathbf{X}}$ and ${\mu_{{{\mathbf{X}}^{*}}}}$ denote the average values, while $\sigma_{\mathbf{X}}$ and ${\sigma_{{{\mathbf{X}}^{*}}}}$ denote the standard deviation of ${\mathbf{X}}$ and ${\mathbf{X}}^{*}$, respectively, and ${{\sigma_{{\mathbf{X}}{{\mathbf{X}}^{*}}}}}$ denote the covariance matrix between ${\mathbf{X}}$ and ${\mathbf{X}}^{*}$.

5.2 Parameter settings

The stopping criterion of the proposed algorithm is:

[TABLE]

where the maximum number of iterations $T$ is fixed at 1000, and the tolerance $tol$ is fixed at $10^{-4}$. Moreover, with regard to the penalty parameter, it is initialized as 0.01, i.e., ${\beta^{0}}=0.01$, the step size is set to 1.1, and ${\beta^{\max}}={10^{5}}$. Additionally, following the settings in [37], $\gamma_{1}=0.1$, and $\rho=2$.

Furthermore, we investigate the effect of the $p$-value selection on the performance of the proposed method. For a given $100\times 100\times 20$ tensor with the tensor tubal rank 5 (which can be generated with the method in synthetic data), we attempt to recover it from $20\%$ sampling observations (i.e., the sampling rate sr = 0.2) to show the performance of our method with different $p$-values. In the experiments, we follow the above parameter settings and vary $p$ in the range $[-2,0.9]$, and the step size is fixed at 0.1. Fig.2 shows the performance of the proposed method with different $p$-values.

As can be seen from Fig.2, the performance tends to be stable when $p\leq-1$. Meanwhile, the optimal performance is achieved. Therefore, we fix $p=-1$ in the rest of experiments.

5.3 Synthetic data

Data sets: As in [6, 34], we assume that the tensor ${\cal X}\in{{\mathbb{R}}^{{I_{1}}\times{I_{2}}\times{I_{3}}}}$ with tubal rank $r$ can be generated by a tensor product ${\cal X}={\cal P}*{\cal Q}$, where ${\cal P}\in{{\mathbb{R}}^{{I_{1}}\times{r}\times{I_{3}}}}$ and ${\cal Q}\in{{\mathbb{R}}^{{r}\times{I_{2}}\times{I_{3}}}}$ are tensors with elements obtained independently from the ${\cal N}\left(0,1\right)$ distribution. Three kinds of data sets are generated:

(i) We fix the tubal rank at 5 and consider the tensor ${\cal X}\in{{\mathbb{R}}^{{I_{1}}\times{I_{2}}\times{I_{3}}}}$ with $I_{1}=I_{2}=I_{3}=n$ that $n$ varies in $\left\{50,100,150,200\right\}$.

(ii) We fix $I=100$ and consider the tensor ${\cal X}\in{{\mathbb{R}}^{{I}\times{I}\times{I_{3}}}}$ with $I_{3}$ varies in $\left\{20,40,60,80,100\right\}$. Moreover, the tensor tubal rank is also fixed at 5.

(iii) We consider the tensor ${\cal X}\in{{\mathbb{R}}^{{I}\times{I}\times{I_{3}}}}$ with $I_{3}$ ($I=100,I_{3}=20$) with different tubal rank $r$ varying in the range $[10,40]$.

Considering the first data set, we investigate the effect of different tensor size on the performance of the proposed method. The sampling rate is varied from 0.05 to 0.5, with a step size equaling 0.05. Results are shown in Fig.3(a). Note that there exists a negative correlation between the tensor size and RSE. In particular, the underlying tensor can be recovered accurately with a much small sampling rate, when its size is large enough. In addition, for a fixed-size tensor, the RSE decreases monotonously with an increasing sampling rate. It is reasonable – the larger the number of observations is, the more information of the underlying tensor is obtained.

Furthermore, we study the influence of different $I_{3}$ on the recovery performance using the second data set. The sampling rate is varied from 0.05 to 0.5, with a step size equaling 0.05. Results are shown in Fig.3(b). At the same sampling rate, it can be seen that a larger $I_{3}$ leads to a smaller RSE. This conclusion is similar to the discussion of small $I_{3}$ and large $I_{3}$ in [38]. Moreover, with an increasing sampling rate, RSEs decreases monotonously.

Considering the third data set, we investigate the effect of tubal rank on the recovery performance of the proposed method. The sampling rate is varied from 0.05 to 0.5, with a step size equaling 0.015. Furthermore, the PSNRs are normalized. Results are shown in Fig.4, where white means the underlying tensor is recovered successfully while black indicates that the recovery failed. As can be seen, the recovery performance can be guaranteed as long as the tubal rank is relatively low and the sampling rate is relatively large.

5.4 Color image data

Data sets: Zhou et al. [30] point out that most natural images have low tubal rank structure, and following the experiment settings in [39] and [40], we use four color images 111https://github.com/Spring-Liu/testimage. to demonstrate the performance of the proposed method in this part. The four images (named lena, facade, baboon, and house, respectively) are resized to $256\times 256\times 3$ tensor. A summary of these four images is reported in Fig.LABEL:Fig.5.

Baseline: We compare the performance of the proposed algorithm with other state-of-the-art methods, including: i) smooth PARAFAC tensor completion (SPC) [40] which is CP-rank based method, ii) the tensor n-rank based methods, i.e., simultaneous tensor decomposition and completion (STDC) [25], iii) the TT rank based methods, i.e., Tmac-TT [26] and tensor train stochastic gradient descent (TT-SGD) [34], and iv) TNN [33] which is the tubal rank based method.

For each color image, we test the above mentioned LRTC methods with sampling rate equaling 0.1, 0.2, 0.3 and 0.4. Results are shown in Table 1. As can be seen from Table 1, the performance improvement is universal for an increasing sampling rate. In particular, our proposed method is generally superior to other state-of-the-arts. More precisely, in addition to the results of facade at sampling rate 0.1, our method achieves highest PSNR, SSIM and smallest RSE. For further visually compare the recovery performance of all methods, we show the recovered results for each color image in Fig.5 (the sampling rate is fixed at 0.2). Note that the recovered color images of our method are closer to the original images. The reason of our method can obtain better details of images than TNN is that the proposed $p$-TNN is a much tighter approximation of the tubal rank than tensor nuclear norm.

5.5 Hyperspectral image data

Data sets: In this part, we use two public hyperspectral image data sets named Washington DC Mall (WDC Mall) and Pavia University (paviaU)222http://www.ehu.eus/ccwintco/index.php?title=Hyperspectral_Remote_Sensing_Scenes.. The whole WDC Mall data set contains $1208\times 307$ pixels and 191 spectral bands, which is from the hyperspectral digital imagery collection experiment. In the following experiments, it is resized to a $200\times 200\times 191$ tensor. The whole paviaU data set is a $610\times 610$ pixels and 103 spectral bands image, which is from Pavia University. In the rest experiments, it is resized to a $200\times 200\times 103$ tensor. A visually summary is shown in Fig. 6.

Baseline: We compare the performance of the proposed algorithm with the same methods mentioned in color image experiments, i.e., SPC, STDC, Tmac-TT, TT-SGD, and TNN.

For the two hyperspectral images, we test the above mentioned LRTC methods with sampling rate equaling 0.1, 0.2, 0.3 and 0.4. Results are shown in Table 2. We use only PSNR and RSE as the evaluation metrics here. Similarly, we can see that our method almost outperforms others. The improvement of average recovery performance is statistically significant. To further prove the superiority of our proposed method for hyperspectral image inpainting, we show the recovered results of all methods from 0.2 sampling observed images in Fig.6, where only three slices of two data sets (the first three slices of WDC Mall and slices [60,61,62] of paviaU) are shown exemplarily. As can be seen from Fig.6, compared with other methods, our proposed algorithm visually performs better recovery results.

5.6 Conclusion

In this paper, we address the problem of low-rank tensor completion. We propose a new definition of tensor $p$-shrinkage nuclear norm ($p$-TNN). The tightness property of $p$-TNN demonstrates that the proposed $p$-TNN is a tighter surrogate of tensor average rank than tensor nuclear norm (TNN). Therefore, we propose a novel LRTC model by employing our proposed $p$-TNN. In particular, the upper bound of recovery error for our LRTC model is further provided to show the underlying tensor can be recovered accurately. Accordingly, we develop an efficient algorithm by incorporating the adaptive momentum scheme. Subsequently, some theoretical analysis are provided from the aspects of complexity and convergence, respectively. The experimental results validate the superiority of our method over the state-of-the-arts.

However, there still several directions in future work. First, the experimental results show that further improvements are needed for the performance of our method at low sampling rate. Moreover, we will focus on the distributed version to apply it to massive data sets.