截断核范数的矩阵回归及在人脸识别中的应用

Matrix Regression Based on Truncated Nuclear Norm and Application in Face Recognition

摘要: 采用回归方法进行人脸识别的过程中对误差的度量采用像素层的F范数，此方法需要假设像素之间是相互独立的，然而在具有连续遮挡的情况下该假设是不成立的。事实上，误差图像在空间上是相关的，采用矩阵的截断核范数可以更好地描述图像的结构信息。对误差图像进行分析后，提出了一种截断核范数的回归模型，并采用交替方向乘子算法(ADMM)求解。与现有的其他回归方法相比，本文的方法将误差检测和误差矫正集成到一个回归模型中，对Extend Yale B人脸数据库的实验也证明了该方法在人脸识别上的优越性。

Abstract: In the process of face recognition using regression method, the error is measured by the F-norm. This method needs to assume that the pixels are independent of each other, but this assumption is not true in the case of continuous occlusion. In fact, the error images are spatially related, and the truncated nuclear norm of the matrix can describe the structural information of the image. By analyzing the error image, a regression model based on truncated nuclear norm is proposed and solved by the Alternating Direction Multiplier Algorithm (ADMM). Compared with other existing regression methods, the method of this paper integrates error detection and error correction into a regression model. The experiment on Extend Yale B and AR face database also proves the superiority of this method in face recognition.

1. 引言

从十九世纪末法国学者提出人脸识别的概念 [1]，到二十世纪末人脸识别的重大突破。人脸识别的发展经历了一个世纪。当前的人脸识别技术正处于飞速发展阶段，许多国家都开展了有关人脸识别的研究。国内多所研究机构及企业成立了人脸识别研究小组。随着人脸识别技术在网络安全、物业管理、银行、军队及计算机等诸多领域的应用。这一人工智能技术会在不断地应用中逐渐完善。然而通常的人脸识别系统都是基于正常光照且无遮挡的情况，对于特殊的情况(带遮挡、光照)未必能起到很好的识别效果。

基于线性回归分析的方法是人脸识别领域的主要方法。I. Naseem等人提出了一个线性回归分类器(LRC) [2]。LRC是专门为解决小样本人脸识别问题而设计的，虽然在实验中能够取得良好的结果，但无法解决大样本人脸识别问题。

通常对线性回归模型对系数采用基于L2范数或L1范数的正则化可以避免过度拟合。J. Wright等人 [3] 提出了一个基于稀疏表示的分类(SRC)方法。SRC方法假设表示系数是稀疏的并使用所有训练样本组成的字典来表示测试样本，该系数中的非零值集中在训练样本与测试样本相同的类别中。J. Yang等人对SRC进行了深入的研究，并为其有效性寻求理论的支持 [4]。研究发现是L1范数决定了表示的效果，而基于L0范数的正则化只能实现稀疏性。A. Wagner等人 [5] 在前人的基础上进一步扩展了SRC模型并将人脸对齐和识别统一到共同的框架之中。L. Zhang等人 [6] 分析了SRC的工作原理，认为协作表示策略比基于L1范数的稀疏约束更重要。因此提出了一种基于岭回归的协作表示分类器(CRC)。CRC可以实现与SRC相似的结果，并且能够提高算法的速度。但是，大量的实践证明基于CRC模型的人脸识别方法不提供消除噪声的机制，所以它不是稳健的人脸识别方法。

目前基于回归的人脸识别方法都是将每个像素上的误差逐个地表征，而这样的做法存在两个方面缺陷。首先，这种建模需要假定误差像素是互相独立的，这一假设对于像素的随机噪声是没有问题的，但是遮挡与光照变化的情况下，这种假设就不成立。因此，使用基于像素误差模型来解决带有遮挡的图像分类在理论上是有问题的。

在回归分析的人脸识别方法中，使用一组训练图像来表示一个测试图像，在理想情况下，误差图像应当是一个低秩矩阵。在更一般的情况下，测试图像中可能存在光照的变化和部分的遮挡。而光照变化，特别是局部照明变化(如阴影)通常导致低秩误差图像。现有的回归方法，单独表征每一个误差像素，不能够有效利用这种结构信息。

