笔记：Sparse and Truncated Nuclear Norm Based Tensor Completion

Han, Z. F., Leung, C. S., Huang, L. T., & So, H. C. (2017). Sparse and Truncated Nuclear Norm Based Tensor Completion. Neural Processing Letters, 45(3), 729-743.
本文是这篇 Neural Processing Letter 期刊论文的笔记，主要是对文中的理论方法进行展开详解。本人学术水平有限，文中如有错误之处，敬请指正。

摘要： 张量（tensor）补全的主要问题是计算其秩（rank）。最近，一种张量核范数（tensor nuclear norm），等于所有张量的展开的矩阵的核范数的加权和，已经被提出来解决这个问题。然而，在此方法中， 所有的奇异值是被同时最小化的。所以，张量的秩并没有被很好地近似。另外，许多现有的方法都忽视了张量的结构信息。此文提出了一种张量补全的算法，基于张量截断核范数（truncated nuclear norm），其优于传统的核范数方法。而且，为了保持结构的信息，一个稀疏约束项，定义于变换域，也被添加到了目标函数中。实验结果表明此文的方法能胜过一些最先进的算法，在张量补全中。

1 简介

在信号处理和机器学习的应用中，比如图像去噪，图像分类和子空间分割，需要根据已知的数据来估丢失的元素。当然，没有对数据信息的先验知识，估计丢失的元素是一个病态的问题。所以，一些假设被广泛地采用。例如，可以使用输入数据的统计或结构的信息，来建立已知数据和未知数据的关系。然而，许多方法都只关注局部联系。

另一种办法是使用输入数据的全局结构信息。许多研究表明真实的信号都存在低维的空间之中。比如，许多的自然图像都包含一些纹理区域。这些纹理的秩通常都很小。所以，许多的矩阵补全方法都研究最小化矩阵的秩。不幸的是，矩阵的秩并不是一个凸函数，意味着低秩最小化不能被有效地求解。

最近，研究表明，核范数，一个凸近似，适合低秩最小化求解。另外，核范数最小化可以通过一些数学方法有效地求解。然而，如 1 所述，使用核范数的最大限制是其不能很好地近似矩阵的秩，因为在优化过程中，所有的奇异值都被同时最小化了。在 2 中，矩阵截断核范数（truncated nunclear norm, MTNN）被提出。不像之前的核范数方法，最小化所有奇异值的和，该 MTNN 方法关注一部分最小的奇异值。经验性的研究表明此 MTNN 方法更优于其他矩阵核范数方法。

最近，作为矩阵补全的扩展，低秩张量补全越来越受关注。Liu et al. 3 首次定义了张量核范数，将矩阵补全扩展到张量补全中。另外，他们构建了两种算法，快速低秩张量补全（FaLRTC）和高精度低秩张量补全（HaLRTC）。由于张量核范数的概念是基于矩阵核范数，所以张量的秩也并没有被很好地近似。另外，这两种算法并没有考虑张量数据地结构信息。

此文定义了张量截断核范数（TTNN）方法。之后，将其应用到张量补全问题中。为了更好地提升性能，此文将 $\ell_1$ 范数约束项加入到目标函数中，保留数据片之间的平滑性质。约束项等于多维离散余弦变换系数的和。

2 背景

2.1 符号

此文中，矩阵用大写加粗字母表示（$\mathbf{X}, \mathbf{Y}$），其元素加下标（$\mathbf{X}_{ij}$）。Frobenius 范数定义为 $|| \mathbf{X} ||_\text{F} = \left( \sum_{ij} |\mathbf{X}_{ij}|^2 \right)^{\frac{1}{2}}$ 。矩阵的内积定义为 $\langle \mathbf{X}, \mathbf{Y} \rangle := \sum_{ij} \mathbf{X}_{ij} \mathbf{Y}_{ij}$ 。张量定义为书法体的大写加粗字母（$\mathcal{X}, \mathcal{Y}$）。一个 $N$ 模的张量，或 $N$ 阶的张量，表示为 $\mathcal{X} \in \mathbb{R}^{I_1 \times I_2 \times \cdots \times I_N}$，其中 $I_m$ 表示第 $m$ 个维度的大小。张量中的一个元素表示为 $\mathcal{X}_{i_1 \cdots i_N}$ 。模 $N$ 矩阵化一个张量 $\mathcal{X}$ 的结果是一个矩阵，表示为 $\mathcal{X}_{(n)}$ 。其过程是将除 $i_n$ 之外的索引全部线性化。矩阵 $\mathcal{X}_{(n)}$ 的大小为 $I_n \times \prod_{k=1,k \neq n}^{N} I_k$ 。在解包过程中，一个元素 $\mathcal{X}_{i_1 \cdots i_N}$ 对应着 $\mathcal{X}_{(n)}$ 的位置 $(i_n, j)$，其中

