CN106056039A - 一种针对视线估计的鲁棒混合回归方法 - Google Patents

Abstract

本发明公开了一种针对视线估计的鲁棒混合回归方法，属于计算机视觉技术领域，涉及鲁棒混合回归方法，主要解决视线估计和视线跟踪等视觉映射问题，可以应用于汽车安全驾驶和感兴趣区域检测等领域。首先对采集到的眼部图像进行大小归一化；之后，在输入眼部图像梯度方向直方图特征和对应视线之间建立混合回归模型；接着，初始化混合回归模型的聚类中心，根据初始化聚类中心划分分组，并求解回归参数；最后，对待估计眼部图像提取梯度方向直方图特征，并利用已学到的鲁棒混合回归模型估计头部姿态。

Description

一种针对视线估计的鲁棒混合回归方法

技术领域

本发明属于计算机视觉技术领域，涉及鲁棒混合回归方法，主要解决视线估计和视线跟踪等视觉映射问题，可以应用于汽车安全驾驶和感兴趣区域检测等领域。

背景技术

在计算机视觉中，视线估计是指根据输入的面部图像，定位人眼区域并根据眼球所在位置自动估计视线方向。现有的视线估计方法包括两类：(1)基于几何的方法和(2)基于回归的方法。详见参考文献：Takahiro Ishikawa,Simon Baker,Iain Matthews,andTakeo Kanade，Passive Driver Gaze Tracking with Active Appearance Models,tech.report CMU-RI-TR-04-08,2004。

基于几何的视线估计方法主要通过定位瞳孔中心、上下眼睑特征点以及眼角位置等来计算视线方向。定位上述眼部区域的特征点主要利用活动外观模型(ActiveAppearance Model)，这是一种通过人脸全局外观定位脸部特征点(眼角点、嘴角点等)的方法。基于几何的视线估计方法的优点在于当眼部特征点定位准确时视线估计准确度高，而缺点在于主动外观模型易受光照、遮挡和姿态的影响，造成眼部特征点定位不准确。详见参考文献：Takahiro Ishikawa,Simon Baker,Iain Matthews,and Takeo Kanade，PassiveDriver Gaze Tracking with Active Appearance Models,tech.report CMU-RI-TR-04-08,2004和Iain Matthews,and Simon Baker,Active Appearance Models Revisited,International Journal of Computer Vision,Vol.60,No.2,pp.135-164,2004.

基于回归的视线估计方法主要通过检测人眼区域，并建立眼部图像特征和视线之间的映射关系。现有方法主要通过支持向量回归、高斯过程回归等建立映射关系。基于回归的视线估计方法的主要优点在于当人眼区域定位准确的前提下简单易实现，而缺点是现有回归方法很难非常准确的描述人眼特征与视线之间的映射关系。详见参考文献：ZhiweiZhu,Qiang Ji,and Kristin P.Bennett,Nonlinear Eye Gaze Mapping FunctionEstimation via Support Vector Regression,The18th International Conference onPattern Recognition,Vol.1,pp.1132-1135,2006.和Oliver Williams,Andrew Blake,and Roberto Cipolla.Sparse and Semi-supervised Visual Mapping with theS3GP.The 2006IEEE Computer Society Conference on Computer Vision and PatternRecognition,Vol.1,pp.230-237,2006.

在基于回归方法的视线估计方法中，最重要的问题是建立从眼部特征到视线的非线性回归模型。在非线性回归模型中，混合回归模型已被证实为解决该问题最好的模型之一。由于其高准确性、低复杂度和强通用性，混合回归模型也是目前被最为广泛应用的模型。近年来，国外很多学者将研究集中在如何克服以上缺陷，提高现有基于混合回归视觉映射模型的准确性，尤其是在有干扰和噪声存在情况下的鲁棒性。然而，现有研究工作大都存在建模和解法上的不足，还需要进行深入的研究和完善。此外，国内在该领域也开始有少量研究工作出现，但多数研究还处于起步阶段。混合回归模型最早由多伦多大学的GeoffreyHinton和加州大学伯克利分校的Micheal Jordan在1991年提出，原始模型假设输出变量服从混合高斯分布，并且认为属于不同高斯的输入和输出数据之间存在局部线性映射关系。近几年，INRIA研究机构的Bill Triggs成功地将混合回归模型应用人体姿态估计问题，并假设输入特征和输出变量联合服从混合高斯分布，从而建立从输入图像特征到人体姿态的概率估计模型。详见参考文献：A.Agarwal and B.Triggs,Monocular Human MotionCapture with a Mixture of Regressors,in IEEE Workshop on Vision for HumanComputer Interaction in Conjunction with IEEE Conference on Computer Visionand Pattern Recognition,pp.72-80,2005.