本文针对带有遮挡和光照的人脸图像提出基于截断核范数的识别方法，对矩阵回归中误差图像采用截断核范数度量，通过ADMM算法求解回归系数，并采用误差矩阵的截断核范数为度量进行图像分类。

2. 基于截断核范数的矩阵回归

首先给出了基于截断核范数的矩阵回归，然后利用了交替方向乘子算法进行模型求解。

2.1. 模型建立

给定一组n个图像矩阵 $A_1,\cdots ,A_n\in R^{p\times q}$ 和一个测试图像矩阵 $B\in R^{p\times q}$，测试图像B可以用矩阵 $A_1,\cdots ,A_n$ 线性表示为

$B=x_1{A}_1+x_2{A}_2+\cdots +x_n{A}_n+E$ (1)

其中， $x_1,x_1,\cdots ,x_n$ 是表示系数， $E$ 表示误差图像矩阵。定义

$A\left(x\right)=x_1{A}_1+x_2{A}_2+\cdots +x_n{A}_n$ (2)

则公式(1)可表示为

$B=A\left(x\right)+E$ (3)

(3)式给出了线性矩阵回归模型的一般形式，误差图像 $E$ 的最优解在许多情况下通常是低秩的，因此直观的想法可以通过秩最小化来计算回归系数

$\underset{x}{\min }\text{rank}\left(A\left(x\right)-B\right)$ (6)

然而，采用秩函数来描述矩阵的低秩结构会使模型求解困难，candés [7] 提出采用矩阵的核范数作为矩阵秩函数的替代。

定义1 给定一个矩阵 $Y\in R^{m\times n}$，截断核范数 ${\Vert Y\Vert }_r$ 被定义为 $\min \left(m,n\right)-r$ 个最小奇异值的和，即 ${\Vert Y\Vert }_r={\displaystyle {\sum }_{i=r+1}^{\min \left(m,n\right)}{\sigma }_i}\left(Y\right)$。

由于截断核范数能够更准确的近似秩函数，上述秩最小化问题可转化为截断核范数最小化问题:

$\underset{x}{\min }{\Vert A\left(x\right)-B\Vert }_r$ (7)

为防止过拟合，加入正则项，得到正则化的截断核范数的矩阵回归模型

$\underset{x}{\min }{\Vert A\left(x\right)-B\Vert }_r+\frac{1}{2}\lambda {\Vert x\Vert }_2^2$ (8)

下面求解这个模型。

2.2. 模型求解

本文采用交替方向乘子算法(ADMM)求解截断核范数的矩阵回归问题 [8] [9] [10]。关于ADMM的更多知识，可以参考文献 [8]。根据A. Hansson等人的研究成果 [10]，(8)式中的模型可以重写为

$\begin{array}{l}\underset{x}{\text{min}}{\Vert E\Vert }_r+\frac{1}{2}\lambda {\Vert x\Vert }_2^2\\ \text{s}\text{.t}\text{.}A\left(x\right)-B=E\end{array}$ (9)

增广拉格朗日函数定义为

$L_{\mu }\left(E,x,Z\right)={\Vert E\Vert }_r+\frac{1}{2}\lambda {\Vert x\Vert }_2^2+Tr\left(Z^{\text{T}}\left(A\left(x\right)-E-B\right)\right)+\frac{\mu }{2}{\Vert A\left(x\right)-E-B\Vert }_F^2$ (10)

$\mu >0$ 是惩罚参数， $Z$ 是拉格朗日乘数法， $\text{Tr}(\cdot )$ 表示迹算子。如果 $\mu =0$，(10)表示标准的拉格朗日函数。

ADMM算法由以下迭代组成

1) 令 $E=E^k$ 且 $Z=Z^k$，更新 $x$ ：

$x^{k+1}=\text{arg}\underset{x}{\min }{L}_{\mu }\left(E,x,Z\right)$ (11)