$$\begin{equation} j = 1 + \sum_{k=1, \, k \neq n}^{N} (i_k - 1) J_k \ \ \text{with} \ \ J_k = \prod_{m=1, \, m \neq n}^{k-1} I_m . \end{equation}$$
其逆操作打包定义为 $\text{Fold}_{(n)} \left( \mathcal{X}_{(n)} \right) := \mathcal{X}$ 。张量 $\mathcal{X}$ 的 Frobenius 范数定义为 $|| \mathcal{X} ||_\text{F} := \left( \sum_{i_1, \cdots, i_N} | \mathcal{X}_{i_1 \cdots i_N} |^2 \right)^{\frac{1}{2}}$，类似地，两个相同大小的 tensor 的内积定义为 $\langle \mathcal{X}, \mathcal{Y} \rangle := \sum_{i_1, \cdots, i_N} \mathcal{X}_{i_1 \cdots i_N} \mathcal{Y}_{i_1 \cdots i_N}$ 。

2.2 矩阵补全

首先回顾以下矩阵补全的概念。令 $\mathbf{M} \in \mathbb{R}^{m \times n}$ 为一个低秩的不完整的矩阵。令 $\Omega$ 为已知数据元素的集合，而 $\Omega^\text{c}$ 是其补集。矩阵补全可以用如下的优化问题表示：

$$\begin{align} \min_{\mathbf{X}} \ & \text{rank} (\mathbf{X}) \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathbf{X}) = {P}_{\Omega} (\mathbf{M}) , \tag{1} \end{align}$$
其中 $P_{\Omega} (\cdot)$ 是映射操作，给出如下
$$\begin{equation} \left( P_{\Omega} (\mathbf{M}) \right)_{ij} = \begin{cases} \mathbf{M}_{ij}, & (i, j) \in \Omega , \\ 0, & (i, j) \in \Omega^\text{c}. \\ \end{cases} \end{equation}$$
在优化问题中 $(1)$，目标是估计 $\mathbf{X}$ 中的缺失的元素，使其的秩尽可能的小。一般的，矩阵的范数最小化问题是一个 NP-hard 问题。

在文献 4 中，介绍了矩阵 $\mathbf{X} \in \mathbb{R}^{m \times n}$ 的核范数，给出如下

$$\begin{equation} || \mathbf{X} ||_* = \sum_{i=1}^{\min(m,n)} \sigma_i , \end{equation}$$
其中 $\sigma_i$ 是矩阵 $\mathbf{X}$ 是第 $i$ 个最大的奇异值。使用矩阵核范数的概念，矩阵补全的问题 $(1)$ 就转变为
$$\begin{align} \min_{\mathbf{X}} \ & || \mathbf{X} ||_* \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathbf{X}) = {P}_{\Omega} (\mathbf{M}) . \tag{2} \end{align}$$
在基于矩阵核范数的方法中，所有的奇异值都是同时被最小化的。所以，矩阵的秩没有被很好地近似 5 。为了解决这个问题，MTNN 方法被提出了。其等价于 $\min(m,n) - r$ 个最小的奇异值之和，即
$$\begin{equation} || \mathbf{X} ||_r = \sum_{i=r+1}^{\min(m,n)} \sigma_i . \end{equation}$$

接着，优化问题 $(2)$ 转变为

$$\begin{align} \min_{\mathbf{X}} \ & || \mathbf{X} ||_r \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathbf{X}) = {P}_{\Omega} (\mathbf{M}) . \tag{3} \end{align}$$
根据文献 [5] 所述，该优化问题 $(3)$ 可以改写为
$$\begin{align} \min_{\mathbf{X}} \ & || \mathbf{X} ||_* - \max_{\mathbf{A} \mathbf{A}^\text{T} = \mathbf{I}, \, \mathbf{B} \mathbf{B}^\text{T} = \mathbf{I}} \text{tr} \left( \mathbf{A} \mathbf{X} \mathbf{B}^\text{T} \right) \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathbf{X}) = {P}_{\Omega} (\mathbf{M}) , \tag{4} \end{align}$$
其中 $\mathbf{A} \in \mathbb{R}^{r \times m}$， $\mathbf{B} \in \mathbb{R}^{r \times n}$， $\mathbf{I}$ 是一个单位矩阵。在文献 [5] 中，作者提出了算法使用交替乘子法（alternating direction method of multipliers, ADMM）求解优化问题。基本思想是设置初始值 $\mathbf{X}^{(1)} = P_{\Omega} (\mathbf{M})$ 。在第 $t$ 次迭代过程中，先固定 $\mathbf{X}^{(t)}$ 和进行奇异值分解（SVD）于 $\mathbf{X}_t$，即 $\mathbf{X}^{(t)} = \mathbf{U} \Sigma \mathbf{V}^\text{T}$，其中 $\mathbf{U} = \left[ \mathbf{u}_1, \cdots, \mathbf{u}_m \right] \in \mathbb{R}^{m \times m}$， $\Sigma \in \mathbb{R}^{m \times n}$， $\mathbf{V} = \left[ \mathbf{v}_1, \cdots, \mathbf{v}_n \right] \in \mathbb{R}^{n \times n}$ 。令 $\mathbf{A}^{(t)} = \left[ \mathbf{u}_1, \cdots, \mathbf{u}_r \right]^\text{T}$ 和 $\mathbf{B}^{(t)} = \left[ \mathbf{v}_1, \cdots, \mathbf{v}_r \right]^\text{T}$ 。之后，固定 $\mathbf{A}^{(t)}$ 和 $\mathbf{B}^{(t)}$，接着更新变量 $\mathbf{X}$ 基于如下形式：
$$\begin{align} \min_{\mathbf{X}} \ & || \mathbf{X} ||_* - \text{tr} \left( \mathbf{A}^{(t)} \mathbf{X} \mathbf{B}^{(t)\text{T}} \right) \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathbf{X}) = {P}_{\Omega} (\mathbf{M}) . \tag{5} \end{align}$$
基于如上两步骤，交替更新矩阵，最终会使算法收敛到局部最小值。

