改进的低秩张量逼近算法

Enhanced Low Rank Tensor Approximation Algorithm

摘要: 本文在低秩矩阵逼近(Low-rank Matrix Approximation, LRMA)问题的基础上将其推广到张量上，提出了利用非凸正则化构造凸优化问题来估计低秩张量的方法，采用参数化的非凸罚函数估计非零的奇异值，并求出了该目标函数的全局最优解。实验结果表明，该方法能很好处理低秩张量逼近问题，实现图像去噪。

Abstract: In this paper, we extended the low rank matrix approximation problem to tensors. A method for estimating low-rank tensors by constructing convex optimization problems with non-convex reg-ularization is proposed. The non-zero singular values are estimated by parameterized non-convex penalty functions, and the global optimal solution of the objective function is obtained. The ex-perimental results show that this method can deal with the low rank tensor approximation problem very well and achieve image denoising.

1. 引言

近年来，低秩逼近广泛应用于信号处理 [1] [2] [3] 和机器学习 [4] [5] 等各个领域，研究人员提出了大量的方法，这些算法大致分为两种，凸优化算法和非凸优化算法。凸优化算法，主要利用核范数极小化。这一思路由文献 [6] 提出，主要包括：文献 [7] 提出的奇异值阈值(Singular Value Thresholding, SVT)算法。这种算法对原始矩阵的不同奇异值的惩罚力度相同，所求出的全局最优解在实际问题中效果差。所以本文提出了一个新的算法，即利用非凸正则化构造凸优化问题来估计低秩张量。该方法不仅能很好处理低秩张量逼近问题，还能实现图像去噪。

2. 预备知识

定义1 (张量的切片)对于一个三阶张量 ${\mathcal{A}}_{ijk}$ ，用 $\mathcal{A}\left(i,:,:\right),\mathcal{A}\left(:,i,:\right),\mathcal{A}\left(:,:,i\right)$ 分别表示张量 $\mathcal{A}$ 第i个水平的，侧面的和正面的切片，通常记正面切片 $A^{\left(i\right)}=\left(:,:,i\right)$ 。

定义2 (离散傅里叶变换)设张量 $\mathcal{A}\in R^{n_1\times n_2\times n_3}$ 沿着第三维度的傅里叶变换得到 $\bar{\mathcal{A}}\in C^{n_1\times n_2\times n_3}$ ，在傅里叶变换域中块循环矩阵可以被映射到傅里叶变换域的一个块对角矩阵。

记 块循环矩阵为

$bcirc\left(\mathcal{A}\right)=\left[\begin{array}{cccc}{A}^{\left(1\right)}& A^{\left(n_3\right)}& \cdots & A^{\left(2\right)}\\ A^{\left(2\right)}& A^{\left(1\right)}& \cdots & A^{\left(3\right)}\\ ⋮& ⋮& \ddots & ⋮\\ A^{\left(n_3\right)}& A^{\left(n_3-1\right)}& \cdots & A^{\left(1\right)}\end{array}\right]\in C^{n_1{n}_3\times n_2{n}_3}$ (1)

记 块对角矩阵为

$\bar{A}=bdiag\left(\bar{\mathcal{A}}\right)=\left[\begin{array}{cccc}{\bar{A}}^{\left(1\right)}& & & \\ & {\bar{A}}^{\left(2\right)}& & \\ & & \ddots & \\ & & & {\bar{A}}^{\left(n_3\right)}\end{array}\right]\in C^{n_1{n}_3\times n_2{n}_3}$ (2)

定义3 [8] (张量的核范数)令 $\mathcal{A}=\mathcal{U}\ast \mathcal{S}\ast {\mathcal{V}}^{\ast }$ ，则有

${\Vert \mathcal{A}\Vert }_{\ast }=\langle \mathcal{S},\mathcal{I}\rangle ={\displaystyle \underset{i=1}{\overset{r}{\sum }}\mathcal{S}\left(i,i,1\right)}$ (3)

其中 $\mathcal{I}$ 表示单位张量， $r=ran{k}_t\left(\mathcal{A}\right)$ 。

