本节将用一种表示方法来建立人脸特征检测器,该方法也许是人们认为最简单的模型,即:线性图像模型。由于该算法需表示一个图象块,因此这种面部特征检测器称为块模型( patch model )。该模型在 patch_model 类中被实现,该类的定义和实现可分别在 patch_model.hpp 和 patch_model.cpp 文件中找到。下面的这段代码给出了 patch_model 类的主要功能:
1 class patch_model{ //correlation-based patch expert 2 public: 3 Mat P; //normalised patch 4 ... 5 Mat //response map (CV_32F) 6 calc_response(const Mat &im, //image to compute response from 7 const bool sum2one=false); //normalize response to sum to one? 8 void 9 train(const vector<Mat> &images, //feature centered training images 10 const Size psize, //desired patch size 11 const float var = 1.0, //variance of annotation error 12 const float lambda = 1e-6, //regularization weight 13 const float mu_init = 1e-3, //initial stoch-grad step size 14 const int nsamples = 1000, //number of stoch-grad samples 15 const bool visi = false); //visualize intermediate results? 16 ... 17 };
用于检测面部特征的块模型储存在矩阵 P 中。该类有两个重要的函数: calc_response 和 train 。calc_response 函数会对搜索区域 im 的每个元素计算块模型的响应值。train 函数用来得到块模型,其大小由 psize 参数决定。
本节将介绍该类的功能。下面先来讨论相关性块(correlation patch)和它的训练过程,相关性块将被用来学习块模型。接下来介绍 patch_models 类,该类保存着每个面部特征的一些块模型,并实现了全局变换的功能。在文件 train_patch_model.cpp 和文件 visualize_patch_model.cpp 中的程序分别用来训练可视化块模型。在本节最后将简单介绍它们的用法。
一、 相关性块模型
学习检测器可用两种主要方法:生成方法(generative)和判断方法(discriminative),这两种方法的思想完全不一样。生成方法会学习一个图象块底层表示,这种表示在各种情况下都能最恰当地生成对象外观。而判断方法根据已有的对象的信息来对新对象做出最好的判断,这样已有样本来源于运行的系统。生成方法的优势在于能对具体对象的属性进行操作,可直观地查看新的对象实例的情况。著名的特征脸(Eigenface)方法就是一种流行的生成方法。判别方法的优势在于所建模型直接针对当前问题,通过已有对象来对新对象做出判别。所有判别方法中最著名的也许是支持向量机。将一个图象块作为面部特征来建模时,判别方法更好一些。
1. 学习基于判别法的块模型
给定一个标注数据集,特征检测器可从这些数据学习得到。判别块模型的学习目标是为了构造这样的图象块:当图象块与含有面部特征的图像区域交叉相关时,对特征区域有一个强的响应,而其他部分响应很弱。该模型用数学化方式可表示为:
这里,P 表示块模型,Ii 表示第 i 个训练图像,Ii(a:b,c:d) 表示左上角位置和右下角位置分别是 (a,c) 和 (b,d),圆点表示内积操作,R 表示理想的响应矩阵这个目标函数的解就是一个块模型,此模型会得到一个响应矩阵,该矩阵通常在最小二乘的度量标准下最靠近理想的响应矩阵。对理想的响应矩阵 R 的一个常见选择(假定面部特征集中在训练图象块的中心)是,除中心外其他地方都设为令。在实践中,由于图像被手工标记,总会有标注误差。为了解决这个问题,通常用衰减函数来刻画 R,该函数从中心开始,随着距离增加,函数值会迅速变小。二维高斯分布可很好地充当这种函数,这也相当于假定标注误差服从高斯分布,这个过程可用下图来描述,该图为人脸左眼角外部。
上面给出的目标函数通常称为线性最小二乘,但该问题的自由度,也就是该方法的变量数与块中像素一样多,因此求解最优块模型的计算代价和所需内存会让人望而却步。
解决该问题的有效替代方法是将其看成线性方程,用随机梯度下降法来求解。该算法将团块模型的解集看成一幅地势图,通过不断迭代获得地势图梯度方向的近似估计,并每次向梯度的反方向前进,直到走到目标函数的极小值(达到阈值或迭代次数上限)。另外,随机梯度法每次随机选取样本,只需要很少的样本就能达到最优解,非常适合实时性要求较高的系统。目标函数的梯度公式为:
具体实现步骤:
1. 生成服从高斯分布的理想反馈图像 F
1 // 生成服从高斯分布的理想反馈图像 F 2 Size wsize = images[0].size(); 3 if((wsize.width < psize.width) || (wsize.height < psize.height)){ 4 cerr << "Invalid image size < patch size!" << endl; throw std::exception(); 5 } 6 // 设置反馈图像大小 7 int dx = wsize.width-psize.width,dy = wsize.height-psize.height; 8 Mat F(dy,dx,CV_32F); 9 for(int y = 0; y < dy; y++){ float vy = (dy-1)/2 - y; 10 for(int x = 0; x < dx; x++){ float vx = (dx-1)/2 - x; 11 // 生成函数 12 F.fl(y,x) = exp(-0.5*(vx*vx+vy*vy)/var); 13 } 14 } 15 // 归一化处理 16 normalize(F,F,0,1,NORM_MINMAX);
2. 随机梯度法求最优的团块模型
1 // 利用随机梯度下降法求最优团块模型 2 RNG rn(getTickCount()); 3 // 给定初始更新速率 4 double mu=mu_init,step=pow(1e-8/mu_init,1.0/nsamples); 5 for(int sample = 0; sample < nsamples; sample++){ 6 int i = rn.uniform(0,N); // i 为随机选中的样本图像标记 7 I = this->convert_image(images[i]); dP = 0.0; // 将图像转换为灰度图 8 for(int y = 0; y < dy; y++){ 9 for(int x = 0; x < dx; x++){ 10 Mat Wi = I(Rect(x,y,psize.width,psize.height)).clone(); 11 Wi -= Wi.dot(O); normalize(Wi,Wi); 12 // 计算目标函数的偏导数 D 13 dP += (F.fl(y,x) - P.dot(Wi))*Wi; 14 } 15 } 16 // 更新团块模型 P 17 P += mu*(dP - lambda*P); mu *= step;
学习过程完整代码 train 函数如下:
1 void 2 patch_model:: 3 train(const vector<Mat> &images, // 包含多个样本图像的矩阵向量 4 const Size psize, // 团块模型窗口的大小 5 const float var, // 手工标注错误的方差(生成理想图像时使用) 6 const float lambda, // 调整的参数(调整上一次得到的团块模型的大小) 7 const float mu_init, // 初始步长(构建梯度下降法求团块模型时的更新速率) 8 const int nsamples, // 随机选取的样本数量(梯度下降算法迭代的次数) 9 const bool visi) // 训练过程是否可观察标志 10 { 11 int N = images.size(),n = psize.width*psize.height; 12 13 // 生成服从高斯分布的理想反馈图像 F 14 Size wsize = images[0].size(); 15 if((wsize.width < psize.width) || (wsize.height < psize.height)){ 16 cerr << "Invalid image size < patch size!" << endl; throw std::exception(); 17 } 18 // 设置反馈图像大小 19 int dx = wsize.width-psize.width,dy = wsize.height-psize.height; 20 Mat F(dy,dx,CV_32F); 21 for(int y = 0; y < dy; y++){ float vy = (dy-1)/2 - y; 22 for(int x = 0; x < dx; x++){ float vx = (dx-1)/2 - x; 23 // 生成函数 24 F.fl(y,x) = exp(-0.5*(vx*vx+vy*vy)/var); 25 } 26 } 27 // 归一化处理 28 normalize(F,F,0,1,NORM_MINMAX); 29 30 //allocate memory 31 Mat I(wsize.height,wsize.width,CV_32F); // 被选中的样本灰度图像 32 Mat dP(psize.height,psize.width,CV_32F); // 目标函数的偏导数,大小同团块模型 33 Mat O = Mat::ones(psize.height,psize.width,CV_32F)/n; // 生成团块模型的归一化模板 34 P = Mat::zeros(psize.height,psize.width,CV_32F); // 团块模型 35 36 // 利用随机梯度下降法求最优团块模型 37 RNG rn(getTickCount()); 38 // 给定初始更新速率 39 double mu=mu_init,step=pow(1e-8/mu_init,1.0/nsamples); 40 for(int sample = 0; sample < nsamples; sample++){ 41 int i = rn.uniform(0,N); // i 为随机选中的样本图像标记 42 I = this->convert_image(images[i]); dP = 0.0; // 将图像转换为灰度图 43 for(int y = 0; y < dy; y++){ 44 for(int x = 0; x < dx; x++){ 45 Mat Wi = I(Rect(x,y,psize.width,psize.height)).clone(); 46 Wi -= Wi.dot(O); normalize(Wi,Wi); 47 // 计算目标函数的偏导数 D 48 dP += (F.fl(y,x) - P.dot(Wi))*Wi; 49 } 50 } 51 // 更新团块模型 P 52 P += mu*(dP - lambda*P); mu *= step; 53 if(visi){ 54 Mat R; 55 matchTemplate(I,P,R,CV_TM_CCOEFF_NORMED); // 在样本图像上寻找与团块模型匹配的区域 56 Mat PP; normalize(P,PP,0,1,NORM_MINMAX); 57 normalize(dP,dP,0,1,NORM_MINMAX); 58 normalize(R,R,0,1,NORM_MINMAX); 59 imshow("P",PP); // 归一化的团块模型 60 imshow("dP",dP); // 归一化的目标函数偏导数 61 imshow("R",R); // 与团块模型匹配的区域 62 if(waitKey(10) == 27)break; 63 } 64 }return; 65 }