2.4 图像特征提取
2.4.1 特征提取算法
(1)特征点
如何选取合适的图像特征一直是图像检测领域的研究热点,不好的图像特征不仅会引入非线性关系,还会引入噪声,影响缺陷判别的实现。对于特征提取算法来说,算法的鲁棒性会直接影响整个缺陷检测系统的鲁棒性,甚至会影响系统的稳定性。目前图像检测与识别领域常用的图像特征有Hu不变矩、Haar特征算子和surf特征,这三种特征都具有平移、旋转、缩放不变性,有良好的鲁棒性。下面介绍前两种特征。
Hu不变矩是由Hu提出的,其利用中心矩构造出7个不变量,可以对区域形状进行描述,而且具有平移、比例、旋转不变性,在机器视觉中是一种十分重要的特征描述方法。
对于分布为f(x,y)的灰度图像,定义其(p+q)阶矩为
且对应的中心矩定义为
式中,,,(x0,y0)为图像质心。
数字图像的分布是不连续的,将矩函数与中心矩函数离散化可得
式中,M和N为图像大小。
定义归一化中心矩为
其中,。
7个不变量分别为
Haar特征是由Papageorgiou提出的,最早用来进行人脸检测与识别的特征算子,随着特征识别算法的普及和传播已经越来越多地应用在手势识别、物体追踪等方面。
首先构造几种不同的特征模板,如图2-28所示。
图2-28 Haar特征模板
a)Haar特征模板1 b)Haar特征模板2 c)Haar特征模板3 d)Haar特征模板4
将模板作为滑动窗口,在目标图像上进行滑动与缩放来求得Haar特征值。对于特征模板1、特征模板2和特征模板3,Haar特征值求解公式为
式中,sum(white)表示模板中白色部分的像素值之和,sum(black)表示模板中黑色部分的像素值之和。
特征模板4的求解方法为
由于Haar特征模板在图像上进行滑动和缩放时产生的特征数量很多,所以之后需要利用主成分分析(PCA)技术筛选和变换特征空间。
为了快速求出如此多的特征,在求解Haar特征时引入了积分图像的概念,公式为
式中,ii(x,y)为待求积分图像;I(i,j)为原始图像在(i,j)像素处的像素值。构建积分图像可以很方便地求任何区域内所有像素值之和。设在原图上有区域A(α,β,γ,ψ),α、β、γ和ψ分别为区域A的4个顶点,则区域A中所有像素和的计算方法为
由式(2-93)和式(2-94)可以看出,求Haar特征的过程就是求出区域像素值和然后作差的过程,这可以映射到积分图像上,使得矩形内的求和运算变成矩形端点的运算,这样,不管矩形尺度如何变换,其计算时间总是一致的,而且只需遍历一次积分图像就可以求得所有子窗口的特征值。
(2)纹理特征
灰度共生矩阵分析方法(GLCM)是建立在图像的二阶组合条件概率密度估计基础上的。它通过计算图像中某一距离和某一方向上两点之间灰度的相关性,来反映图像在方向、间隔、变化快慢及幅度上的综合信息,从而准确描述纹理的不同特性。
灰度共生矩阵是一个联合概率矩阵,它描述了图像中满足一定方向和一定距离的两个灰度值出现的概率,具体定义为:灰度值分别为i和j的一对像素点,位置方向为θ,像素距离为d时的概率,记作p(i,j,d,θ)。通常,对θ=0°,45°,90°,135°,d=1的数字图像而言,其灰度共生矩阵计算公式为
式中,i,j=0,1,…L-1,L是灰度等级,取L=256,是图像中像素的坐标。
由于d、θ选取不同,灰度共生矩阵中向量的意义和范围也不同,因此有必要对p(i,j,d,θ)进行归一化处理
式中,Num为归一化常数,这里取相邻像素对的个数。
为简便起见,后文中忽略了对d和θ的讨论,将归一化后的图像灰度共生矩阵简化为Pij。作为图像纹理分析的特征量,灰度共生矩阵不能直接用于图像特征的分析,而是需要在灰度共生矩阵的基础上计算图像的二阶统计特征参数。Haralick提出了多种基于灰度共生矩阵的纹理特征统计参数,这里采用了常用的7种,并给出了详细的描述。
1)反差,即主对角线的惯性矩:
惯性矩度量灰度共生矩阵的值分布情况和图像中局部变化的大小反映了图像的清晰度和纹理的粗细。粗纹理的Pij值较集中于主对角线附近,i-j较小,所以反差较小,图像较模糊;细纹理反差较大,图像较清晰。
2)熵:
熵度量图像纹理的不规则性。当图像中像素灰度分布非常杂乱、随机时,灰度共生矩阵中的像素值很小,熵值很大;反之,当图像中像素分布井然有序时,熵值很小。
3)逆差距:
逆差矩度量图像纹理局部变化的大小。当图像纹理的不同区域间缺少变化时,其局部灰度非常均匀,图像像素对的灰度差值较小,逆差矩较大。
4)灰度相关:
其中
式中,μ表示均值,σ表示标准差。
灰度相关用来描述矩阵中行或列元素之间的灰度相似度。相关值大表明矩阵中元素均匀相等;反之,表明矩阵中元素相差很大。当图像中相似的纹理区域有某种方向性时,相关值较大。
5)能量(角二阶距):
能量反映图形灰度分布的均匀性和纹理粗细度。当Pij数值分布较集中时,能量较大;当Pij中所有值均相等时,能量较小。如果一幅图像的灰度值均相等,则其灰度共生矩阵Pij只有一个值(等于图像的像素总数),其能量值最大。所以,能量值大就表明图像灰度分布较均匀,图像纹理较规则。
6)集群荫:
7)集群突出:
2.4.2 主成分分析
从原始图像中提取到的特征还不能直接用于分类器的训练:一方面原始特征会带有噪声,直接用于分类器训练容易造成分类器的过拟合;另一方面原始特征往往维度较高而有用数据所占比例较少,直接用于分类器训练会导致训练需要更多的迭代次数。而且有时候由于数据噪声比较多,还会使数据呈非线性分布,这样使用一般的分类器就无法获得很高的判断精度。因此需要对数据进行降维、压缩、去噪及数据空间的变换。目前,常用的数据处理方法有主成分分析(Principal Component Analysis,PCA)、基于核函数的主成分分析(Kernel Principal Component Analysis,KPCA)、线性判别分析(Linear Discriminant Analysis,LDA)。
主成分分析是一种无监督数据特征空间线性转换技术,广泛用于数据降维、去噪。在大数据领域,它是数据挖掘必不可少的一步;在生物工程学领域,它主要应用于基因表达与基因分析;在金融领域,它应用在股票交易与信号去噪中。
主成分分析的主要工作是将原始高维数据映射到新的数据子空间,通过特征之间的相关性可确定数据中存在的模式,还可以利用原始数据中的最大方差方向。通常,子空间的维度不大于原始特征空间。新的数据子空间正交轴即原始高维数据的主成分,也是待求的原始空间最大方差方向。以二维数据空间为例,如图2-29所示。其中,X1、X2为原始数据空间的特征轴;PC1、PC2为变换后新的数据子空间特征轴,也就是原始数据的主成分。
主成分分析的数据压缩过程其实就是求解变换矩阵W的过程,W为k×d的矩阵。其中,d为原始数据空间的维度;k为映射后新的数据子空间维度,且k≤d。数据变换过程为
式中,X为原始数据空间中的数据;Z为变换后新数据子空间中的数据。图2-30为主成分分析的系统流程图。
图2-29 主成分分析的数据子空间示意图
图2-30 主成分分析的系统流程图
数据标准化是机器学习中常用的数据缩放算法,经过特定的缩放过程,使数据变成均值为0、方差为1的分布,从而消除特征范围对主成分分析方向的影响。
令xj为数据X第j维度的特征;μj和σj分别为第j维度数据的平均值与标准差。数据标准化过程如下
式中,为标准化后第j维度的数据。
对于第k维度与第j维度的协方差可以表示为它们与各自期望之差的乘积和后求期望,公式如下
式中,xij表示第i个数据在j特征维度上的特征;xik表示第i个数据在k特征维度上的特征;μj和μk为对应特征维度的平均值。由于数据在计算协方差前进行了标准化,由标准化的定义可知μj和μk均等于0,故公式(2-106)可改写为
由协方差计算公式可知,协方差矩阵是d×d的对称矩阵。协方差矩阵的特征向量代表了主成分(最大方差的方向),对应的特征值决定了特征向量绝对值的大小。因此求解其特征值特征向量,计算方法如下
式中,τ为协方差矩阵;λ为待求特征值;v为特征向量。因为协方差矩阵是实对称矩阵,因此可以采用雅可比迭代法编程实现。
计算出特征值与特征向量后,通过筛选特征值构建变换矩阵W。在降维压缩过程中,需要利用包含原始数据空间中最多信息的特征向量,由主成分分析的定义可知,协方差矩阵的特征值大小反映了特征向量包含信息的多少,因此将特征值与对应的特征向量从大到小排序并提取前k个特征相量构建变换矩阵。同时引入方差解释率
式中,δi为待求第i个特征值的方差解释率;为所有特征值之和。
在实际应用中通常选择k个特征向量,使得满足以下条件
这样在数据压缩的同时可以尽量多地保留原始数据空间中的主成分。最后利用k个特征向量组成的变换矩阵W,根据式(2-104)完成数据转换。
2.4.3 SIFT特征点
SIFT(Scale Invariant Feature Transform,尺度不变特征变换)方法是David Lowe于1999年提出的一种基于尺度空间的图像局部特征表示方法,它具有图像缩放、旋转甚至仿射变换不变的特性,并于2004年得到了更深入的发展和完善。SIFT本质上是一种在不同尺度空间下检测关键点(特征点),并对关键点的方向进行计算的算法。SIFT被广泛应用在机器视觉、三维重建等领域。
一般的SIFT算法分为以下几个步骤完成。
1. 尺度空间的生成
尺度空间理论即采用高斯核理论对初始图片进行尺度变换运算,得到图片在多个尺度下的尺度空间描述序列,最后在尺度空间下对得到的序列进行特征提取。图像尺度的变换通常由高斯卷积核进行唯一确定,不同尺度下的目标图像和高斯卷积核进行卷积的结果就是图像的尺度空间,可以表示为
其中的高斯函数G(x,y,σ)可以实现尺度的变换,其表达式为
式中,σ为高斯尺度因子,随着σ的增大,图像平滑程度慢慢变大,图像变得越模糊;反之,图像保留的细节越丰富,图像变得越清晰。
高斯差分尺度空间公式为
式中,k是两个相邻尺度空间的尺度因子在变换时的倍数,它发生在建立尺度金字塔的过程中。
2. DOG极值点检测与定位
DOG(Difference of Guassian,高斯差分)算子局部极值点就是SIFT算子下图像特征点的子集,在进行极值点检测时,为找到极值点,要将每一个像素点和其三维领域内的26个点进行比较,如果它是这些点中的最大值或最小值,则被保存下来,作为目标图像在这个标准下的特征点,即为候选特征点。
为了消除对比度较低的点和DOG算子产生的不稳定边缘点,可通过拟合三维二次函数来计算出特征点的位置和尺度,进而增强后续图像匹配的稳定性和抗噪能力。将图像的尺度空间函数通过泰勒公式进行展开可得
对其进行求导,选定特征点的精准位置
将式(2-115)代入式(2-114)并取其前两项得
如果,则该特征点被保留下来,否则舍掉。特征点的位置和尺度可以表示为
由于高斯差分后的算子极值点其主曲率在x方向上的值较大,在y方向的数值较小,经过2×2黑塞(Hessian)矩阵的计算得到的结果为
H的最大特征值为α,最小特征值为β,α=γβ,其中γ为比例系数,如果
则去除了边缘相应的较大极值点。
3. 特征点方向分配
上一步中得到了图像的特征点,然后利用特征点邻域像素的梯度方向分布特性对每个特征点添加一个方向,使得这些特征点具有旋转不变性,公式为
式中,m(x,y)和θ(x,y)分别是特征点(x,y)处梯度的模值和方向,L表示各个特征点所在的尺度。
首先以各个特征点为中心创建邻域窗口,然后对创建的邻域窗口进行采样处理,最后将每一个像素梯度方向的次数用直方图来表示,如图2-31所示。
图2-31 梯度方向直方图
至此,特征点检测完毕,每个特征点都可以由3个信息表示:二维位置信息(x,y)、尺度空间信息σ和主方向信息θ。
4. 特征点描述子的生成
为了使图像具有旋转不变性,首先将坐标轴旋转至与特征点主方向一致,然后以特征点为中心,将其周围的8×8邻域窗口(见图2-32)分为4个子块,接着计算每个子块8个方向的梯度方向直方图,最后绘制出每个梯度方向的累加值,从而形成种子点,每个种子点有8个矢量信息。
图2-32 特征点描述子的生成
SIFT特征点的匹配可参考以下程序(https://xmfbit.github.io/2017/01/30/cs131-sift/):
2.4.4 SURF特征点
Bay提出的SURF(Speeded Up Robust Feature,加速稳健特征)算法是一个速度较快、鲁棒性能较好的方法。它是SIFT算法的改进,融合了哈里斯(Harris)特征和积分图像,加快了程序的运行速度。具体来说,该算法可分为建立积分图像、构建黑塞(Hessian)矩阵和高斯金字塔尺度空间、定位特征点、确定主方向、生成特征点描述子等几步完成。
1. 建立积分图像
由于SURF算法的积分图用于加速图像卷积,所以加快了SURF算法的计算速度。对于一个灰度图像I(i,j),图像中的像素积分为
2. 构建黑塞矩阵和高斯金字塔尺度空间
(x,y)为图像中的任意一点,在(x,y)处,尺度为σ的黑塞矩阵H(x,y,σ)可以表示为
式中,Lxx(x,y,σ)是高斯函数与二阶微分在点(x,y)处与图像I(x,y)的卷积,Lxy(x,y,σ)和Lyy(x,y,σ)与此类似。
SURF算法选用DOG算子D(x,y,σ)代替LOG算子来近似地表达,得到类似的黑塞矩阵为
式中,ω=0.9,为矩阵的权重;Dxx、Dyy、Dxy表示箱式滤波和图像卷积的值,取代了Lxx、Lyy、Lxy的值。在进行极值点判断时,如果det(Happrox)的符号为正,则该点为极值点。
图2-33上面的图为先使用高斯平滑滤波,然后在y方向上进行二阶求导的结果;下面的图为滤波后在x和y方向上进行二阶求导的结果。
图2-33 SURF箱式滤波器
3. 定位特征点
得到各像素点的黑塞矩阵后,根据其行列式的正负判断是否为极值点,并使用非极大值抑制法在3×3×3立体邻域检测极值点,只将比它所在尺度层的周围8个点和上下两层对应的9个点都大或者都小的极值点作为候选特征点。
4. 确定主方向
以每个候选特征点为中心,6s(s为特征点尺度)为特征点尺度的半径,以Harr小波统计总响应的60°扇区和x在y方向的所有特征点(Harr小波尺寸4S),高斯分配权重系数的响应,然后以中心角60°扇区模板遍历整个圆形区域,如图2-34所示,将最长的向量作为特征点的方向。
5. 生成特征点描述子
确定主方向后,需要生成特征点描述子。将20s×20s正方形区域中感兴趣的区域分割成4×4的正方形子区域(每个子区域的大小是5s×5s)。如图2-35所示,对于被计算的每一个子区域,将Harr小波响应的水平分量表示为dx,垂直分量表示为dy,然后响应区域dx、dy的和以及绝对值|dx|、|dy|被计算出来,每个子区域形成一个四维的描述向量
图2-34 选取特征点的主方向
这样,最终生成的每一个特征点描述子都是一个4×(4×4)=64维的特征向量,比SIFT算法减小了很多,所以提高了匹配的速度。
图2-35 构造64维的特征点描述子
图2-36所示为SURF特征点匹配结果。特征点的匹配采用了MATLAB工具箱中的matchFeatures函数。
图2-36 SURF特征点匹配
a)参考图 b)匹配图 c)SURF特征点匹配结果