2.3 张量补全

张量补全是矩阵补全的一般化问题。给定一个低秩的不完整的张量 $\mathcal{T} \in \mathbb{R}^{I_1 \times \cdots \times I_N}$，令 $\Omega$ 对应已观测到的元素的集合，而 $\Omega^{\text{c}}$ 是其相应的补集。在张量补全中，需要估计一个低秩的 $\mathcal{X} \in \mathbb{R}^{I_1 \times \cdots \times I_N}$ 并满足 $P_{\Omega} (\mathcal{X}) = P_{\Omega} (\mathcal{T})$ 。在文献 [3] 中，低秩张量补全问题给出如下

$$\begin{align} \min_{\mathcal{X}} \ & \text{rank}(\mathcal{X}) \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathcal{X}) = {P}_{\Omega} (\mathcal{T}) . \tag{6} \end{align}$$
不像矩阵的情况，很难定义一个张量的秩。许多的张量的秩的版本提出 6 7 。通过将张量转变为矩阵使用矩阵化， Liu et al. [3] 提出了基于张量补全的张量核范数，如下
$$\begin{equation} || \mathcal{X} ||_* := \sum_{i=1}^n \alpha_i || \mathcal{X}_{(i)} ||_* , \end{equation}$$
其中 $\alpha_i \geq 0$ 且 $\sum_{i=1}^n \alpha_i = 1$ 。低秩张量补全问题 $(6)$ 定义如下
$$\begin{align} \min_{\mathcal{X}} \ & || \mathcal{X} ||_* \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathcal{X}) = {P}_{\Omega} (\mathcal{T}) . \tag{7} \end{align}$$
另外，Liu et al. [3] 提出了 3 种方法求解张量补全问题，包括简单版（SiLRTC），快速版（FaLRTC），和高精度版（HaLRTC）。由于矩阵核范数的内部依赖性，介绍了虚拟的矩阵作为中间变量， $\mathbf{M}_i, \ i = 1,2, \cdots, n$ 。将优化问题 $(7)$ 进行转化
$$\begin{align} \min_{\mathcal{X}, \ \mathbf{M}_i} \ & \sum_{i=1}^n \alpha_i || \mathbf{M}_{(i)} ||_* \\ \text{s.t.} \ \ & {P}_{\Omega} (\mathcal{X}) = {P}_{\Omega} (\mathcal{T}) , \tag{8} \\ & \mathcal{X}_{(i)} = \mathbf{M}_i, \ i = 1, 2, \cdots, n . \end{align}$$
然而，在 FaLRTC 中，所有的奇异值都是被同时最小化的。所以秩也没有被很好地近似。另外，结构信息通常也被忽略了。

2.4 ADMM

交替乘子法 8 是一个经典的算法，用于求解带约束的最小化问题。考虑如下的约束最小化问题：