假设1 [9] ：非凸性惩罚函数 $\varphi :R\to R$ 满足以下条件：

1) $\varphi$ 在R上连续并且 $\varphi$ 在 $R/\left\{0\right\}$ 上是连续可微的且为偶函数即 $\varphi \left(-x;a\right)=\varphi \left(x;a\right)$

2) ${\varphi }^{\prime }\left(x;a\right)\ge 0,\forall a>0$

3) ${\varphi }^{″}\left(x;a\right)\le 0,\forall a>0$ (4)

4) ${\varphi }^{\prime }\left(0^{-};a\right)=1$

5) $\underset{X\ne 0}{\inf }{\varphi }^{″}\left(x;a\right)={\varphi }^{″}\left(0^{+};a\right)=-a$

定理1 [8] (张量奇异值分解，T-SVD)令 $\mathcal{A}\in R^{n_1\times n_2\times n_3}$ ，可以将其分解为

$\mathcal{A}=\mathcal{U}\ast \mathcal{S}\ast {\mathcal{V}}^{\ast }$ (5)

其中 $\mathcal{U}\in R^{n_1\times n_1\times n_3}$ 和 $\mathcal{V}\in R^{n_2\times n_2\times n_3}$ 都是正交张量， $\mathcal{S}\in R^{n_1\times n_2\times n_3}$ 是F-对角张量。

定理 2 [9] 令 $Y=U\Sigma V^{\text{T}}$ 为Y的奇异值SVD分解，假设 $\varphi$ 是 $R\to R$ 上的满足假设1的非凸惩罚函数，当 $0\le a<\frac{1}{\lambda }$ 时， $\frac{1}{2}{\Vert Y-X\Vert }_F^2+\lambda {\displaystyle \underset{i=1}{\overset{k}{\sum }}\varphi }\left({\sigma }_i\left(X\right),a\right)$ 的全局最小值为

$\bar{X}=U\cdot \Theta \left(\Sigma ;\lambda ,a\right)\cdot V^{\text{T}}$ (6)

其中 $\Theta \left(y;\lambda ,a\right):=\arg \underset{x\in R}{\min }\left\{\frac{1}{2}{\left(y-x\right)}^2+\lambda \varphi \left(x;a\right)\right\}$ 是函数 $\varphi$ 的近端算子。

3. 改进的低秩张量逼近算法

3.1. 模型的建立

设 $\mathcal{Y}$ 为受到噪声污染后的低秩张量， $\mathcal{X}$ 为潜在的目标低秩张量， $\mathcal{E}$ 为高斯噪声张量。则有 $\mathcal{Y}=\mathcal{X}+\mathcal{E}$ 。所以低秩张量逼近模型如下：

$\underset{\mathcal{X}}{\arg \min }\left\{\Psi \left(\mathcal{X}\right)=\frac{1}{2}{\Vert \mathcal{Y}-\mathcal{X}\Vert }_F^2+\lambda {\displaystyle \underset{i=1}{\overset{k}{\sum }}\varphi \left(\mathcal{S}\left(i,i,1\right),a\right)}\right\}$ (7)

本文假设所有张量都是三阶张量即 $\mathcal{X}\in R^{n_1\times n_2\times n_3},\mathcal{Y}\in R^{n_1\times n_2\times n_3}$ 。其中 $k=\min \left(n_1{n}_3,n_2{n}_3\right)$ ， $\varphi$ 是稀疏产生的正则化函数，所以它可能是非凸的。 $\mathcal{S}\left(i,i,1\right)$ 表示张量 $\mathcal{X}$ 奇异值分解后F-对角张量的第一个正面切面。

3.2. 模型的求解

由式(1)可得 ${\displaystyle \underset{i=1}{\overset{r}{\sum }}\mathcal{S}\left(i,i,1\right)}={\Vert \mathcal{A}\Vert }_{\ast }$ 则有 $\lambda {\displaystyle \underset{i=1}{\overset{k}{\sum }}\varphi \left(\mathcal{S}\left(i,i,1\right),a\right)}=\lambda \varphi \left({\Vert \mathcal{A}\Vert }_{\ast },a\right)$ (8)