发明内容

本发明提供了一种基于鲁棒混合回归的视线估计方法。首先对采集到的眼部图像进行大小归一化；之后，在输入眼部图像梯度方向直方图特征和对应视线之间建立混合回归模型；接着，初始化混合回归模型的聚类中心，根据初始化聚类中心划分分组，并求解回归参数；最后，对待估计眼部图像提取梯度方向直方图特征，并利用已学到的鲁棒混合回归模型估计头部姿态。算法示意图参见图3。

为了方便地描述本发明内容，首先对一些术语进行定义。

定义1：视线方向。在二维空间中眼球注视方向的角度通常由一个向量表示，该向量由二个元素构成，第一个元素为水平角度，第二个元素为垂直角度。

定义2：梯度方向直方图特征。利用像素强度梯度或边缘的方向分布描述一幅图像中的物体的表象和形状的视觉特征提取方法。其实现方法先将图像分成小的叫做方格单元的连通区域；然后采集方格单元中各像素点的梯度方向或边缘方向直方图；最后把这些直方图组合起来就可以构成特征描述子。为了提高精确度，还可以把这些局部直方图在图像的更大的区间(block)中进行对比度归一化(contrast-normalized)，此方法通过先计算各直方图在这个区间(block)中的密度，然后根据这个密度值对区间中的各个方格单元做归一化。通过该归一化能对光照变化和阴影有更强的鲁棒性。

定义3：高斯分布。又名正态分布，是一个在数学、物理及工程等领域都非常重要的概率分布。具体是指：对于任意随机变量x，若其概率密度函数满足其中μ为正态分布的中心，σ为正态分布的方差，则称其满足高斯分布，常用符号表示。若随机变量为向量x，则对应高维高斯分布μ为高维正态分布的中心，Σ为高维正态分布的协方差，或记为

定义4：混合高斯分布。是指随机变量服从由多个高斯函数混合组成的分布函数，即：其中K为高斯分布的个数，z为指示变量x属于哪一个高斯分布的隐变量。z为一个K维的状态向量，该向量只有一个元素为1，其余元素为0，当第k个元素z_k为1时，则表示x属于第k个高斯分布。p(z)表示隐变量z的出现概率，该概率通常用符号π_k表示，即：p(z_k＝1)＝π_k。

定义5：隐变量。无法被观测到的变量，例如在混合高斯分布 中，z指示变量x属于哪一个高斯分布，无法被观测即为隐变量。

定义6：混合系数。在混合高斯分布中代表变量x属于某一高斯分布的概率，即π_k为混合系数。

定义7：线性回归。是指输出变量(向量)y和输入变量(向量)x之间存在线性映射关系，即：y＝A^Tx+b，其中A^T代表回归矩阵，b代表偏差向量。

定义8：混合回归关系。是指输出变量(向量)y和输入变量(向量)x之间存在非线性映射关系，即：其中代表第k个回归矩阵，b_k代表第k个偏差向量。π_k(x)表示x属于某一线性回归模型的概率。

定义9：似然函数。似然函数是一种关于统计模型中的参数的函数，表示模型参数的似然性或可能性。

定义10：期望最大值算法。一种求解混合高斯分布参数似然函数最大值的优化方法，它通常包含两个基本的步骤，即：期望步骤(E步)和最大化步骤(M步)。