$$\begin{align} \min_{\mathbf{x}, \ \mathbf{z}} \ & f(\mathbf{x}) + g (\mathbf{z}) \\ \text{s.t.} \ \ & \mathbf{E} \mathbf{x} + \mathbf{F} \mathbf{z} = \mathbf{c} , \tag{9} \end{align}$$
其中 $\mathbf{x} \in \mathbb{R}^{n_x}$， $\mathbf{z} \in \mathbb{R}^{n_z}$， $\mathbf{E} \in \mathbb{R}^{n_c \times n_x}$， $\mathbf{F} \in \mathbb{R}^{n_c \times n_z}$， $f(\cdot)$ 和 $g(\cdot)$ 均是凸函数。ADMM 算法考虑一个增广 Lagrangian 函数，如下：
$$\begin{equation} L(\mathbf{x}, \mathbf{z}, \mathbf{y}) = f(\mathbf{x}) + g(\mathbf{z}) + \langle \mathbf{y} , \mathbf{E} \mathbf{x} + \mathbf{F} \mathbf{z} - \mathbf{c} \rangle + \frac{\beta}{2} || \mathbf{E} \mathbf{x} + \mathbf{F} \mathbf{z} - \mathbf{c} ||^2 , \end{equation}$$
其中 $\mathbf{y}$ 是 Lagrangian 乘子向量。ADMM 迭代机制是：
$$\begin{align} \mathbf{x}^{(k+1)} &= \arg\min_{\mathbf{x}} \ L(\mathbf{x}, \mathbf{z}^{(k)}, \mathbf{y}^{(k)}), \\ \mathbf{z}^{(k+1)} &= \arg\min_{\mathbf{z}} \ L(\mathbf{x}^{(k+1)}, \mathbf{z}, \mathbf{y}^{(k)}), \\ \mathbf{y}^{(k+1)} &= \mathbf{y}^{(k)} + \beta \left( \mathbf{E} \mathbf{x}^{(k+1)} + \mathbf{F} ( \mathbf{y}^{(k+1)} ) - \mathbf{c} \right) , \end{align}$$
其中 $\beta$ 是一个正的惩罚参数。

3 张量截断核范数补全

由于 Matrix TNN 优于矩阵的核范数，此文中将 Matrix TNN 扩展到 Tensor TNN，即

$$\begin{equation} || \mathcal{X} ||_r := \sum_{i=1}^n \alpha_i || \mathcal{X}_{(i)} ||_r . \end{equation}$$
另外，空间结构信息可以用变换系数的稀疏性表征。所以，此文添加一个稀疏约束项（定义于变换域）到目标函数中，保持空间的结构信息。如此，根据之前的 TNN 的定义，可以将目标优化问题 $(4)$ 转化为
$$\begin{align} \min_{\mathcal{X},\ \mathcal{W}} \ & \left( \sum_{i=1}^n \alpha_i || \mathcal{X}_{(i)} ||_* - \alpha_i \max_{\mathbf{A}_i \mathbf{A}_i^\text{T} = \mathbf{I}, \ \mathbf{B}_i \mathbf{B}_i^\text{T} = \mathbf{I} } \ \text{tr} \left( \mathbf{A}_i \mathcal{X}_{(i)} \mathbf{B}_i^\text{T} \right) \right) + \lambda || \mathcal{W} ||_0 \\ \text{s.t.} \ \ & P_\Omega (\mathcal{X}) = P_\Omega (\mathcal{T}), \ \mathcal{W} = \mathcal{F} (\mathcal{X}) , \tag{10} \end{align}$$
其中 $\mathcal{F}(\cdot)$ 是一个正向的 $n$ 维离散余弦变换，$\lambda > 0$ 是一个加权的常系数。由于 $\ell_0$ 范数是一个非凸的函数，此文使用 $\ell_1$ 范数来近似测量 $\mathcal{W}$ 的稀疏性。该约束优化问题 $(10)$ 可以转化为
$$\begin{align} \min_{\mathcal{X},\ \mathcal{W}} \ & \left( \sum_{i=1}^n \alpha_i || \mathcal{X}_{(i)} ||_* - \alpha_i \max_{\mathbf{A}_i \mathbf{A}_i^\text{T} = \mathbf{I}, \ \mathbf{B}_i \mathbf{B}_i^\text{T} = \mathbf{I} } \ \text{tr} \left( \mathbf{A}_i \mathcal{X}_{(i)} \mathbf{B}_i^\text{T} \right) \right) + \lambda || \mathcal{W} ||_1 \\ \text{s.t.} \ \ & P_\Omega (\mathcal{X}) = P_\Omega (\mathcal{T}), \ \mathcal{W} = \mathcal{F} (\mathcal{X}) . \tag{11} \end{align}$$
因为 $\mathcal{X}_{(i)}$ 之间有内部依赖性，所以这里需要引入一个虚拟的变量 $\mathcal{M}_i$ 加入目标函数 $(11)$ 中，来分离这些有依赖的项，即
$$\begin{align} \min_{\mathcal{X},\ \mathcal{W}} \ & \left( \sum_{i=1}^n \alpha_i || \mathcal{M}_{i,(i)} ||_* - \alpha_i \max_{\mathbf{A}_i \mathbf{A}_i^\text{T} = \mathbf{I}, \ \mathbf{B}_i \mathbf{B}_i^\text{T} = \mathbf{I} } \ \text{tr} \left( \mathbf{A}_i \mathcal{M}_{i,(i)} \mathbf{B}_i^\text{T} \right) \right) + \lambda || \mathcal{W} ||_1 \\ \text{s.t.} \ \ & P_\Omega (\mathcal{X}) = P_\Omega (\mathcal{T}), \ \mathcal{W} = \mathcal{F} (\mathcal{X}), \ \mathcal{X} = \mathcal{M}_i, \ i =1, 2, \cdots, n, \tag{12} \end{align}$$
其中 $\mathcal{M}_{i,(i)}$ 是 $\mathcal{M}_i$ 的模 $i$ 的矩阵化的结果。这里定义此文的方法 $(12)$ 记为 TTNNL1。当忽视约束项的时候，即 $\lambda = 0$，将方法记为 TTNN。