2) 令 $x=x^{k+1}$ 且 $Z=Z^k$，更新 $E$ ：

$E^{k+1}=\arg \underset{E}{\min }{L}_{\mu }\left(E,x,Z\right)$ (12)

3) 令 $x=x^{k+1}$ 且 $E=E^{k+1}$，更新 $Z$ :

$Z^{k+1}=Z^k+\mu \left(E^{k+1}-E^{k+1}\right)$ (13)

算法的关键是解决上述步骤中的优化问题(11)和(12)。因为

$\begin{array}{l}Tr\left(Z^{\text{T}}\left(A\left(x\right)-E-B\right)\right)+\frac{\mu }{2}{\Vert A\left(x\right)-E-B\Vert }_F^2\\ =\frac{\mu }{2}{\Vert A\left(x\right)-\left(B+E-\frac{1}{\mu }Z\right)\Vert }_F^2-\frac{1}{2\mu }{\Vert Z\Vert }_F^2\end{array}$ (14)

得到

$L_{\mu }=\left(E,x,Z\right)={\Vert E\Vert }_r+\frac{1}{2}\lambda {\Vert x\Vert }_2^2+\frac{\mu }{2}{\Vert A\left(x\right)-\left(B+E-\frac{1}{\mu }Z\right)\Vert }_F^2-\frac{1}{2\mu }{\Vert Z\Vert }_F^2$ (15)

$x$ 的解可以通过下式求：

$x^{k+1}=\arg \underset{x}{\min }{L}_{\mu }\left(E,x,Z\right)=\arg \underset{x}{\min }\left({\Vert A\left(x\right)-\left(B+E-\frac{1}{\mu }Z\right)\Vert }_F^2+\frac{\lambda }{\mu }{\Vert x\Vert }_2^2\right)$ (16)

令 $H=\left\{\text{Vec}\left(A_1\right),\cdots ,\text{Vec}\left(A_n\right)\right\}$，则 $A\left(x\right)={\displaystyle \underset{i=1}{\overset{n}{\sum }}{x}_i}{A}_i$ 可以写成矩阵形式 $Hx$.令

$g=\text{Vec}\left(B+E-\frac{1}{\mu }Z\right)$ (17)

因此(16)等同于

$x^{k+1}=\arg \underset{x}{\min }\left({\Vert Hx-g\Vert }_2^2+\frac{\lambda }{\mu }{\Vert x\Vert }_2^2\right)$ (18)

(18)式为标准的岭回归模型，可以得到精确解为：

$x^{k+1}={\left(H^{\text{T}}H+\frac{\lambda }{\mu }I\right)}^{-1}{H}^{\text{T}}g$ (19)

下面考虑问题(12)的解，由于 ${\Vert E\Vert }_r$ 是矩阵奇异值的非凸函数，直接求解并不容易。由文献 [8]，可知

${\Vert E\Vert }_r={\displaystyle \underset{i=1}{\overset{\min \left(m,n\right)}{\sum }}{\sigma }_i}\left(E\right)-{\displaystyle \underset{i=1}{\overset{r}{\sum }}{\sigma }_i}\left(E\right)={\Vert E\Vert }_{*}-\underset{A{A}^{\text{T}}=I,B{B}^{\text{T}}=I}{\max }Tr\left(AE{B}^{\text{T}}\right)$ (20)

因此，求解优化问题 $\underset{E}{\min }{\Vert E\Vert }_r$ 可转化为：

$\underset{E}{\text{min}}{\Vert E\Vert }_{*}-\underset{A{A}^{\text{T}}=I,B{B}^{\text{T}}=I}{\max }Tr\left(AE{B}^{\text{T}}\right)$ (21)

这里 $A\in R^{r\times m},B\in R^{r\times n}$。

定义2 关于矩阵 $Y\in R^{m\times n}$，秩为r的奇异值分解(SVD)，

$Y=U\sum V^{\text{T}},\sum =\text{diag}\left({\left\{{\sigma }_i\right\}}_{1\le i\le r}\right)$ (22)