可以得到目标函数的如下形式：

$\frac{1}{2}\left[\frac{1}{n_3}{\Vert \bar{Y}\Vert }_F^2+{\Vert \bar{X}\Vert }_F^2-2\cdot \frac{1}{n_3}\langle \bar{Y},\bar{X}\rangle \right]+\lambda \varphi \left({\Vert \mathcal{X}\Vert }_{\ast },a\right)$ (9)

由张量核范数与矩阵核范数之间的关系，所以目标函数可以写成如下形式：

$\begin{array}{l}\frac{1}{2n_3}{\Vert \bar{Y}-\bar{X}\Vert }_F^2+\frac{1}{n_3}\lambda \varphi \left({\Vert \bar{X}\Vert }_{\ast },a\right)\\ =\frac{1}{n_3}\left[\frac{1}{2}{\Vert \bar{Y}-\bar{X}\Vert }_F^2+\lambda {\displaystyle \underset{i=1}{\overset{k}{\sum }}\varphi }\left({\sigma }_i\left(\bar{X}\right),a\right)\right]\end{array}$ (10)

式中 $\bar{X},\bar{Y}\in R^{n_1{n}_3\times n_2{n}_3}$ ， ${\sigma }_i$ 表示矩阵 $\bar{X}$ 的第i个奇异值， $k=\min \left(n_1{n}_3,n_2{n}_3\right)$ 。

根据定理2可以得到目标函数的最优解为

$\frac{1}{n_3}\left(U\cdot \Theta \left(\Sigma ;\lambda ,a\right)\cdot V^{\text{T}}\right)$ (11)

其中 $\Theta \left(\Sigma ;\lambda ,a\right):=\arg \underset{X\in R}{\min }\left\{\frac{1}{2}{\left(\Sigma -X\right)}^2+\lambda \varphi \left(X;a\right)\right\}$

4. 实验

对于一个噪声图像 $\mathcal{Y}$ ，将其进行奇异值分解后得到的F-对角张量的第一个正面切片记为Y，通过块匹配方法可以得到非局部自相似块矩阵 $Y_j$ ，则有 $Y_j=X_j+W_j$ ，其中 $X_j$ 和 $W_j$ 分别表示清晰图像和高斯噪声。通过叠加这些非局部的相似块形成张量，为了从 ${\mathcal{Y}}_j$ 中提取出 ${\mathcal{X}}_j$ ，我们构造了如下目标函数：

${\hat{\mathcal{X}}}_j:=\underset{\mathcal{X}}{\arg \min }\left\{\frac{1}{2}{\Vert {\mathcal{Y}}_j-\mathcal{X}\Vert }_F^2+\lambda {\displaystyle \underset{i=1}{\overset{m}{\sum }}\varphi \left(S\left(i,i,1\right);a\right)}\right\}$ (10)

根据定理2可得到(10)的最优解为 ${\hat{\mathcal{X}}}_j=\frac{1}{n_3}\left(U_j\cdot \Theta \left({\Sigma }_j;\lambda ,a\right)\cdot V_j^{\text{T}}\right)$ 。

Figure 1. Comparison of denoising effects of different methods

图1. 不同方法去噪效果对比

图1显示了使用不同的方法对噪声图像进行去噪的结果。表1列出了三种测试图像的平均PSNR值和运行时间。所提出的LRTA方法与BM3D、PS和WNNM方法相比可以获得更高的PSNR。在图像质量方面，所提出的LRTA方法也更为清晰一些。

Table 1. PSNR values and run time data

表1. PSNR值和运行时间数据

5. 结论

本文提出了一种新的方法来估计低秩张量。尽管正则化项不是凸的，但所提出的目标函数是凸的。通过对输入张量进行离散的傅里叶变换后再进行奇异值分解，得到F-对角张量的第一个正面切面形成的矩阵的奇异值进行阈值化，得到目标函数的非迭代解。在基于非局部自相似的图像去噪算法中，证明了该方法的有效性。