此算法的基本思想是给定一个当前的估计 $\mathcal{X}$，基于 SVD，可以获得 $\{ \mathbf{A}_i, \mathbf{B}_i \}, \ i =1, 2, \cdots, N$ 。之后，可以更新 $\mathcal{X}$ 基于如下的优化问题

$$\begin{align} \min_{\mathcal{X},\ \mathcal{W}} \ & \left( \sum_{i=1}^n \alpha_i || \mathcal{M}_{i,(i)} ||_* - \alpha_i \, \text{tr} \left( \mathbf{A}_i \mathcal{M}_{i,(i)} \mathbf{B}_i^\text{T} \right) \right) + \lambda || \mathcal{W} ||_1 \\ \text{s.t.} \ \ & P_\Omega (\mathcal{X}) = P_\Omega (\mathcal{T}), \ \mathcal{W} = \mathcal{F} (\mathcal{X}), \ \mathcal{X} = \mathcal{M}_i, \ i =1, 2, \cdots, n. \tag{13} \end{align}$$

4 张量补全算法

4.1 算法

基于 $(13)$，可以定义如下的增广 Lagrangian 函数

$$\begin{align} L(\mathcal{X}, \mathcal{M}_i, \mathcal{W}, \mathcal{Y}_i, \mathcal{B}) = & \sum_{i=1}^n \alpha_i || \mathcal{M}_{i,(i)} ||_* - \alpha_i \, \text{tr} \left( \mathbf{A}_i \mathcal{M}_{i,(i)} \mathbf{B}_i^\text{T} \right) \\ & + \langle \mathcal{X} - \mathcal{M}_i, \mathcal{Y}_i \rangle + \frac{\beta}{2} || \mathcal{X} - \mathcal{M}_i ||_\text{F}^2 + \langle \mathcal{W} - \mathcal{F}(\mathcal{X}), \mathcal{B} \rangle \\ & + \frac{\beta}{2} || \mathcal{W} - \mathcal{F} (\mathcal{X}) ||_\text{F}^2 + \lambda || \mathcal{W} ||_1. \tag{14} \end{align}$$
其中 $\mathcal{B}$ 和 $\mathcal{Y}_i$ 是 Lagrange 乘子张量，其维度与 $\mathcal{X}$ 一致， $\beta$ 是一个惩罚参数。使用 ADMM 的迭代机制，优化问题 $(13)$ 按照如下的迭代机制求解：
1. 计算 $\mathcal{X}^{(k+1)}$：
$$\begin{align} & \mathcal{X}' = \arg\min_{\mathcal{X}} \ L \left( \mathcal{X}, \, \mathcal{M}_{i,(i)}^{(k)}, \, \mathcal{W}^{(k)}, \, \mathcal{Y}_i^{(k)}, \, \mathcal{B}^{(k)} \right), \\ & \mathcal{X}_{\Omega^\text{c}}^{(k+1)} = \mathcal{X}_{\Omega^\text{c}} ', \tag{15} \\ & \mathcal{X}_{\Omega}^{(k+1)} = \mathcal{T}_{\Omega} . \end{align}$$
$\quad$ 这里后两步表示：需要保持已知的元素不变。
2. 计算 $\mathcal{M}_i^{(k+1)}$：
$$\begin{equation} \mathcal{M}_i^{(k+1)} = \arg\min_{\mathcal{M}_i} \ L \left( \mathcal{X}^{(k+1)}, \, \mathcal{M}_{i,(i)}, \, \mathcal{W}^{(k)}, \, \mathcal{Y}_i^{(k)}, \, \mathcal{B}^{(k)} \right). \tag{16} \end{equation}$$
3. 计算 $\mathcal{W}^{(k+1)}$：
$$\begin{equation} \mathcal{W}^{(k+1)} = \arg\min_{\mathcal{W}} \ L \left( \mathcal{X}^{(k+1)}, \, \mathcal{M}_{i,(i)}^{(k+1)}, \, \mathcal{W}, \, \mathcal{Y}_i^{(k)}, \, \mathcal{B}^{(k)} \right). \tag{17} \end{equation}$$
4. 计算 $\mathcal{Y}_i^{(k+1)}$：
$$\begin{equation} \mathcal{Y}_i^{(k+1)} = \mathcal{Y}_i^{(k)} + \beta \left( \mathcal{X}^{(k+1)} - \mathcal{M}_i^{(k+1)} \right). \tag{18} \end{equation}$$
5. 计算 $\mathcal{B}^{(k+1)}$：
$$\begin{equation} \mathcal{B}^{(k+1)} = \mathcal{B}^{(k)} + \beta \left( \mathcal{W}^{(k+1)} - \mathcal{F} ( \mathcal{X}^{(k+1)} ) \right). \tag{19} \end{equation}$$
算法总结于 Algrithm 1 中。具体的更新步骤 1，2，3 后面会给出。