定义奇异值阈值算子 $D_{\tau }$ 如下： $D_{\tau }\left(Y\right)=U{D}_{\tau }(\sum )V^{\text{T}}$

$D_{\tau }(\sum )=\text{diag}\left({\left\{{\sigma }_i-\tau \right\}}_{+}\right)$ (23)

定理1 对于每个 $\tau \ge 0$ 和 $X\in R^{m\times n}$，有

$D_{\tau }\left(X\right)=\underset{X}{\arg \min }\frac{1}{2}{\Vert Y-X\Vert }_F^2+\tau {\Vert Y\Vert }_{*}$ (24)

因此，问题(12)可转化为

$E^{\text{k}+1}=\underset{E}{\arg \min }{\Vert E\Vert }_{*}-\text{Tr}\left(A_l{T}_k{B}_l^{\text{T}}\right)+\frac{\mu }{2}{\Vert E-T_k\Vert }_F^2-\text{Tr}\left(Z_k^{\text{T}}\left(E-T_k\right)\right)$ (25)

其中 $T_k$ 为求解 $E$ 而引进的辅助变量。

忽略常量项，可得

$E^{k+1}=\underset{E}{\arg \min }{\Vert E\Vert }_{*}+\frac{\mu }{2}{\Vert E-\left(T_k+\frac{1}{\mu }{Z}_k^{\text{T}}\right)\Vert }_F^2$ (26)

由定理1可得

$E^{k+1}=D_{\frac{1}{\mu }}\left(T_k+\frac{1}{\mu }{Z}_k^{\text{T}}\right)$ (27)

计算 $E^{k+1}$ 通过修正 $T_k$ 和 $Z_k$

设置 $E_1=T_1$ 作为初始化，在第l次迭代中，首先基于 $E_l$ 奇异值分解确定了 $E_l$ 并计算 $A_l$ 和 $B_l$。

这里

$A={\left(u_1,\cdots ,u_r\right)}^{\text{T}},B={\left(v_1,\cdots ,v_r\right)}^{\text{T}}$

$\underset{A{A}^{\text{T}}=I,B{B}^{\text{T}}=I}{\max }\text{Tr}\left(AE{B}^{\text{T}}\right)={\displaystyle \underset{i=1}{\overset{r}{\sum }}{\sigma }_i}(\; E\; )$

$T_k$ 的更新过程为

$T_{k+1}=\underset{T}{\arg \min }L\left(E^{k+1},Z_k,T_k,\mu \right)$ (28)

它是T的二次函数，容易求得

$T_{k+1}=E^{k+1}+\frac{1}{\mu }\left(A_l^{\text{T}}{B}_l-Z_k\right)$ (29)

最后计算 $Z_{k+1}$，通过修正 $E^{k+1}$ 和 $T_{k+1}$

$Z_{k+1}=Z_k+\mu \left(E^{k+1}-T_{k+1}\right)$ (30)

S. Boyd [9] 给出了ADMM算法的最优条件和停止准则。

$\max \left({\Vert x^{k+1}-x^k\Vert }_F,{\Vert E^{k+1}-E^k\Vert }_F\right)\le \epsilon$ (32)

在本算法中，由(19)式可知， ${\left(H^{\text{T}}H+\frac{\lambda }{\mu }I\right)}^{-1}{H}^{\text{T}}$ 在每次迭代中都是固定的，因而可以提前计算并保存。令

$M=\left(H^{\text{T}}H+\frac{\lambda }{\mu }I\right)H^{\text{T}}$ (33)

然后，在每个迭代中更新 $x$ 、 $g$ 和矩 $Mg$。另一方面，更新 $E$ 时，主要计算用于执行矩阵 $Q=T_k+\frac{1}{\mu }{Z}_k^{\text{T}}$ 的奇异值分解。

TNMR算法步骤：

1) 输入图像矩阵 $Y$ 和训练图像 $A_1,A_2,\cdots ,A_n$，模型中的参数 $\lambda ,\mu$，及终止条件 $\epsilon$。计算 $H=\left[\text{Vec}\left(A_1\right),\cdots ,\text{Vec}\left(A_n\right)\right]$， $M=\left(H^{\text{T}}H+\frac{\lambda }{\mu }I\right)H^{\text{T}}$