定义11：Q函数。在求解混合高斯分布似然函数最大值的期望最大值算法中，Q函数为似然函数的下界函数。

定义12：对角矩阵。是一个主对角线之外的元素皆为0的矩阵。对角线上的元素可以为0或其他值。diag{}是对角矩阵的符号，花括号中的元素代表对角矩阵对角线上的元素。

定义13：流形。流形(manifold)是局部具有欧几里得空间性质的空间，是欧几里得空间中的曲线、曲面等概念的推广。欧几里得空间就是最简单的流形的实例。地球表面这样的球面则是一个稍微复杂的例子。一般的流形可以通过把许多平直的片折弯并粘连而成。

本发明详细技术方案如下：

步骤1：采集N幅包含不同视线的眼部图像(见图1)，并记录采集每幅图像时对应的视线方向y_n的第一维表示水平方向，第二维表示垂直方向，下标n表示第n幅图像对应的视线方向；

步骤2：将步骤1中采集的眼部图像归一化为大小为100×60像素，并提取梯度方向直方图特征(Histogram of Oriented Gradient,HOG)。在梯度方向直方图特征计算的过程中，区域个数的参数设置为2×2，每一区域中图像单元的个数参数设置为8×8，方向柱的个数设置为9，最后得到任意一幅图像对应的梯度方向直方图特征的维数为4086，并记任意第n幅图像对应的梯度方向直方图特征向量

步骤3：将所有N幅图像对应的梯度方向直方图特征向量按顺序排列可以得到输入数据矩阵X，即X＝[x₁,x₂,…,x_N]。

步骤4：将N幅图像对应的视线方向向量按顺序排列为数据矩阵Y，即Y＝[y₁,y₂,…,y_N]；

步骤5：假设在视线估计问题中，所有输入特征和输出视线方向联合服从混合高斯分布，因此任一输入输出观测对{x_n,y_n}的概率分布为：

各高斯分布中心、协方差和混合系数，分别用符号μ_k，Σ_k和π_k表示，其中下角标k表示第k个高斯分布，k＝1,2,…K。为隐变量，决定第n个观测样本{x_n,y_n}属于哪一个高斯，Θ表示此混合回归模型中的所有参数。同时，假设已知{x_n,y_n}属于某一高斯的条件下，输入与输出之间存在线性回归关系，回归参数为参数矩阵A_k，偏差向量b_k和噪声方差水平γ_k，即：

其中为单位矩阵。根据上述假设，该关系的示意图见图2。最后令Θ＝{μ_k,∑_k,π_k,A_k,b_k,γ_k}代表参数集合。另外，由上述假设可以推断输入梯度方向直方图特征也分别服从混合高斯分布，即：

其中表示输入特征的各高斯分布中心和协方差。由输入和输出之间的线性关系，可知：

${\mathrm{\μ}}_k=\left(\begin{array}{c}{\mathrm{\μ}}_k^x\\ A_k^T{\mathrm{\μ}}_k^x+b_k\end{array}\right),{\mathrm{\Σ}}_k=\left(\begin{array}{cc}{\mathrm{\Σ}}_k^x& {\mathrm{\Σ}}_k^x{A}_k\\ A_k^T{\mathrm{\Σ}}_k^x& A_k^T{\mathrm{\Σ}}_k^x{A}_k+{\mathrm{\γ}}_k{I}_2\end{array}\right)$

步骤6：为了提高混合回归模型的鲁棒性，假设上述聚类可以归为不同的L个分组,即这里代表属于第l个组的聚类标号集合。并假设同组聚类对应的回归参数A_k有相同的先验假设，即：

其中表示高斯先验的均值，τ_lI_2D代表协方差矩阵。这样假设的依据是同一组的聚类为数据流形分布上邻近的聚类，具有相似的输入输出映射关系。同时，我们假设分组在优化的过程中学习得到，邻近的聚类归为一个组。

步骤7：根据步骤5和步骤6中得到鲁棒混合回归模型的所有参数Θ，给定输入数据矩阵X和输出数据矩阵Y时的联合后验概率分布为：