---

Algorithm 1 基于 ADMM 优化
Input: 张量数据 $\mathcal{T}$, 已知元素的位置集合 $\Omega$ 。
Initialization: $P_\Omega (\mathcal{X}) = P_\Omega (\mathcal{T}), \ \mathcal{Y}_i = \mathbf{0}, \ \mathcal{M}_i = \mathbf{0}, \ \mathbf{A}_i = \mathbf{0}, \ \mathbf{B}_i = \mathbf{0}, \ \mathcal{W} = \mathbf{0}, \ \mathcal{B} = \mathbf{0}$, $\lambda, \ \beta, \ \alpha$ 。
while not converged do
   Step 1: 更新 $\mathcal{X}^{(k+1)}$ 用 $(15)$；
   Step 2: 更新 $\mathcal{M}_i^{(k+1)}$ 用 $(16)$；
   Step 3: 更新 $\mathcal{W}^{(k+1)}$ 用 $(17)$；
   $\mathcal{Y}_i^{(k+1)} = \mathcal{Y}_i^{(k)} + \beta \left( \mathcal{X}^{(k+1)} - \mathcal{M}_i^{(k+1)} \right)$；
   $\mathcal{B}^{(k+1)} = \mathcal{B}^{(k)} + \beta \left( \mathcal{W}^{(k+1)} - \mathcal{F} ( \mathcal{X}^{(k+1)} ) \right)$；
end while
Output: $\mathcal{X}$ 。

---

4.2 计算 $\mathcal{X}^{(k+1)}$

在计算 $(15)$ 中，需要求解如下的无约束优化问题：

$$\begin{align} \mathcal{X}' &= \arg\min_{\mathcal{X}} \ \sum_{i=1}^{n} \left\langle \mathcal{X} - \mathcal{M}_i^{(k)}, \mathcal{Y}_i^{(k)} \right\rangle + \frac{\beta}{2} \left|\left| \mathcal{X} - \mathcal{M}_i^{(k)} \right|\right|_\text{F}^2 \\ & \quad + \left\langle \mathcal{W}^{(k)} - \mathcal{F} (\mathcal{X}), \mathcal{B}^{(k)} \right\rangle + \frac{\beta}{2} \left|\left| \mathcal{W}^{(k)} - \mathcal{F} (\mathcal{X}) \right|\right|_\text{F}^2 \\ &= \arg\min_{\mathcal{X}} \ \sum_{i=1}^n \left|\left| \mathcal{X} - \mathcal{M}_i^{(k)} + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right|\right|_\text{F}^2 + \left|\left| \mathcal{W}^{(k)} - \mathcal{F} (\mathcal{X}) + \frac{\mathcal{B}^{(k)}}{\beta} \right|\right|_\text{F}^2 . \tag{20} \end{align}$$
然而，该问题并没有直接的闭式解。因为不能很容易地从 $\mathcal{F} (\mathcal{X})$ 中分离出 $\mathcal{X}$ 。所以，该问题需要被迭代地求解。这需要消耗很多时间，意味着需要一个有效的求解方法。Parseval 定理 9 指出信号的能量是不变的，如果是在单位变换下。根据这个定理和离散余弦变换的单位性质， $(20)$ 中的最后一项可以改写为

$$\begin{equation} \left|\left| \mathcal{W}^{(k)} - \mathcal{F} (\mathcal{X}) + \frac{\mathcal{B}^{(k)}}{\beta} \right|\right|_\text{F}^2 = \left|\left| g \left( \mathcal{W}^{(k)} - \mathcal{F} (\mathcal{X}) + \frac{\mathcal{B}^{(k)}}{\beta} \right) \right|\right|_\text{F}^2 = \left|\left| \mathcal{X} - g \left( \mathcal{W}^{(k)} + \frac{\mathcal{B}^{(k)}}{\beta} \right) \right|\right|_\text{F}^2 , \tag{21} \end{equation}$$
其中 $g(\cdot)$ 表示 $n$ 维的反离散傅里叶变换。注意的是 Frobenius 范数的平方即为测量张量数据的能量。将 $(21)$ 带入 $(20)$ 中，可以得到
$$\begin{equation} \mathcal{X} ' = \arg\min_{\mathcal{X}} \ \sum_{i=1}^n \left|\left| \mathcal{X} - \mathcal{M}_i^{(k)} + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right|\right|_\text{F}^2 + \left|\left| \mathcal{X} - g \left( \mathcal{W}^{(k)} + \frac{\mathcal{B}^{(k)}}{\beta} \right) \right|\right|_\text{F}^2 . \tag{22} \end{equation}$$
该问题是一个最小二乘问题，并且有闭式解：
$$\begin{equation} \mathcal{X} ' = \frac{1}{n+1} \left\{ \sum_{i=1}^n \left( \mathcal{M}_i - \frac{\mathcal{Y}_i^{(k)}}{\beta} \right) + g \left( \mathcal{W}^{(k)} + \frac{\mathcal{B}^{(k)}}{\beta} \right) \right\} . \tag{23} \end{equation}$$