2) while not converged do

a) 更新 $x^{k+1}$， $g=\text{Vec}\left(B+E-\frac{1}{\mu }{Z}^k\right)$， $x^{k+1}=Mg$

b) 更新 $E^{k+1}$， $E^{k+1}=D_{\frac{1}{\mu }}\left(T_k+\frac{1}{\mu }{Z}^k\right)$

通过对 $E^{k+1}$ 的奇异值分解得到 $A_{k+1}$， $B_{k+1}$， $T_{k+1}=E_{k+1}+\frac{1}{\mu }\left(A_{k+1}^{\text{T}}{B}_{k+1}-Z_k\right)$

c) 更新 $Z_{k+1}$， $Z_{k+1}=Z_k+\mu \left(E^{k+1}-T_{k+1}\right)$

d) 判断收敛条件

end while

输出：最优回归系数向量 $x^{*}$

算法可以在两步迭代策略中进行解释，以获得稳健人脸识别策略。更新 $x$ 实际上是用于确定表示系数，更新 $E$ 实际上是一个误差检测步骤，用于确定真正损坏的部分。因此，TNMR方法提供了一个统一的框架，将误差检测和误差校正集成到一个简单的模型中。本文的算法收敛性证明类似于文献 [11] 中所提供的核范数的ADMM算法收敛性证明。

2.3. 算法的复杂度

给定n个训练样本和图像大小 $p\times q$，令 $m=p\times q$，由矩阵乘法决定的计算复杂度为 $O\left(mn\right)$。更新 $E$ 的计算复杂度是 $O\left(\min \left(p^{2}q,p{q}^2\right)\right)$，主要由 $p\times q$ 矩阵 $Q=T_k+\frac{1}{\mu }{Z}_k^{\text{T}}$ 的奇异值分解决定的。因此，算法总的计算复杂度是 $O\left(k\left(\min \left(mp,mq\right)+mn\right)\right)$，其中k是迭代次数。

3. 分类器设计

以给定的一组图像 $A_1,A_2,\cdots ,A_{\text{n}}$ 作为训练样本，与SRC的策略类似，使用所有的训练样本来计算表示系数。对于测试图像 $B$，使用所有训练样本来表示它，并通过求解下面的TNMR模型得到表示系数向量

$x^{*}=\underset{x}{\arg \min }{\Vert A\left(x\right)-B\Vert }_r+\frac{1}{2}\lambda {\Vert x\Vert }_2^2$ (34)

在 $x$ 的最优解基础上，得到重构图像 $\hat{B}$ 和误差图像 $E=B-\hat{B}$。定义

$e_i\left(B\right)={\Vert B-{\hat{B}}_i\Vert }_r$ (35)

${\hat{B}}_i$ 表示第i类训练样本的重构图像。

分类规则定义为：

$e_l\left(B\right)=\underset{i}{\min }{e}_i\left(B\right)$，

则 $B$ 被判定为第l类。

4. 实验

本文采用Extend Yale B数据库 [12] 和AR人脸数据库 [13] 以验证光照变化和带遮挡情况下的人脸识别效果。TNMR方法将与下列方法比较：线性表示分类器(LRC) [2]、协作表示分类器(CRC) [6]、稀疏表示分类器(SRC) [3]、相关熵的稀疏表示(CESR) [14]、稳健稀疏编码(RSC) [15]、HQ_A和HQ_M [16]、核范数矩阵回归(NMR)。除了矩阵回归方法，其他向量回归的方法都采用欧氏距离为分类准则。另外，所有实验都是在原始图像进行。关于实验参数，CRC、SRC都是选择最优的参数使识别率最高。其他方法的参数设置参考原文作者的建议。

4.1. 带光照变化的人脸识别实验