ln p(Θ|x,Y)∝ln p(x,Y|Θ)+ln p(Θ)

混合组回归模型对应的图模型见图3。将参数步骤6中的先验分布p(Θ)和步骤5中的似然函数ln p(x,Y|Θ)的具体表达式代入到上式中，可以得到参数Θ的后验分布为：

步骤8：为了求解步骤8中的参数Θ，我们采用EM算法进行参数优化。EM算法包括E步骤和M步骤。

E步骤：在E步中利用上一步的参数估计结果Θ^old计算任意第n个样本属于第k个聚类的后验概率

符号w_kn代表此后验概率。同时建立关于参数Θ的Q函数，

$\begin{array}{c}Q(\mathrm{\Θ},{\mathrm{\Θ}}^{old})=\underset{Z}{\mathrm{\Σ}}p\left(Z\right|X,Y,{\mathrm{\Θ}}^{old})\ln p\left(Y\right|X,Z,\mathrm{\Θ})+\ln p\left(\mathrm{\Θ}\right)\\ +\underset{Z}{\mathrm{\Σ}}p\left(Z\right|X,Y,{\mathrm{\Θ}}^{old})\ln p(X,Z|\mathrm{\Θ})\end{array}$

其中为隐含变量矩阵，其第k行n列的元素为z_kn。

M步骤：最大化Q函数Q(Θ,Θ^eld)来更新参数Θ，根据上式的第一行得到关于回归参数的优化目标函数

λ>0为回归噪声方差水平γ_k和参数τ_l之间的比值，这里假设同一组的子回归器有相同的回归噪声水平γ_k。上述目标式可以进一步转化为：

$\underset{\tilde{A},G:G^{T}G=I_L}{min}\underset{n=1}{\overset{N}{\Σ}}\underset{k=1}{\overset{K}{\Σ}}{w}_{kn}\left|\right|{\hat{y}}_{nk}-(I_2\⊗x_n^T){\tilde{a}}_k|{|}_2^2+\mathrm{\λ}tr\left(\tilde{A}\right(I_K-{\mathrm{GG}}^T\left){\tilde{A}}^T\right)$

这里为单位矩阵，为所有回归参数构成的参数矩阵，为vec(A_k)， 为偏差向量b_k前一次的取值，分组矩阵的第k行l列元素g_kl指示第k个子回归器是否属于第l个分组，表示否则g_kl＝0。在上述目标函数中，将对分组矩阵G的约束松弛为G:G^TG＝I_L。令M＝GG^T，并固定得到求解M的目标函数：

其中符号代表对称的半正定矩阵集合，通过特征值分解可以得到M的解。

根据上面求得的M＝GG^T的值，得到{A_k}和{b_k}的解分别为：

$vec\left(\tilde{A}\right)={\[\underset{n=1}{\overset{N}{\Σ}}{X}_n^T{W}_n{X}_n+\mathrm{\λ}(I_{2KD}-{\mathrm{GG}}^T\⊗I_{2D})\]}^{-1}\underset{n=1}{\overset{N}{\Σ}}{X}_n^T{W}_n{\hat{y}}_n$

$b_k=\frac{1}{N_k}\underset{n=1}{\overset{N}{\Σ}}{w}_{kn}(y_n-A_k^T{x}_n)$

其中和为单位矩阵，由训练样本输入特征构成的矩阵而

根据Q函数第二行可以得到关于聚类参数{π_k}的目标函数

最后得到解为：

${\mathrm{\μ}}_k^x=\frac{1}{N_k}\underset{n=1}{\overset{N}{\Σ}}{w}_{kn}{x}_n$

${\mathrm{\Σ}}_k^x=\frac{1}{N_k}\underset{n=1}{\overset{N}{\Σ}}{w}_{kn}(x_n-{\mathrm{\μ}}_k^x){(x_n-{\mathrm{\μ}}_k^x)}^T$

${\mathrm{\π}}_k=\frac{N_k}{N}$

步骤9：一直重复步骤9中的E步和M步，直到收敛，即：前后两次参数的值不再改变。

附图说明

图1是视线估计示意图。

