On the tensor nuclear norm and the total variation regularization for image and video completion
aa r X i v : . [ m a t h . NA ] F e b ON THE TENSOR NUCLEAR NORM AND THE TOTAL VARIATIONREGULARIZATION FOR IMAGE AND VIDEO COMPLETION
A.H BENTBIB ∗ , A. EL HACHIMI § , K. JBILOU §† , AND
A. RATNANI ‡ Abstract.
In the present paper we propose two new algorithms of tensor completion for three-order tensors. The proposed methods consist in minimizing the average rank of the underlyingtensor using its approximate function namely the tensor nuclear norm and then the recovered datawill be obtained by using the total variation regularisation technique. We will adopt the AlternatingDirection Method of Multipliers (ADM), using the tensor T-product, to solve the main optimizationproblems associated to the two algorithms. In the last section, we present some numerical experimentsand comparisons with the most known image completion methods.
Key words.
ADM, Tensor completion, Tensor nuclear norm, T-product, T-SVD.
1. Introduction.
Tensors become an important notion that treat the high di-mensional data and it plays an important role in a wide range of real-world applica-tions. In this work, we will be interested in the problem of tensor completion, withthe aim of estimating the missing values from an observed data, e.g., inpainting colorimages [3, 10, 19, 27], hyperspectral image recovery [11, 17, 14, 28], magnetic reso-nance image recovery [32, 16], and higher order web site link analysis [2].The main idea behind the completion problem is to find a low-rank data containingthe main information of the original data. For matrices, the problem of completion ismathematically formulated as follows min A rank ( A ) s.t P Ω ( A ) = P Ω ( M ) , (1.1)where A ∈ R n × n is the underlying matrix, M ∈ R n × n is the observed matrix, Ω isthe set of known pixels and P Ω is the projection operator onto Ω. This optimizationproblem is not easy to solve because of the non-convexity of the rank function . Forthat reason, Fazel [6], and Kurucz [13] proposed to solve the problem (1.1) by usingthe convex surrogate of the rank and the SVD decomposition. In [9], Candes andRecht proved theoretically that under certain conditions, the following optimizationproblem min A k A k ∗ s.t P Ω ( A ) = P Ω ( M ) (1.2)recovers well the data, where k A k ∗ is the nuclear norm of A which will be defined later.Since tensors are the generalization of matrices, the problem of tensor completion canbe formulated as follows min A rank ( A ) s.t P Ω ( A ) = P Ω ( M ) (1.3) ∗ Facult´e des Sciences et Techniques-Gueliz, Laboratoire de Math´ematiques Appliqu´ees et Infor-matique, Marrakech, Morocco † Universit´e du Littoral Cote d’Opale, LMPA, 50 rue F. Buisson, 62228 Calais-Cedex, France ‡ laboratory MSDA, Mohammed VI Polytechnic University, Green City, Morocco.1 here A is the underlying tensor of order 3, M is the observed tensor, Ω is the set ofthe known data and P Ω is the projection operator defined by P Ω ( A ) i,j,k = (cid:26) A i,j,k , ( i, j, k ) ∈ Ω0 otherwise. As stated in [22], the optimization problem (1.3) is N-P hard and then one shouldstudy the tensor version of the problem (1.2) and this is the main subject of the presentwork. The new tensor-rank optimisation problem will be solved using the tensor T-product which is based on the Fast Fourier Transform (FFT). Notice that the notionof tensor rank is complicated as compared to the matrix rank and many tensor rankdefinitions and procedures such as Tucker-rank [4], CP-rank [4], TT-rank [29] and thetensor tubal rank [7], have been introduced the last years; see also [15, 25, 5, 20].The outline of this paper is as follows: In Section 2 we give some notations and prelim-inaries that will be used in the paper. Section 3 is devoted to the development of ourproposed tensor completion approaches. We will show how to use the tensor nuclearnorm in combination with the TV regularisation procedure to derive the new comple-tion algorithms. In the two approaches, we will use the well known tensor T-product.Some numerical experiments with comparisons to the most well known methods arepresented in Section 5, showing the effectiveness of the presented approaches.
2. Notations and preliminaries.
In this paper we denote tensors by calli-graphic letters, e.g., A . Matrices are denoted by capital letters, e.g., A and vectorsare denoted by lower case letters, e.g., a .Let A ∈ R n × n × n be an 3-order tensor. We define its Frobenius-norm by kAk F = vuut n X i =1 n X i =1 n X i =1 a i ,i ,i . The inner product between the two tensors A and B in R n × n × n is given by hA , Bi = n X i =1 n X i =1 n X i =1 a i ,i ,i b i ,i ,i . Let
A ∈ R n × n × n , then the i th frontal slice of the tensor A is the matrix A (: , : , i ) andwill be denoted by A ( i ) . For two matrices A ∈ R n × m and B ∈ R p × q , the Kroneckerproduct is the np × mq matrix given as A ⊗ B = [ a ij B ] i =1: n ; j =1: m ∈ R np × mq . Let v ∈ R n , we denote its DiscreteFourier Transform (DFT) by ˆ v and it is defined byˆ v = F n v where F n denotes the DFT matrix and it is defined by F n = . . . . . . ω ω ω . . . ω n − ω n ... ... ... ... ... ... ...1 ω n − ω n − . . . . . . ω ( n − n − ω ( n − ∈ R n × n (2.1) otice that F n √ n is a unitary matrix, i.e., F ∗ n F n = F n F ∗ n = nI n . The Fast Fourier Transform (FFT) allows us to compute the matrix-vector F n v in avery economical way. It computes this product with a cost of O ( nlogn ) instead of O (cid:0) n (cid:1) and it is represented in Matlab by the command f f t and its inverse by if f t :ˆ v = fft ( v ) and v = ifft (ˆ v ).The circulant matrix associated to the vector v is given as circ ( v ) = v v n v n − . . . v v v v n . . . v ... ... . . . . . . ... v n v n − . . . v v ∈ R n × n which it can be diagonalized by using the DFT and we get F n circ ( v ) F − n = diag (ˆ v )with diag (ˆ v ) denotes the diagonal matrix, where the i th element of its diagonal is ˆ v i . Lemma 1. [31] Given a real vector v ∈ R n , the associated ˆ v = F n v satisfies ˆ v ∈ R and conj (ˆ v i ) = ˆ v n − i +2 , i = 2 , ..., (cid:20) n + 12 (cid:21) . (2.2) Conversely, for any given complex ˆ v ∈ C n satisfying (2.2) , there exists a real circulantmatrix circ ( v ) such that (2.1) holds. Let
A ∈ R n × n × n be a 3-order tensor, we denote its DFT along each tubes ˆ A ∈ C n × n × n , . This operation can be done in Matlab by using the following commandˆ A = fft ( A , [ ] , . Conversely, we can obtain A from ˆ A using the Matlab command A = ifft (cid:16) ˆ A , [ ] , (cid:17) . Thanks to Lemma 1, we haveˆ A (1) ∈ R n × n and conj (cid:16) ˆ A ( i ) (cid:17) = ˆ A ( n − i +2) f or i = 2 , ..., (cid:20) n + 12 (cid:21) . We have also kAk F = 1 √ n (cid:13)(cid:13)(cid:13) ˆ A (cid:13)(cid:13)(cid:13) F and hA , Bi = 1 n D ˆ A , ˆ B E . (2.3)We define the block diagonal matrix associated to the tensor A ∈ R n × n × n as follows bdiag ( A ) = A (1) A (2) A (3) . . . A ( n ) (2.4) lso, we define its block circulant matrix by bcirc ( A ) = A (1) A ( n ) . . . . . . A (2) A (2) A (1) . . . . . . A (3) ... . . . . . . . . . ...... . . . . . . . . . ... A ( n ) A ( n − . . . . . . A (1) . (2.5)As bcirc ( A ) is a block circulant matrix, it can be block diagonalized using the DFT[8]. Then we get ( F n ⊗ I n ) bcirc ( A ) (cid:0) F ∗ n ⊗ I n (cid:1) = bdiag ( ˆ A ) . (2.6) Let
A ∈ R n × n × n a third-order tensor, wedefine the following operators unfold ( A ) = h A (1) , A (2) , . . . , A ( n ) i T , fold ( unfold ( A )) = A . Definition 1.
T-product [26]Let
A ∈ R n × n × n and B ∈ R n × n × n , we define the t-product between A and B by A ∗ B = f old ( bcirc ( A ) unfold ( B )) . (2.7)We notice that from (2.6), we can compute the T-product between two tensors A and B of appropriate sizes using the following property D = A ∗ B ⇐⇒ bdiag ( ˆ D ) = bdiag ( ˆ A ) bdiag ( ˆ B ) . The following algorithm summarises the different steps for the T-product
Algorithm 1
The T-product via the FFT. Inputs:
A ∈ R n × n × n and B ∈ R n × n × n Output: C = A ∗ B ∈ R n × n × n Compute ˆ A = fft ( A , [ ] ,
3) and ˆ B = fft ( B , [ ] , for i = 1 to (cid:20) n + 12 (cid:21) do ˆ C (: , : , i ) = ˆ A (: , : , i ) ˆ B (: , : , i ) end for for i = (cid:20) n + 12 (cid:21) + 1 to n do ˆ C (: , : , i ) = conj (cid:16) ˆ C (: , : , n + 2 − i ) (cid:17) end for ˆ C (: , : , i ) = ˆ A (: , : , i ) ˆ B (: , : , i ) C = ifft ( ˆ C , [ ] , A ∈ R n × n × n and B ∈ R n × n × n is O (cid:16) nn n n (cid:17) instead of O (cid:0) n n n n (cid:1) if we use directly the relation(2.7). .3. The tensor SVD. In the sequel, we need the following definitions.
Definitions 1. [26] • Conjugate transpose:
The conjugate transpose of a tensor
A ∈ C n × n × n is the tensor A ∗ ∈ C n × n × n obtained by conjugate transposing each of itsfrontal slices and then reversing the order of transposed frontal slices through n . • Identity tensor:
The identity tensor
I ∈ R l × l × n is the tensor with its firstfrontal slice being the l × l identity matrix, and other frontal slices being allzeros. • F-diagonal tensor:
A tensor is called f-diagonal if each of its frontal slicesis a diagonal matrix. • Orthogonal tensor:
A tensor
Q ∈ R l × l × n is orthogonal if it satisfies Q ∗ Q ∗ = Q ∗ ∗ Q = I . The Singular Value Decomposition (SVD) for matrices, was generalized to thetensor case using the T-product [26] as is stated in the following theorem.
Theorem 2.
Let
A ∈ R n × n × n real valued tensor, then A can be factored as A = U ∗ S ∗ V T with U ∈ R n × n × n and V ∈ R n × n × n are orthogonal tensors, and S ∈ R n × n × n is an f -diagonal tensor. The process called T-SVD of decomposing a 3-order tensor via the tensor T-productis summarized in the following algorithm
Algorithm 2
T-SVD Impute
A ∈ R n × n × n Output: t-SVD components U , S and V ˆ A = fft ( A , [ ] , for i = 1 to (cid:20) n + 12 (cid:21) do [ ˆ U (: , : , i ) , ˆ S (: , : , i ) , ˆ V (: , : , i )] = SV D ( ˆ A (: , : , i )) end for for i = (cid:20) n + 12 (cid:21) + 1 to n do ˆ U (: , : , i ) = conj ( ˆ U (: , : , n + 2 − i )) ˆ S (: , : , i ) = conj ( ˆ S (: , : , n + 2 − i )) ˆ V (: , : , i ) = conj ( ˆ V (: , : , n + 2 − i )) end for U = ifft ( ˆ U , [ ] , S = ifft ( ˆ S , [ ] ,
3) and V = ifft ( ˆ V , [ ] , Definition 3.
For
A ∈ R n × n × n , the tensor tubal rank, denoted as rank t ( A ) ,is the number of nonzero singular tubes of S , where S is from the t-SVD of A = U ∗ S ∗ V ∗ . We can write rank t ( A ) = card ( { i/ S ( i, i, :) = 0 } ) . (2.8) efinition 4. For
A ∈ R n × n × n , the tensor average rank, denoted as rank a ( A ) ,is defined as rank a ( A ) = 1 n rank ( bcirc ( A )) . (2.9) We first recall the matrix nuclear norm. Let A ∈ R n × m be a matrix, then the nuclear norm denoted by k A k ∗ , is defined as thedual norm of the matrix spectral norm, i.e. k A k ∗ = arg min k B k≤ |h A, B i| (2.10)where k B k denotes the matrix spectral norm. We notice that we also have k A k ∗ = r X i =1 σ i (2.11)where { σ i } ri =1 are the singular values of A and r is the tank of A . Definition 5.
Let
A ∈ R n × n × n be a 3-order tensor. Then the tensor spectralnorm of A is defined as kAk = k bcirc ( A ) k (2.12)Also, from (2.4) and (2.5), we get kAk = (cid:13)(cid:13)(cid:13) bdiag ( ˆ A ) (cid:13)(cid:13)(cid:13) . (2.13)The tensor nuclear norm is an extension of the matrix nuclear norm to tensors and isdefined in the following definition Definition 6.
Let
A ∈ R n × n × n be a three mode-tensor. The tensor nuclearnorm of A is defined as follows kAk ∗ = arg min kBk≤ |hA , Bi| . (2.14)Using (2.3) and (2.14), we get the following relations kAk ∗ = 1 n k bcirc ( A ) k ∗ = 1 n (cid:13)(cid:13)(cid:13) bdiag ( ˆ A ) (cid:13)(cid:13)(cid:13) ∗ , (2.15)and kAk ∗ = r X i =1 S ( i, i,
1) (2.16)where r is the tubal rank of A , and S is given from the T-SVD of A . Theorem 7. [12]On the set {A ∈ R n × n × n / kAk ≤ } the convex envelope of the average rank of A is its tensor nuclear norm kAk ∗ . .5. Tensor singular value thresholding . In this subsection, we first recallthe Tensor singular value thresholding [12] and give an algorithm summarizing thewhole process that will be used later.
Definition 8.
Let
A ∈ R n × n × n be a tensor and consider its tensor SVDdecomposition as A = U ∗ S ∗ V T . The tensor singular value thresholding of A with a given parameter τ is defined by D τ ( A ) = U ∗ S τ ∗ V T (2.17) where S τ = if f t (cid:16) max (cid:16) ˆ S − τ, (cid:17) , [ ] , (cid:17) . Theorem 9. [12] For any τ > and A ∈ R n × n × n , the tensor singular valuethresholding (2.17) is connected to the nuclear norm via the following relation D τ ( A ) = arg min X ∈ R n × n × n τ kX k ∗ + kA − X k F . (2.18)The tensor singular value thresholding process is summarized in the following algo-rithm Algorithm 3
Tensor singular value thresholding algorithm Impute
A ∈ R n × n × n Output: D τ ( A ) ˆ A = fft ( A , [ ] , for i = 1 to (cid:20) n + 12 (cid:21) do [ ˆ U (: , : , i ) , ˆ S (: , : , i ) , ˆ V (: , : , i )] = SV D ( ˆ A (: , : , i )) ˆ S ( i ) τ = (cid:16) ˆ S ( i ) − τ (cid:17) + D τ (cid:16) ˆ A (cid:17) ( i ) = ˆ U (: , : , i ) ˆ S (: , : , i ) τ ˆ V (: , : , i ) end for for i = (cid:20) n + 12 (cid:21) + 1 to n do D τ (cid:16) ˆ A (cid:17) ( i ) = D τ (cid:16) ˆ A (cid:17) ( n +2 − i ) end for D τ ( A ) = ifft ( D τ (cid:16) ˆ A (cid:17) , [ ] ,
3. The proposed methods.
Our proposed approaches are based on the min-imization of the average rank of a three-order tensor
A ∈ R n × n × n . The problemcan be formulated as follows min A rank a ( A ) s.t P Ω ( A ) = P Ω ( M ) . (3.1)where A ∈ R n × n × n is the underlying tensor, M ∈ R n × n × n is the observed tensorand Ω is the set of the known pixels. As stated in [1] and thanks to Theorem 7, we an replace the problem (3.1) by the following one min A kAk ∗ s.t P Ω ( A ) = P Ω ( M ) (3.2)where k . k ∗ is the tensor nuclear norm defined above. It is known that in real world,the problem (3.2) can be very ill-conditioned and one needs regularisation techniquessuch as the ( TV )-regularisation which we will consider in the present paper. As inthe matrix case, other regularization procedures are also possible. For our first approach, we consider the TV-regularized problem min A kAk ∗ + λ TV ( A ) s.t P Ω ( A ) = P Ω ( M ) (3.3)where λ is a regularization parameter, and TV ( A ) = h T V ( A (1) ) | T V ( A (2) | ... | T V ( A ( n ) ) i ∈ R n × n × n , with T V ( A ( n ) ) = n X i =1 n X j =1 q(cid:0) D A ( n ) (cid:1) i,j + (cid:0) D A ( n ) (cid:1) i,j , n ∈ { , , ..., n } and D and D are the derivative operators in the first and the second direction,respectively, with D A ( n ) = A ( n ) C n , D A ( n ) = C n A ( n ) where C and C are the matrices defined as C m = − . . . − . . . . . .
00 1 − . . . ...... . . . . . . . . . . . . ...... . . . . . . . . . . . . ...0 0 0 . . . − ∈ R m × m , C p = − . . . − . . . . . .
00 0 − . . . ...... . . . . . . . . . . . . ...... . . . . . . . . . . . . 11 0 0 . . . − ∈ R p × p . To solve the constrained optimization problem (3.3), we have to go through the fol-lowing intermediate optimization problem min A , Z , Y kZk ∗ + λ n X n =1 n X j =1 n X i =1 (cid:13)(cid:13)(cid:13) Y ( n ) i,j (cid:13)(cid:13)(cid:13) (3.4) s.t P Ω ( A ) = P Ω ( M ) , Z = A , Y = D A and Y = D A with Y ( n ) i,j = h ( Y ) ( n ) i,j , ( Y ) ( n ) i,j i for n ∈ { , , ..., n } , i ∈ { , , ..., n } and j ∈{ , , ..., n } , D X = (cid:2) D X (1) | D X (2) | ... | D X ( n ) | (cid:3) and D X = (cid:2) D X (1) | D X (2) | ... | D X ( n ) | (cid:3) . he constrained optimization problem (3.4) can be written as min A , Z , W [ F ( Z ) + G ( W )] (3.5) s.t P Ω ( A ) = P Ω ( M ) , Z = A , W = D A where F ( Z ) = kZk ∗ , G ( W ) = λ n X n =1 n X j =1 n X i =1 (cid:13)(cid:13)(cid:13) Y ( n ) i,j (cid:13)(cid:13)(cid:13) , D = (cid:18) D D (cid:19) and W = (cid:18) W W (cid:19) = (cid:18) Y Y (cid:19) ∈ R n × n × n . To solve the regularized optimisation problem (3.5), we can use the well knownADM method [21, 33]. The augmented Lagrangian associated to the problem (3.5) isgiven by L ( A , Z , W , Q , B ) = F ( Z ) + G ( W ) + hA − Z , Qi + β kA − Zk F + hD A − W , Bi + β kD A − Wk F (3.6)with Q ∈ R n × n × n , B = (cid:18) B B (cid:19) ∈ R n × n × n are the Lagrangian multipliers and β , β > (cid:0) A k , Z k , W k (cid:1) = arg min X , Z , W L ( A , Z , W , Q k , B k ) , (3.7) Q k +1 = Q k + β ( A k − Z k ) , (3.8) B k +1 = B k + β (cid:0) D A k − W k (cid:1) . (3.9)Let us see now how to solve each of those sub-problems. • Solving the A -problem : From (3.7) for a given Z , W , we compute theapproximation A k by solving for A the minimization problem A k = arg min A β (cid:13)(cid:13)(cid:13)(cid:13) A − Z + Q k β (cid:13)(cid:13)(cid:13)(cid:13) F + β (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) D A − W + (cid:0) B k (cid:1) β (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) F Then, the optimal value A k satisfies the following equation β A k + β (cid:0) D (cid:1) T D A k + β (cid:0) D (cid:1) T D A k = R (3.10)with R = β Z − Q k + β (cid:0) D (cid:1) T W − (cid:0) D (cid:1) T B k + β (cid:0) D (cid:1) T W − (cid:0) D (cid:1) T B k . To solve the tensor equation (3.10), we can transform it to a matrix equationby considering the frontal slices of the tensor A k . Then, for n ∈ { , , ..., n } ,the matrix (cid:0) A k (cid:1) ( n ) satisfies the following matrix equation β (cid:0) A k (cid:1) ( n ) + β ( D ) T D (cid:0) A k (cid:1) ( n ) + β ( D ) T D (cid:0) A k (cid:1) ( n ) = R ( n ) (3.11) hich can be also written as β (cid:0) A k (cid:1) ( n ) + β (cid:0) A k (cid:1) ( n ) (cid:0) C n (cid:1) T C n + β (cid:0) C n (cid:1) T C n (cid:0) A k (cid:1) ( n ) = R ( n ) (3.12)Since C n and C n are circulant matrices, they are diagonalizable using thediscrete Fourier transformation, i.e., there exist Λ and Λ diagonal matricessuch that C n = F ∗ n Λ F n , C n = F ∗ n Λ F n where F n and F n are the respectively the matrices representing the discreteFourier transformation of size n × n and n × n , then (cid:0) C n (cid:1) T C n = F ∗ n Λ F n , (cid:0) C n (cid:1) T C n = F ∗ n Λ F n By referring to [25], we can rewrite (3.12), for each n ∈ { , , ..., n } , as (cid:0) F ∗ n ⊗ F ∗ n (cid:1) (cid:0) β I ⊗ I + β Λ ⊗ I + β I ⊗ Λ (cid:1) ( F n ⊗ F n ) vect (cid:16)(cid:0) A k (cid:1) ( n ) (cid:17) = vect (cid:16) R ( n ) (cid:17) and since β I ⊗ I + β Λ ⊗ I + β I ⊗ Λ is an invertible matrix we get that vect (cid:16)(cid:0) A k (cid:1) ( n ) (cid:17) = (cid:0) F ∗ n ⊗ F ∗ n (cid:1) (cid:0) β I ⊗ I + β Λ ⊗ I + β I ⊗ Λ (cid:1) − ( F n ⊗ F n ) vect ( R ( n ) ) . (3.13)As the parameters β and β are strictly positive numbers, this shows thatfor each n ∈ { , , ..., n } , the equation (3.12) has a unique solution. • Solving the Z -problem : Given X and W , the value of Z k satisfies thefollowing optimization problem Z k = arg min Z F ( Z ) + β (cid:13)(cid:13)(cid:13)(cid:13) Z − A − Q k β (cid:13)(cid:13)(cid:13)(cid:13) F . = arg min Z kZk ∗ + β (cid:13)(cid:13)(cid:13)(cid:13) Z − A − Q k β (cid:13)(cid:13)(cid:13)(cid:13) F . Then, from the result of Theorem 9, we get Z k = D τ (cid:18) A + Q k β (cid:19) (3.14)with τ = 1 β . • Solving the W -problem : For a given X and Z , W k is obtained by solvingthe following sub-problems: for n ∈ { , , ..., n } the n th sub-problem is givenby (cid:0) W k (cid:1) ( n ) = arg min W n λ n X j =1 n X i =1 (cid:13)(cid:13) Y ni,j (cid:13)(cid:13) + β (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) Y ( n )1 − D A ( n ) − (cid:0) B k (cid:1) n β (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) F + β (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) Y ( n )2 − D A ( n ) − (cid:0) B k (cid:1) n β (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) F . hich is equivalent to solve the following Y -problem (cid:0) Y k (cid:1) ( n ) i,j = arg min Y ni,j λ (cid:13)(cid:13)(cid:13) Y ( n ) i,j (cid:13)(cid:13)(cid:13) + β (cid:16) Y ( n )1 (cid:17) i,j − (cid:16) D A ( n ) (cid:17) i,j − (cid:0) B k (cid:1) ( n ) i,j β + β (cid:16) Y ( n )2 (cid:17) i,j − (cid:16) D A ( n ) (cid:17) i,j − (cid:0) B k (cid:1) ( n ) i,j β . By using the 2-D shrinkage formula, we will get, for 1 ≤ i ≤ n and 1 ≤ j ≤ n , the following expression (cid:0) Y k (cid:1) ( n ) i,j = max { (cid:13)(cid:13)(cid:13) S ( n ) i,j (cid:13)(cid:13)(cid:13) − λβ , } S ( n ) i,j (cid:13)(cid:13)(cid:13) S ( n ) i,j (cid:13)(cid:13)(cid:13) (3.15)where S ( n ) i,j = (cid:0) D A ( n ) (cid:1) i,j + (cid:0) B k (cid:1) ( n ) i,j β , (cid:0) D A ( n ) (cid:1) i,j + (cid:0) B k (cid:1) ( n ) i,j β and we set0 (cid:18) (cid:19) = 0 . The tensor completion procedure using the tensor nuclear norm and first order totalvariation (TNN-TV ) is summarized in Algorithm 4. Algorithm 4
Tensor completion using the tensor nuclear norm and the first ordertotal variation (TNN-TV ). Initialize Z , W , Q , B , λ , β and β . while not converged do Update A k from (3.13) ; Update Z k from (3.14); Update W k from (3.15); Update Q k +1 from (3.8); Update B k +1 from(3.9); end while3.2. Tensor completion using tensor nuclear norm and the second ordertotal variation. In this part we apply the total variation with the second orderderivative. Then second proposed model is formulated as min A kAk ∗ + λ TV ( A ) s.t P Ω ( A ) = P Ω ( M ) (3.16)where A , M k . k ∗ and Ω play the same role as in the first model (3.3). The expressionof TV is given by TV ( A ) = h T V ( A (1) ) | T V ( A (2) ) | ... | T V ( A ( n ) ) i ∈ R n × n × n where T V (cid:16) A ( n ) (cid:17) = n X i =1 n X j =1 q(cid:0) D A ( n ) (cid:1) i,j + (cid:0) D A ( n ) (cid:1) i,j or each n ∈ { , , ..., n } . The matrices D and D are the second derivative operatorsin the first and in the second direction, respectively , satisfying for n ∈ { , , ..., n } D A ( n ) = A ( n ) C n and D A ( n ) = C n A ( n ) with C i = 12 − ... − ...
00 . . . . . . . . . . . . ...... . . . . . . . . . . . . ...0 ... ... − ... ... − ∈ R i × i ; i = n , n . The optimization problem (3.16) is equivalent to the following one min A , Z , Y kZk ∗ + λ n X n =1 n X j =1 n X i =1 (cid:13)(cid:13)(cid:13) Y ( n ) i,j (cid:13)(cid:13)(cid:13) s.t P Ω ( A ) = P Ω ( M ) , Z = A , Y = D A and Y = D A (3.17)with Y ( n ) i,j = h ( Y ) ni,j , ( Y ) ( n ) i,j i for n ∈ { , , ..., n } , i ∈ { , , ..., n } and j ∈ { , , ..., n } , D A = (cid:2) D A (1) | D A (2) | .... | D A ( n ) (cid:3) and D A = (cid:2) D A (1) | D A (2) | .... | D A ( n ) (cid:3) .The constrained optimization problem (3.17) is transformed to the following one arg min A , Z , W [ F ( Z ) + G ( W )] s.t P Ω ( A ) = P Ω ( M ) , A = Z , D A = W (3.18)where F ( Z ) = kZk ∗ , G ( W ) = λ n X n =1 n X j =1 n X i =1 (cid:13)(cid:13)(cid:13) Y ( n ) i,j (cid:13)(cid:13)(cid:13) , D = (cid:18) D D (cid:19) and W = (cid:18) W W (cid:19) = (cid:18) Y Y (cid:19) .We observe that (3.18) is similar to (3.5). Thus, we can use the same procedure ofsolving (3.5) to solve (3.18). The augmented Lagrangian associated to the optimiza-tion problem (3.18) is given by L ( A , Z , W, Q , B ) = F ( Z ) + G ( W ) + hA − Z , Qi + β kA − Zk F + hD A − W , Bi + β kD A − Wk F . (3.19)Therefore, using ADM, we have to solve the following sub-problems (cid:0) A k , Z k , W k (cid:1) = arg min A , Z , W L ( A , Z , W , Q k , B ) (3.20) Q k +1 = Q k + β (cid:0) A k − Z k (cid:1) (3.21) B k +1 = B k + β (cid:0) D A k − W k (cid:1) . (3.22)Let us see now how to solve each of those sub-problems. The A -problem: For fixed Z and W , each frontal slice of the approximation A k satisfies the Sylvester matrix equation β (cid:0) A k (cid:1) ( n ) + β ( D ) T D (cid:0) A k (cid:1) ( n ) + β ( D ) T D (cid:0) A k (cid:1) ( n ) = R ( n ) (3.23)where R = β Z − Q k + β (cid:0) D (cid:1) T W − (cid:0) D (cid:1) T B k + β (cid:0) D (cid:1) T W − (cid:0) D (cid:1) T B k ,which can be written as β (cid:0) A k (cid:1) ( n ) + β (cid:0) A k (cid:1) ( n ) C Tn C n + β C Tn C n (cid:0) A k (cid:1) ( n ) = R ( n ) . Using the same idea as for 3.12, each frontal slice of A k satisfies vect (cid:16)(cid:0) A k (cid:1) ( n ) (cid:17) = (cid:0) F ∗ n ⊗ F ∗ n (cid:1) (cid:0) β I ⊗ I + β Λ ⊗ I + β I ⊗ Λ (cid:1) − ( F n ⊗ F n ) vect ( R ( n ) ) . (3.24)with Λ = F n C n F ∗ n , and Λ = F n C n F ∗ n , where F n i is the Fourier matrix of size n i × n i for i = 1 , • The Z -problem: For τ = 1 β and for a given A and W we get Z k = D τ (cid:18) A + Q k β (cid:19) . (3.25) • The W -problem: By applying 2D shrinkage formula on each frontal sliceof Y for a given A and Z , we get (cid:0) Y k (cid:1) ( n ) i,j = max (cid:26)(cid:13)(cid:13)(cid:13) S ( n ) i,j (cid:13)(cid:13)(cid:13) − λβ , (cid:27) S ( n ) i,j (cid:13)(cid:13)(cid:13) S ( n ) i,j (cid:13)(cid:13)(cid:13) (3.26)where S ( n ) i,j = (cid:0) D A ( n ) (cid:1) i,j + (cid:0) B k (cid:1) ( n ) i,j β , (cid:0) D A ( n ) (cid:1) i,j + (cid:0) B k (cid:1) ( n ) i,j β with 0 (cid:18) (cid:19) =0 . The different steps of the tensor completion using the tensor nuclear norm and totalvariation (TNN-TV ) is summarized in the following algorithm. Algorithm 5
Tensor completion using the tensor nuclear norm and the second ordertotal variation (TNN-TV ). Initialize Z , W , Q , B , λ , β and β . while not converged do Update A k from (3.24) ; Update Z k from (3.25); Update W k from (3.26); Update Q k +1 from (3.21); Update B k +1 from (3.22); end while Now we discuss the complexity of the two algorithms TNN-TV and TNN-TV . As we used the fast Fourier transform, the cost of computing the A -sub-problem is O ( n n n log ( n n )). The cost of computing Z in (3.14) and (3.25) is O ( n n n + n n )). Computing W in (3.15) and (3.26) requires O ( n n n ) arithmetic operations. . Convergence analysis. In this section we study the convergence of the pro-posed approaches. As the two methods are similar, we will give theoretical resultsonly for sequences obtained by Algorithm 4. Notice first that the functions F and G defined earlier are closed, proper and convex. Then, thanks to [18][30], the optimiza-tion problem (3.5) is solvable, i.e., there exist Z ∗ and W ∗ not necessarily unique thatminimize (3.5).Let us define the space E = R n × n × n × R n × n × n × R n ×× n × n × R n ×× n × n × R n ×× n × n which is closed and nonempty. We first recall the following theorem. Theorem 10. [34] A ∗ is a solution of (3.3) if and only if there exist ( Z ∗ , W ∗ ) ∈ R n × n × n × R N × n × n and ( Q ∗ , B ∗ )) ∈ R n × n × n × R N × n × n such that ( A ∗ , Z ∗ , W ∗ , Q ∗ , B ∗ ) ∈ E is a saddle point of L ,i.e. L ( A ∗ , Z ∗ , W ∗ , Q , B ) ≤ L ( A ∗ , Z ∗ , W ∗ , Q ∗ , B ∗ ) ≤ L ( A , Z , W , Q ∗ , B ∗ ) ∀ ( A , Z , W , Q , B ) ∈ E . (4.1)The next theorem gives some convergence results on the sequences obtained fromAlgorithm 4. Theorem 11.
Assume that ( A ∗ , Z ∗ , W ∗ , Q ∗ , B ∗ ) is a saddle point of L . Thesequence ( X k , Z k , W k , Q k , B k ) generated by Algorithm 4 satisfies:1. lim k → + ∞ F ( Z k ) + G ( W k ) = F ( Z ∗ ) + G ( W ∗ ) . lim k → + ∞ (cid:13)(cid:13) A k − Z k (cid:13)(cid:13) = 0 . lim k → + ∞ (cid:13)(cid:13) D A k − W k (cid:13)(cid:13) = 0 . Proof 1.
From the first inequality of (4.1) we get hA ∗ − Z ∗ , Qi + hD A ∗ − W ∗ , Bi ≤ hA ∗ − Z ∗ , Q ∗ i + hD A ∗ − W ∗ , B ∗ i ∀ ( Q , W )(4.2) which gives (cid:26) A ∗ = Z ∗ D A ∗ = W ∗ . Let us define the following quantities ¯ X k = A k − A ∗ , ¯ Z k = Z k − Z ∗ , ¯ W k = W k − W ∗ , ¯ Q k = Q k − Q ∗ , ¯ B k = B k − B ∗ . The main idea of the proof is to show that the sequence (cid:16) β (cid:13)(cid:13) ¯ Q k (cid:13)(cid:13) F + β (cid:13)(cid:13) ¯ B k (cid:13)(cid:13) F (cid:17) k ≥ is decreasing. Notice that ¯ Q k +1 = ¯ Q k − β (cid:0) ¯ A k − ¯ Y k (cid:1) , ¯ W k +1 = ¯ W k − β (cid:0) D ¯ A k − ¯ B k (cid:1) . We have (cid:16)(cid:13)(cid:13) ¯ Q k (cid:13)(cid:13) F + (cid:13)(cid:13) ¯ B k (cid:13)(cid:13) F (cid:17) − (cid:16)(cid:13)(cid:13) ¯ Q k +1 (cid:13)(cid:13) F + (cid:13)(cid:13) ¯ B k +1 (cid:13)(cid:13) F (cid:17) = − β (cid:10) ¯ Q k , ¯ A k − ¯ Y k (cid:11) − β (cid:10) ¯ B k , D ¯ A k − ¯ W k (cid:11) − β (cid:13)(cid:13) ¯ A k − ¯ Y k (cid:13)(cid:13) F − β (cid:13)(cid:13) D ¯ A k − ¯ W k (cid:13)(cid:13) F . (4.3) From the second inequality of (4.1) we obtain the following inequalities for ( A , Z , W ) =( A k , Z k , W k ) (cid:10) A k − Z ∗ , Q ∗ (cid:11) + (cid:10) D A k − W ∗ , B ∗ (cid:11) + β (cid:10) A k − Z ∗ , A k − A ∗ (cid:11) + β (cid:10) D A k − W ∗ , D ( A k − A ∗ ) (cid:11) ≥ ( Z k ) − F ( Z ∗ ) + (cid:10) A ∗ − Y k , Q ∗ (cid:11) + β (cid:10) A ∗ − Y k , Y ∗ − Y k (cid:11) ≥ .G ( W k ) − F ( W ∗ ) + (cid:10) D A ∗ − W k , B ∗ (cid:11) + β (cid:10) D A ∗ − W k , W ∗ − W k (cid:11) ≥ . Using (3.7) , we get for ( A , Z , W ) = ( A k , Z k , W k ) that (cid:10) A ∗ − Z k , Q k (cid:11) + (cid:10) D A ∗ − W k , B k (cid:11) + β (cid:10) A ∗ − Z k , Z ∗ − Z k (cid:11) + β (cid:10) D A ∗ − W k , W ∗ − W k (cid:11) ≥ F ( Z ∗ ) − F ( Z k ) + (cid:10) A k − Y ∗ , Q k (cid:11) + β (cid:10) X k − Y ∗ , A k − A ∗ (cid:11) ≥ .G ( W ∗ ) − G ( W k ) + (cid:10) D A k − W ∗ , B k (cid:11) + β (cid:10) D A k − W ∗ , D ( A k − A ∗ ) (cid:11) ≥ . By regrouping terms we get − (cid:10) ¯ Q k , ¯ A k − ¯ Y k (cid:11) − (cid:10) ¯ B k , D ¯ A k − ¯ W k (cid:11) ≥ β (cid:13)(cid:13) ¯ A k − ¯ Y k (cid:13)(cid:13) F + β (cid:13)(cid:13) D ¯ A k − ¯ W k (cid:13)(cid:13) F then (cid:16) β (cid:13)(cid:13) ¯ Q k (cid:13)(cid:13) F + β (cid:13)(cid:13) ¯ B k (cid:13)(cid:13) F (cid:17) − (cid:16) β (cid:13)(cid:13) ¯ Q k +1 (cid:13)(cid:13) F + β (cid:13)(cid:13) ¯ B k +1 (cid:13)(cid:13) F (cid:17) ≥ β β (cid:13)(cid:13) ¯ A k − ¯ Y k (cid:13)(cid:13) F + β β (cid:13)(cid:13) D ¯ A k − ¯ W k (cid:13)(cid:13) F . Thus, the sequence (cid:16) β (cid:13)(cid:13) ¯ Q k (cid:13)(cid:13) F + β (cid:13)(cid:13) ¯ B k (cid:13)(cid:13) F (cid:17) k ≥ is decreasing, which gives + ∞ X k =0 (cid:16) β β (cid:13)(cid:13) ¯ A k − ¯ Y k (cid:13)(cid:13) F + β β (cid:13)(cid:13) D ¯ A k − ¯ W k (cid:13)(cid:13) F (cid:17) ≤ β (cid:13)(cid:13) ¯ Q (cid:13)(cid:13) F + β (cid:13)(cid:13) ¯ W (cid:13)(cid:13) F . Therefore (cid:0) Q k (cid:1) k ≥ and (cid:0) B k (cid:1) k ≥ are boundedlim k → + ∞ (cid:13)(cid:13) A k − Z k (cid:13)(cid:13) = 0 .lim k → + ∞ (cid:13)(cid:13) D A k − W k (cid:13)(cid:13) = 0 . In addition, by using again the second inequality of (4.1) for ( A , Z , W ) = ( A k , Z k , W k ) we obtain F ( Z ∗ ) + G ( W ∗ ) ≤ F ( Z k ) + G ( W k ) + (cid:10) A k − Z k , Q k (cid:11) + β (cid:13)(cid:13) A k − Z k (cid:13)(cid:13) F + (cid:10) D A k − W k , B k (cid:11) + β (cid:13)(cid:13) D A k − W k (cid:13)(cid:13) F , and F ( Z ∗ ) + G ( W ∗ ) ≥ F ( Z k ) + G ( W k ) + (cid:10) A k − Z k , Q k (cid:11) + β (cid:13)(cid:13) A k − Z k (cid:13)(cid:13) F + (cid:10) D A k − W k , B k (cid:11) + β (cid:13)(cid:13) D A k − W k (cid:13)(cid:13) F . Hence lim sup k → + ∞ F ( Z k ) + G ( W k ) ≤ F ( A ∗ ) + G ( W ∗ ) ≤ lim inf k → + ∞ F ( Z k ) + G ( W k ) . . Numerical experiments. In this section, we give some numerical tests toshow the performance of our proposed algorithms TNN-TV and TNN-TV and com-pare them with the results obtained by other known algorithms for image and videocompletion, such as TNN [1], SiLRTC-TT [5] and MF-TV [25].The quality of the recovered images is measured by computing the relative squarederror (RSE), and the peak signal-to-noise-ration (PSNR), defined by RSE = kA ori − Ak F kAk F , and P SN R = 10 log M ax A kA − A ori k F , where A ori is the original tensor, A is the recovered tensor and M ax A is the maximumpixel value of A . The convergence stopping criterion is defined by computing therelative error of A between two successive iterations as follows (cid:13)(cid:13) A k +1 − A k (cid:13)(cid:13) F kA k k F ≤ − . (5.1)In all the experiments we used fixed values of regularization and penalty parameters.For the algorithm TNN-TV we used λ = 0 . , β = 0 .
01 and β = 0 . we set λ = 1, β = 0 .
01 and β = 0 . For this example, we used color images of size 256 × ×
3. InFigure 5.1, we reported the obtained visual results of TNN, SiLRTC-TT, MF-TV,TNN-TV and TNN-TV , with SR = 0 . SR represents the percentage ofthe data remained in the image. In Table 5.1 we compared the efficiency of our twoalgorithms with TNN, SiLRTC-TT and MF-TV by comparing the values of RSE and
P SN R . igure 5.1: The results of the algorithms TNN, MF-TV, ,SiLRTC-TT, TNN-TV andTNN-TV for the images ”Lena”, ”Barbara”, ”airplane” and ”house” for SR = 0 . gives better results thanTV .In Table 5.2, we reported the execution times needed to achieve the convergencecriterion for each method. As can be seen from this table, the results obtained byTV and TV are faster as compared to the ones obtained by the other three methods. SR algorithms RSE PSNR RSE PSNR RSE PSNR RSE PSNR0.1 TNN 0.1942 19.2585 0.2362 18.8644 0.1415 18.2747 0.1596 20.2903SiLRTC-TT 0.1521 21.4342 0.2020 20.2838 0.1338 19.4385 0.1451 20.9854MF-TV 0.3345 14.6340 0.3867 14.4866 0.1946 16.1082 0.1867 18.8811TNN-TV TNN-TV Table 5.1: The values of RSE and PSNR for TNN, SiLRTC-TT, MF-TV, TNN-TV and TNN-TV with the images ”Lena”, ”Barbara”, ”airplane” and ”house” using SR = 0 . , . . SR 0.2 0.3RSE 0.13 0.11 0.09 0.08 0.1 0.09 0.08 0.07TNN 3.0392 - - - 1.6736 - - -TNN-TV , TNN-TV , MF-TV andSiLRTC-TT for SR = 0 . , . In this part we test the performance of our algorithms on somevideos. In our example, we used the video of ”Suzie” of size 128 × × SR = 0 . SR = 0 .
05 for the second one. igure 5.2: The results of the algorithms TNN, SiLRTC-TT, MF-TV, TNN-TV andTNN-TV for a frame of the video ”Suzie” with SR = 0 . SR = 0 .
05 for the second lineFigure 5.3: The values of RSE by the algorithms SiLRTC-TT, TNN, MF-TV, TNN-TV and TNN-TV for each frame of the video of ”Suzie” for SR = 0 . SR = 0 . and TNN-TV . As can be seen from this figure, TNN-TV and TNN-TV return the best results. In this subsection we test our methods on the MRI data of the frontdirection. In this example, we used a video of MRI of front direction of size 181 × × SR = 0 . igure 5.4: The results of the algorithms TNN, MF-TV, TNN-TV and TNN-TV for the video of the front direction with SR = 0 . and TNN-TV .In Figure 5.5 we plotted the values of RSE for each frame of the video MRI obtainedwith TNN, MF-TV, TNN-TV and TNN-TV . As one can see from this figure, TNN-TV returns the best result.
6. Conclusion.