4.3 计算 $\mathcal{M}_i^{(k+1)}$

在计算 $(16)$ 中，需要求解如下的无约束优化问题：

$$\begin{align} \mathcal{M}_i^{(k+1)} &= \arg\min_{\mathcal{M}_i} \ \alpha_i || \mathcal{M}_{i,(i)} ||_* - \alpha_i \, \text{tr} \left( \mathbf{A}_i \mathcal{M}_{i,(i)} \mathbf{B}_i^\text{T} \right) + \left\langle \mathcal{X}^{(k+1)} - \mathcal{M}_i, \mathcal{Y}_i^{(k)} \right\rangle + \frac{\beta}{2} \left|\left| \mathcal{X}^{(k+1)} - \mathcal{M}_i \right|\right|_\text{F}^2 \\ &= \arg\min_{\mathcal{M}_i} \ \alpha_i || \mathcal{M}_{i,(i)} ||_* + \frac{\beta}{2} \left|\left| \mathcal{X}^{(k+1)} - \mathcal{M}_i + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right|\right|_\text{F}^2 + \left\langle \mathcal{M}_{i,(i)} - \alpha_i \mathbf{A}_i \mathbf{B}_i^\text{T} \right\rangle \\ &= \arg\min_{\mathcal{M}_i} \ \alpha_i || \mathcal{M}_{i,(i)} ||_* + \frac{\beta}{2} \left|\left| \mathcal{M}_i - \left( \mathcal{X}^{(k+1)} + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right)_{(i)} - \frac{\alpha_i}{\beta} \mathbf{A}_i \mathbf{B}_i^\text{T} \right|\right|_\text{F}^2 .\tag{24} \end{align}$$
该问题是一般的核范数最小化问题。可以用如下的定理解决
Theorem 1 对于一个 $\epsilon > 0$ 和 $\mathbf{Y} \in \mathbb{R}^{m \times n}$，矩阵的核范数可以如下求解
$$\begin{equation} \mathcal{D}_{\frac{1}{\epsilon}} (\mathbf{Y}) = \arg\min_{\mathbf{X}} \ || \mathbf{X} ||_* + \frac{\epsilon}{2} || \mathbf{X} - \mathbf{Y} ||_\text{F}^2 , \end{equation}$$
其中 $\mathcal{D}_\tau$ 是一种奇异值收缩操作，基于矩阵 $\mathbf{Y}$ 的奇异值分解。令 $\mathbf{U} \mathbf{\Sigma} \mathbf{V}^\text{T}$ 为 $\mathbf{Y}$ 的奇异值分解，其中 $\mathbf{\Sigma} = \text{diag}(\{ \sigma_i \}_{1 \leq i \leq \min(m, n)})$ 。该 $\mathcal{D}_\tau$ 定义为
$$\begin{align} \mathcal{D}_\tau (\mathbf{Y}) &= \mathbf{U} \mathcal{D}_\tau (\mathbf{\Sigma}) \mathbf{V}^\text{T} , \\ \mathcal{D}_\tau (\mathbf{\Sigma}) &= \text{diag} (\{ \max (\sigma_i - \tau, 0) \}_{1 \leq i \leq \min(m,n) }). \end{align}$$
所以，从定理 1 中，可以得到 $\mathcal{M}_{i, (i)}^{(k+1)}$ 的闭式解
$$\begin{equation} \mathcal{M}_{i, (i)}^{(k+1)} = \mathcal{D}_{\frac{\alpha_i}{\beta}} \left( \left( \mathcal{X}^{(k+1)} + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right)_{(i)} + \frac{\alpha_i}{\beta} \mathbf{A}_i \mathbf{B}_i^\text{T} \right) . \tag{25} \end{equation}$$
于是有
$$\begin{equation} \mathcal{M}_i^{(k+1)} = \text{Fold}_{(i)} \left\{ \mathcal{D}_{\frac{\alpha_i}{\beta}} \left( \left( \mathcal{X}^{(k+1)} + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right)_{(i)} + \frac{\alpha_i}{\beta} \mathbf{A}_i \mathbf{B}_i^\text{T} \right) \right\} . \tag{26} \end{equation}$$
注意的是计算 $\mathcal{M}_i^{(k+1)}$ 的过程中，生成 $\left( \mathcal{X}^{(k+1)} + \frac{\mathcal{Y}_i^{(k)}}{\beta} \right)_{(i)} + \frac{\alpha_i}{\beta} \mathbf{A}_i \mathbf{B}_i^\text{T}$ 的左、右奇异向量。他们会被矩阵更新 $\mathbf{A}_i$ 和 $\mathbf{B}_i$ 在下一轮迭代中。使用 $r$ 来表示 $\frac{\alpha_i}{\beta}$，即 MTNN 是奇异值之和，其小于 $|| \mathbf{X} ||_r = \sum_{\sigma_i < \frac{\alpha_i}{\beta}} \sigma_i$ 。