图2是针对视线估计的鲁棒混合回归模型。

图3是混合组回归的概率图模型。

具体实施方式

实现语言：Matlab,C/C++

硬件平台：Intel core2E7400+4G DDR RAM

软件平台：Matlab2012a,VisualStdio2010

根据本发明的方法，首先采集一定数量的包含眼部的面部图像并记录这些图像对应的视线方向，根据本发明专利，利用Matlab或者C语言编写基于鲁棒混合回归的视线估计程序，并在采集数据上训练本发明对应模型的参数；接着在各种需要视线估计的场景安装采集摄像头采集原始图像；接着提取采集到的图像眼部区域的梯度方向直方图特征；根据之前训练的参数可以估计该图像对应的视线方向。本发明的方法，可以用于各种场景下视线估计与跟踪。

Claims (1)

1.一种针对视线估计的鲁棒混合回归方法，包括以下步骤：

步骤1：采集N幅包含不同视线的眼部图像，并记录采集每幅图像时对应的视线方向y_n的第一维表示水平方向，第二维表示垂直方向，下标n表示第n幅图像对应的视线方向；

步骤2：将步骤1中采集的眼部图像归一化为大小为100×60像素，并提取梯度方向直方图特征；在梯度方向直方图特征计算的过程中，区域个数的参数设置为2×2，每一区域中图像单元的个数参数设置为8×8，方向柱的个数设置为9，最后得到任意一幅图像对应的梯度方向直方图特征的维数为4086，并记任意第n幅图像对应的梯度方向直方图特征向量

步骤3：将所有N幅图像对应的梯度方向直方图特征向量按顺序排列可以得到输入数据矩阵X，即X＝[x₁,x₂,…,x_N]；

步骤4：将N幅图像对应的视线方向向量按顺序排列为数据矩阵Y，即Y＝[y₁,y₂,…,y_N]；

步骤5：假设在视线估计问题中，所有输入特征和输出视线方向联合服从混合高斯分布，因此任一输入输出观测对{x_n,y_n}的概率分布为：

各高斯分布中心、协方差和混合系数，分别用符号μ_k，∑_k和π_k表示，其中下角标k表示第k个高斯分布，k＝1,2,…K；为隐变量，决定第n个观测样本{x_n,y_n}属于哪一个高斯，Θ表示此混合回归模型中的所有参数；同时，假设已知{x_n,y_n}属于某一高斯的条件下，输入与输出之间存在线性回归关系，回归参数为参数矩阵A_k，偏差向量b_k和噪声方差水平γ_k，即：

其中为单位矩阵。根据上述假设，最后令Θ＝{μ_k,∑_k,π_k,A_k,b_k,γ_k}代表参数集合。另外，由上述假设可以推断输入梯度方向直方图特征也分别服从混合高斯分布，即：

其中表示输入特征的各高斯分布中心和协方差。由输入和输出之间的线性关系，可知：

${\mathrm{\μ}}_k=\left(\begin{array}{c}{\mathrm{\μ}}_k^x\\ A_k^T{\mathrm{\μ}}_k^x+b_k\end{array}\right),{\Σ}_k=\left(\begin{array}{cc}{\Σ}_k^x& {\Σ}_k^x{A}_k\\ A_k^T{\Σ}_k^x& A_k^T{\Σ}_k^x{A}_k+{\mathrm{\γ}}_k{I}_2\end{array}\right)$

步骤6：为了提高混合回归模型的鲁棒性，假设上述聚类可以归为不同的L个分组,即这里代表属于第l个组的聚类标号集合；并假设同组聚类对应的回归参数A_k有相同的先验假设，即：

其中表示高斯先验的均值，τ_lI_2D代表协方差矩阵；这样假设的依据是同一组的聚类为数据流形分布上邻近的聚类，具有相似的输入输出映射关系；同时，我们假设分组在优化的过程中学习得到，邻近的聚类归为一个组；

步骤7：根据步骤5和步骤6中得到鲁棒混合回归模型的所有参数Θ，给定输入数据矩阵X和输出数据矩阵Y时的联合后验概率分布为：