Extend Yale B人脸数据库由耶鲁大学计算机视觉与控制中心创建，包含38人的9个姿势和64种光照条件下的人脸图像，在相机帧率为30帧/秒的情况下，获得的一个人的64张正面人脸图像，所以64张图片的头部姿势和面部表情只有微小的变化。本文的实验采用图片的大小为96 × 84像素。某人的部分图像如图1所示。

Figure 1. Face images under different conditions in Extend Yale B

图1. Extend Yale B中不同光照条件下的人脸图像

本文将分别使用子集1和子集2作为训练集，对应的子集4和子集5 (极端光照条件下的人脸图像)作为测试(如图2，图3所示)，得出的各个方法的识别率如下表1和表2。

Table 1. Subset 2 as the training group and subset 4 as the test group (%)

表1. 子集2作为训练组子集4作为测试组得出的各方法识别率(%)

Table 2. Subset 1 as the training group and subset5 as the test group (%)

表2. 子集1作为训练组子集5作为测试组得出的各方法识别率(%)

Figure 2. Partial face images of subset 2 and subset 4

图2. 子集2和子集4的部分人脸图像

Figure 3. Partial face images of subset 1 and subset 5

图3. 子集1和子集5的部分人脸图像

从表1实验结果来看，除了CESR，SSEC，HQ_A算法外，提出的大部分算法都得到了较好的结果。而本文所提出的TNMR算法达到了最佳效果。一些稳健的稀疏表示方法，如CESR、HQ_A、HQ_M，在极端的光照下不太可靠。SSEC作为一种专为连续遮挡而设计的方法，也不适合极端的光照条件。传统的线性回归方法LRC对光照变化的敏感性较稳健稀疏表示方法要小。而本文TNMR算法得出的结果远优于其他算法结果，并再次达到了最佳效果。上述的实验结果证明了所提出的方法(TNMR)在正常情况下与极端光照条件下的优越性能。

4.2. 带遮挡的人脸识别实验

AR人脸数据库由巴塞罗那计算机视觉中心创建，该数据库包含126人在不同光照、表情下的4000多张人脸图像。图像在两个不同的时间段采集，由于人脸图像具有光照、表情和遮挡变化，因此该数据库是公认的人脸识别数据库。

为了验证本文方法对遮挡的人脸识别的有效性，采用AR数据中两个时间段所拍摄的不同光照和表情的图像作为训练数据(如图4所示)，分别用戴围巾和墨镜的图像作为测试数据(如图5所示)。

Figure 4. Face image under different illumination and expression in AR

图4. AR中不同光照和表情下的人脸图像

Figure 5. Face images with different occlusion in AR

图5. AR中带有不同遮挡的人脸图像

从实验结果可以看出，戴围巾遮挡的各方法识别率比较低(如表3所示)，带有墨镜遮挡的测试数据对识别的影响不是很大(如表4所示)，本文的方法相对于其他方法也没有显著性的提高。而对于带有墨镜遮挡的测试数据，LRC、CESR方法就表现的很不好，只有NMR方法能达到和我们接近的结果。所以，TNMR方法对于带有“污染”的数据仍就表现的很稳健。

Table 3. Recognition rate of each method of wearing scarf to cover (%)

表3. 戴围巾遮挡的各方法识别率(%)

Table 4. Recognition rate of wearing sunglasses (%)

表4. 戴墨镜的各方法识别率(%)

5. 结束语

本文提出了一种基于截断核范数的矩阵回归模型，并利用交替方向乘子算法计算回归系数，提出的TNMR分类器在Extend Yale B人脸图像数据库中进行的研究结果表明：

1) 与常规的基于回归的方法对光照变化和遮挡情况下的人脸识别技术相比TNMR更加稳健。

2) TNMR比结构化稀疏有着更强大的误差编码模型，并且结果对光照的变化不敏感。

3) TNMR可以在极端光照甚至部分遮挡的情况下获得重构后的满意的一般人脸图像。

尽管TNMR比许多稳健的回归方法要好，并且在实验中获得了良好的效果。但本文的模型将截断核范数作为捕获低秩结构噪声的度量，该模型是否对更复杂的噪声有效，或者如何扩展一般噪声模型需要进一步的研究。