4.4 计算 $\mathcal{W}^{(k+1)}$

为了更新 $\mathcal{W}$，需要计算如下优化问题：

$$\begin{align} \mathcal{W}^{(k+1)} &= \arg\min_{\mathcal{W}} \ \lambda || \mathcal{W} ||_1 + \left\langle \mathcal{W} - \mathcal{F} (\mathcal{X}^{(k+1)}), \mathcal{B}^{(k)} \right\rangle + \frac{\beta}{2} \Big|\Big| \mathcal{W} - \mathcal{F} (\mathcal{X}^{(k+1)}) \Big|\Big|_\text{F}^2 \\ &= \arg\min_{\mathcal{W}} \ \lambda || \mathcal{W} ||_1 + \frac{\beta}{2} \left|\left| \mathcal{W} - \mathcal{F} (\mathcal{X}^{(k+1)}) + \frac{\mathcal{B}^{(k)}}{\beta} \right|\right|_\text{F}^2 . \tag{27} \end{align}$$
这是一个 $\ell_1$ 范数最小化问题，有闭式解。给定一个矩阵 $\mathbf{Y}$，对于一个 $\epsilon > 0$，令 $\mathbf{Z}$ 等于
$$\begin{equation} \mathbf{Z} = \arg\min_{\mathbf{X}} \ || \mathbf{X} ||_1 + \frac{\epsilon}{2} || \mathbf{X} - \mathbf{Y} ||_\text{F}^2 . \end{equation}$$
$\mathbf{Z}$ 可以如下的软阈值操作求解：
$$\begin{equation} \mathbf{Z} = S_{\frac{1}{\epsilon}} (\mathbf{Y}) , \end{equation}$$
其中 $S_{\frac{1}{\epsilon}} (\cdot)$ 是一个按元素的操作，定义为
$$\begin{equation} S_{\frac{1}{\epsilon}} (y_i) = \text{sgn} (y_i) \times \max \left\{ |y_i| - \frac{1}{\epsilon}, 0 \right\} . \end{equation}$$
根据这个定义，可以将 $(27)$ 按如下方式求解
$$\begin{equation} \mathcal{W}^{(k+1)} = S_{\frac{\lambda}{\beta}} \left( \mathcal{F} (\mathcal{X}^{k+1}) - \frac{\mathcal{B}^{(k)}}{\beta} \right) . \tag{28} \end{equation}$$

5 实验结果

略，详见原文。

---

1. Hu Y, Zhang D, Ye J, Li X, He X (2013) Fast and accurate matrix completion via truncated nuclear norm regularization. IEEE Trans Pattern Anal Mach Intell 35(9):2117–2130 ↩
2. Hu Y, Zhang D, Ye J, Li X, He X (2013) Fast and accurate matrix completion via truncated nuclear norm regularization. IEEE Trans Pattern Anal Mach Intell 35(9):2117–2130 ↩
3. Liu J,Musialski P,Wonka P, Ye J (2013b) Tensor completion for estimating missing values in visual data. IEEE Trans Pattern Anal Mach Intell 35(1):208–220 ↩
4. Candès E, Recht B (2009) Exact matrix completion via convex optimization. Found Comput Math 9(6):717–772 ↩
5. Hu Y, Zhang D, Ye J, Li X, He X (2013) Fast and accurate matrix completion via truncated nuclear norm regularization. IEEE Trans Pattern Anal Mach Intell 35(9):2117–2130 ↩
6. Kolda TG, Bader BW (2009) Tensor decompositions and applications. SIAM Rev 51(3):455–500 ↩
7. Yang L, Huang Z-H, Shi X (2013) A fixed point iterative method for low n-rank tensor pursuit. IEEE Trans Signal Process 61(11):2952–2962 ↩
8. Boyd S, Parikh N, Chu E, Peleato B, Eckstein J (2011) Distributed optimization and statistical learning via the alternating direction method of multipliers. Found Trends Mach Learn 3(1):1–122 ↩
9. Merhav N, Kresch R (1998) Approximate convolution using DCT coefficient multipliers. IEEE Trans Circuits Syst Video Technol 8(4):378–385 ↩