lnp(Θ|X,Y)∝lnp(X,Y|Θ)+lnp(Θ)

将参数步骤6中的先验分布p(Θ)和步骤5中的似然函数lnp(X,Y|Θ)的具体表达式代入到上式中，可以得到参数Θ的后验分布为：

步骤8：为了求解步骤8中的参数Θ，我们采用EM算法进行参数优化。EM算法包括E步骤和M步骤；

E步骤：在E步中利用上一步的参数估计结果Θ^old计算任意第n个样本属于第k个聚类的后验概率

符号w_kn代表此后验概率。同时建立关于参数Θ的Q函数，

$\begin{array}{c}Q(\mathrm{\Θ},{\mathrm{\Θ}}^{\mathrm{old}})=\underset{Z}{\mathrm{\Σ}}p\left(Z\right|X,Y,{\mathrm{\Θ}}^{\mathrm{old}})\ln p\left(Y\right|X,Z,\mathrm{\Θ})+\ln p\left(\mathrm{\Θ}\right)\\ +\underset{Z}{\mathrm{\Σ}}p\left(Z\right|X,Y,{\mathrm{\Θ}}^{\mathrm{old}})\ln p(X,Z|\mathrm{\Θ})\end{array}$

其中为隐含变量矩阵，其第k行n列的元素为z_kn；

M步骤：最大化Q函数Q(Θ,Θ^old)来更新参数Θ，根据上式的第一行得到关于回归参数的优化目标函数

λ>0为回归噪声方差水平γ_k和参数τ_l之间的比值，这里假设同一组的子回归器有相同的回归噪声水平γ_k；上述目标式可以进一步转化为：

$\underset{\tilde{A},G:G^{T}G=I_L}{\min }\underset{n=1}{\overset{N}{\mathrm{\Σ}}}\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{w}_{\mathrm{kn}}{\left|\right|{\hat{y}}_{\mathrm{nk}}-(I_2\⊗x_n^T){\tilde{a}}_k\left|\right|}_2^2+\mathrm{\λtr}\left(\tilde{A}(I_K-{\mathrm{GG}}^T){\tilde{A}}^T\right)$

这里为单位矩阵，为所有回归参数构成的参数矩阵，为 为偏差向量b_k前一次的取值，分组矩阵的第k行l列元素g_kl指示第k个子回归器是否属于第l个分组，表示否则g_kl＝0；在上述目标函数中，将对分组矩阵G的约束松弛为G:G^TG＝I_L；令M＝GG^T，并固定得到求解M的目标函数：

其中符号代表对称的半正定矩阵集合，通过特征值分解可以得到M的解。

根据上面求得的M＝GG^T的值，得到{A_k}和{b_k}的解分别为：

$vec\left(\tilde{A}\right)={\[\underset{n=1}{\overset{N}{\Σ}}{X}_n^T{W}_n{X}_n+\mathrm{\λ}(I_{2KD}-{\mathrm{GG}}^T\⊗I_{2D})\]}^{-1}\underset{n=1}{\overset{N}{\Σ}}{X}_n^T{W}_n{\hat{y}}_n$

$b_k=\frac{1}{N_k}\underset{n=1}{\overset{N}{\Σ}}{w}_{kn}(y_n-A_k^T{x}_n)$

其中和为单位矩阵，由训练样本输入特征构成的矩阵而

根据Q函数第二行可以得到关于聚类参数的目标函数

最后得到解为：

${\mathrm{\μ}}_k^x=\frac{1}{N_k}\underset{n=1}{\overset{N}{\Σ}}{w}_{kn}{x}_n$

${\Σ}_k^x=\frac{1}{N_k}\underset{n=1}{\overset{N}{\Σ}}{w}_{kn}(x_n-{\mathrm{\μ}}_k^x){(x_n-{\mathrm{\μ}}_k^x)}^T$

${\mathrm{\π}}_k=\frac{N_k}{N}$

步骤9：一直重复步骤9中的E步和M步，直到收敛，即：前后两次参数的值不再改变。