一种开放环境下中医舌象图像分割方法与流程
本发明涉及医学图像处理，尤其是涉及一种开放环境下中医舌象图像分割方法。
背景技术：
传统的中医舌诊，主要依靠医生的肉眼观察进行分析判断，诊断结果往往与医生的知识水平和经验积累有关，近年来，信息技术的发展推动了中医舌象分析客观化、数字化和自动化的进程。国内外学者对此进行了许多有益的探索，并开发了一些舌象分析系统，取得了较好的效果(余兴龙,竺子民,金国潘,等.中医舌诊自动识别系统[j].仪器仪表学报,1994,15(1):67-71；蒋依吾.电脑化中医舌诊系统[j].中国中西医结合杂志,2000,20(2):145-147；caiy.anovelimagingsystemfortongueinspection[c].ieeeimtc,alaska.ieee,2002:159–163；zhanghz,wangkq,zhangd,pangb,huangb.computeraidedtonguediagnosissystem[j].confprocieeeengmedbiosoc,2005,7:6754-6757)，但这些数字化舌诊分析系统的拍摄环境通常是固定的，即在密闭、光照稳定的环境下拍摄采集，这样能够得到高质量的舌象图像，但存在一定局限性，表现在采集设备的特殊性与不易携带。在舌象图像分割方面，香港理工大学david.zhang等人在这一领域取得了较新的研究成果，进行了许多有益的探索(pangb,zhangd,wangk.thebi-ellipticaldeformablecontouranditsapplicationtoautomatedtonguesegmentationinchinesemedicine[j].ieeetransactionsonmedicalimaging,2005,24(8):946；ningj,zhangd,wuc,etal.automatictongueimagesegmentationbasedongradientvectorflowandregionmerging[j].neuralcomputingandapplications,2012,21(8):1-8；cuiz,zuow,zhangh,etal.automatedtonguesegmentationbasedon2dgaborfiltersandfastmarching[m]//intelligencescienceandbigdataengineering.2013:328-335；wuk,zhangd.robusttonguesegmentationbyfusingregion-basedandedge-basedapproaches[j].expertsystemswithapplications,2015,42(21):8027-8038.)。
随着智能手机、平板电脑等移动设备的普及，通过移动设备在开放的自然环境下进行舌象采集，获得个人健康信息逐渐成为一个发展方向。但随之而来的问题是，由于在开放环境下采集舌象，存在光源色温、光线强弱、拍摄角度、设备差异等诸多不确定因素的影响，使得最终获取的图像同固定环境相比，采集的舌象图像往往存在较大差异，对后续的分析带来困难。不可避免的是，在采集舌象的过程中，难免有一些非舌象区域被采集，例如脸、牙齿、嘴唇和背景中的其他物体。然而，舌象特征提取与诊断结果完全依赖于舌体，因此，在开放环境下采集舌象进行分析时，对舌体进行检测与分割显得十分重要。在对舌象进行分析前，对图像进行舌象图像分割有助于提高系统后续的分析速度和准确率，也有助于中医医生直观、准确、方便的观察舌象，提高诊断的速度。舌象图像分割的目的是将舌体图像从图像背景中分离出来，舌象图像分割的准确性直接影响整个系统的容错率与健壮性。
目前，针对开放环境下中医舌象图像分割展开的研究并不多，与该发明方法最接近的技术为中国专利cn203970354u(申请人为厦门强本科技有限公司，发明人王博亮)公开一种基于移动终端的中医舌象分析系统，对移动设备采集的舌象进行处理、分析。该系统虽然实现了舌象图像分割功能，但也存在一些不足之处：
1、采用穷举的方式对待分割目标切分小块，转换到霍夫空间投票，算法复杂，处理速度慢；而在开放环境下采用移动设备进行舌象采集、诊断，往往对处理过程的实时性要求较高；
2、处理过程中需要用较多的参数来描述舌象模型，比较复杂，实现比较困难；
3、需要有监督的训练大量数据，需要人工手动标记出前景和背景，从而建立检测模型，训练过程难度大；
4、进行数据训练时，要求样本尽可能多，对数据集的依赖性高且容易过拟合，陷入局部极值。
此外，中国专利cn106023151a(本申请的发明人黄晓阳,郑丰,王博亮,王彦晖)公开一种开放环境下中医舌象目标检测方法，但是方法采用单一纹理特征进行检测，对舌体边缘的细节部分分割效果不够精确，该方法对舌体分割的准确率和鲁棒性有待提高。
技术实现要素：
本发明的目的在于针对现有的在开放环境下采集舌象过程中存在诸多不足之处，如：光源色温、光线强弱、拍摄角度、设备差异等诸多不确定因素的影响，使得最终获取的图像同固定环境相比，采集的舌象图像往往存在较大差异，对后续的分析带来困难等问题，提供可快速并准确地分割出可供后续分析的舌象，以便于进行后续舌象分析的一种开放环境下中医舌象图像分割方法。
本发明包括以下步骤：
1)输入在开放环境下采集的图像a；
2)对采集的图像a进行颜色校正，得到校正图像b；
3)对步骤2)得到的校正图像b进行图像初步分割；
4)对步骤3)得到的图像c中的各个连通域进行特征提取；
5)对步骤4)提取到的各个凸包的特征，用随机森林分类器计算舌体位置；
6)以步骤5)得到的舌体区域d为基础，对图像b进行开放环境下中医舌象图像分割方法。
在步骤2)中，所述对采集的图像a进行颜色校正的具体方法可为：
(1)在标准光照环境下采集舌象图像s1，计算舌象图像s1的rgb三个颜色通道均值分别与舌象图像s1的整体均值ks的比值αr、αg、αb：
k＝(ravg+gavg+bavg)/3
其中，k为图像a的整体均值，ravg、gavg、bavg分别为图像a的rgb三个颜色通道的均值；rd、gd、bd为校正图像b每个像素点的rgb三个颜色通道的值，rs、gs、bs为图像a每个像素点的rgb三个颜色通道的值；
在步骤3)中，所述对步骤2)得到的校正图像b进行图像初步分割的具体方法可为：
(1)将校正图像b转换到灰度空间图像fb1，按照最大类间方差法(参见文献：ohtsun.athresholdselectionmethodfromgray-levelhistograms[j].systemman&cyberneticsieeetransactionson,1979,9(1):62-66.)对灰度空间图像fb1进行阈值分割，得到分割图像b1’，并对分割图像b1’运用形态学运算平滑连通域，得到图像b1；
(2)将校正图像b转换到hsv颜色空间图像fb2，对图像fb2的h通道进行阈值分割，得到分割图像b2’，并对分割图像b2’运用形态学运算平滑连通域，得到图像b2；
(3)采用rgb三色分量方差法(参见文献：蒋依吾.电脑化中医舌诊系统[j].中国中西医结合杂志,2000,20(2):145-147)对校正图像b进行阈值分割,得到分割图像b3’，并对分割图像b3’运用形态学运算平滑连通域，得到图像b3。
(4)对图像b1、图像b2、图像b3三个图像进行逻辑“与”运算，得到图像c。
在步骤3)第(1)、(2)、(3)部分中，所述将分割结果图像运用形态学运算平滑连通域，其具体步骤为对分割结果图像b(x,y)进行先腐蚀再膨胀：
b(x,y)＝dilate(erode(b(x,y),element))
其中，element定义为形态学运算中的结构元素；dilate定义为形态学运算中的膨胀操作；erode定义为形态学运算中的腐蚀操作。
在步骤4)中，所述对步骤3)得到的图像c中的各个连通域进行特征提取的具体方法可为：
(1)对于图像c中的每一个连通域，计算该连通域的凸包si；
(2)提取凸包区域的空间位置特征；
在步骤4)第(2)部分中，所述提取凸包区域的空间位置特征具体步骤如下：
a)计算凸包si与图像c的面积比
b)计算连通域质心ci离图像中心c0的欧式距离
c)计算连通域最小右边界矩形的长宽比rscale＝w/h，其中w表示最小外接矩形的宽，h表示最小外接矩形的长。
(3)提取凸包区域的颜色特征；
在步骤4)第(3)部分中，所述提取凸包区域的颜色特征具体步骤如下：
a)将校正图像b转换到hsv颜色空间图像c’；
b)凸包si对应于图像c’中的区域，计算区域的非均匀量化颜色特征。
所述非均匀量化颜色特征定义如下：
根据人体的感知特性对h，s，v三个颜色通道进行非均匀量化，将色调h分为8份，饱和度s分为2份，亮度v分为1份。其中色调h的取值范围为[0,360]，饱和度s的取值范围为[0,1]，亮度v的取值范围为[0,1]。
v＝{0v∈(0.15,1]
按照以上的量化等级，将h、s、v三个个颜色分量合成为一维特征矢量：
l＝hqsqv+sqv+v
其中qs和qv分别是分量s和v的量化级数，这里qs＝2，qv＝1，因此上式可表示为：
根据上式，l取值范围为[0,1,…,19]，计算l获得20柄的一维直方图，通过统计并归一化我们能够得到20个颜色特征。
(4)提取凸包区域的形状特征；
在步骤4)第(4)部分中，所述提取凸包区域的形状特征具体为计算凸包si的7个hu不变矩mi(i∈[1,7])(参见文献：hum.visualpatternrecognitionbymomentinvariants[j].informationtheoryiretransactionson,1962,8(2):179-187.)。
所述hu不变矩定义如下：
对于大小为m×n的图像f(x,y)，f(x,y)的二维(p+q)阶矩定义为：
相应的(p+q)阶中心矩定义为：
由ηpq表示的归一化中心矩定义为：
hu用归一化中心矩的线性组合构造具有平移、伸缩、旋转不变的7个hu不变矩
(5)提取凸包区域的纹理特征；
在步骤4)第(5)部分中，所述提取凸包区域的纹理特征具体具体步骤如下：
a)将校正图像b转换到灰度空间图像c”；
b)凸包si对应于图像c”中的区域，分别计算区域的水平、垂直、对角偏移量三种灰度
共生矩阵的联合概率密度分布；
设灰度共生矩阵的联合概率密度分布记为[pmn]l×l，其中l为灰度取值范围，m＝[0,l-1],n＝[0,l-1]；根据联合概率密度分布计算得到6个纹理特征(参见文献：haralickrm,shanmugamk,dinsteini.texturalfeaturesforimageclassification[j].systemsman&cyberneticsieeetransactionson,2010,smc-3(6):610–621.)，所述6个纹理特征包括角二阶矩pcon、对比度pasm、熵pent、逆差矩pidm、中值pmean和灰度相关pcor。
在步骤5)中，所述随机森林分类器的定义为：
随机森林(参见文献：surhonelm,tennoemt,henssonowsf,etal.randomforest[j].machinelearning，2001,45(1):5-32)是一种以决策树作为基础分类器的分类算法。随机森林分类器由多棵决策树构成，这些决策树彼此独立，而且在训练样本的选择、特征子集的选择和树的生长过程中引入了随机性，即每棵决策树的构造都是从原数据集中随机抽取出一部分作为样本子集，然后再从样本子集中随机选取部分特征子集进行处理，并且所有决策树都自然生长，不进行剪枝。因此，多颗决策树可以并行生长，这样随机森林分类器相对于其他分类器具有更短的训练和识别时间，并且由于引入了随机性，分类器对噪声也不敏感。随机森林的最终决策是通过所有决策树投票的方法而得，因此，具有更高的检测率与鲁棒性。
所述随机森林采用bagging方法选取每个训练子集，即采用有放回的抽取方式，这样会导致原数据集中的部分样本可能重复出现在新的训练子集中，而另外一部分样本一次也不出现，这些没有被选择到的数据叫做“outofbag(oob)”数据，一般占所有数据的1/3。
对于模型中的每个节点，决策树的构造通过随机抽取特征子集来分裂生长，随机抽取的数量为特征总数的开方。对于分类问题，决策树生长过程中节点最佳分裂的度量选择基尼指数，对于一个含有m类样本的集合s，集合s的纯度用基尼值度量为：
其中pk为第k类样本占总样本数量的比例，k＝1,2,...,m。集合s的纯度与计算得到的基尼值成反比。对于节点的分裂，假定用属性v对集合s进行分裂，其中属性v为具有n个可能取值的离散属性{v1,v2,...,vn}，则分裂后的集合s会产生n个分支，其中每个分支对应于原始集合s中在属性v上取值相同的所有样本，记为sv。则根据公式可计算出sv的基尼纯度，再对每个分支赋予权重|sv|/|s|，其中|sv|为第v个分支的样本数量，|s|表示所有样本数量。于是，可计算出用属性v对集合s进行分裂所获得的基尼指数：
从上式中可以看出，基尼指数值越小，意味着所有分支节点的加权纯度越小，对集合s纯度的提升越明显。因此，在候选属性集合a中，选择使得分裂后gindex值最小的属性作为最佳分裂属性，即：
训练随机森林分类器的具体方法可为：
a)对于训练集中的所有图像，对每个图像按照步骤2)进行颜色校正然后按照步骤3)对校正图像进行图像分割，最后对分割图像的各个凸包按照步骤4)提取特征并人工标记样本；
b)在舌象训练样本中用bagging方法形成n个子样本；
c)对于每个子样本，随机选择p个属性作为节点分裂的候选属性，其中p取值为所有属性的开方；
d)在候选属性中计算基尼指数，选择吉尼指数最小的属性对决策树进行分裂；
f)重复步骤d)直到吉尼指数小于某一规定阈值；
g)重复步骤c)到f)，直到生成n棵决策树；
h)对于未知样本的分类，随机森林输出其决策树多数投票结果。
所述对步骤4)提取到的各个凸包的特征，用步骤5)第(1)部分训练好的随机森林分类器输出最终舌体位置，其具体步骤如下：
a)对于各个凸包，随机森林分类器输出其为舌体的可能值pi；
b)所有pi值里最大值对应的凸包即为最终舌体区域d。
在步骤6)中，所述以步骤5)得到的舌体区域d为基础，对图像b进行开放环境下中医舌象图像分割方法的具体方法可为：
(1)将图像b转换到灰度空间图像b’；
(2)截取图像b’中的舌象图像分割区域，具体步骤如下：
a)计算舌体区域d对应的最小右边界矩形r；
b)将最小右边界矩形r向外扩散20个像素点，得到舌象图像分割区域r’。
(3)对舌象图像分割区域r’进行边缘增强，具体步骤如下：
a)计算舌象图像分割区域r’水平与垂直的梯度矢量场，分别为u和v；
b)将水平与垂直梯度矢量场叠加到舌象图像分割区域r’，得到边缘增强后的最终舌象图像分割区域r”，即区域r”的最终灰度取值为r’+u+v。
(4)以舌体区域d对应的凸包为初始轮廓，采用gvfsnake模型(参见文献：xuc,princejl.snakes,shapes,andgradientvectorflow[j].ieeetransactionsonimageprocessingapublicationoftheieeesignalprocessingsociety,1998,7(3):359-69.)对最终舌象图像分割区域r”进行分割，得到舌象分割图像。
所述gvfsnake模型定义为：
在传统主动轮廓模型(snake)的基础上，采用梯度矢量流作为曲线演化的外部能量。传统snake模型是在图像区域内定义的一条轮廓曲线，其参数表达形式为x(s)＝[x(s),y(s)],s∈[0,1]，其中s是用傅里叶变换形式描述轮廓边界的变量，x(s)和y(s)表示曲线在图像区域内的坐标。该轮廓曲线的运动形式是通过最小化能量函数：
其中，α和β分别为曲线弹性和刚性的加权系数，相应的，x′(s)和x″(s)分别为曲线x(s)关于s的一阶导数和二阶导数。
要使能量函数e最小化，曲线x(s)必须满足欧拉(euler)方程：
为了求解该方程，引入时间变量t，将x(s)视为s和时间t的函数x(s,t)，则上述动态演化方程可表示为对t求偏微分：
xt(s,t)＝αx″(s,t)-βx″″(s,t)-▽eext
当xt(s,t)稳定时，xt(s,t)＝0，满足欧拉方程，曲线演化完毕。
在上述公式的基础上，gvfsnake模型定义一个新的静态外部能量
xt(s,t)＝αx″(s,t)-βx″″(s,t)+v
对于新的静态外部能量gvf，相应的，将梯度矢量场定义为v(x,y)＝(u(x,y),v(x,y))，通过最小化下述能量方程得到：
其中，下标x和y表示的是图像在x和y方向上的偏导数。参数μ是控制能量函数在第一项和第二项之间权衡的正则化参数，根据图像中的噪声而定，当图像中的噪声较大时，μ取大值。
相应地，最小化能量方程可以转化为求解欧拉方程组：
其中，▽2为拉普拉斯算子。
进一步的，类似于传统snake模型的求解，欧拉方程组的求解可将分量u和v看作是时间t的函数：
当方程组趋于稳定的解即为欧拉方程组的期望解。由于方程组是解耦的，因此可以作为单独的标量偏导数求解方程。为了方便求解，改写方程组为：
其中，
本发明针对开放环境下中医舌象图像的特点对图像进行舌象图像分割，首先对图像进行颜色校正的预处理，减少因外界光源色温带来的影响；然后为了确定初始轮廓，先对图像进行分割，得到多个连通区域；并对各连通区域进行特征提取，最终通过训练好的随机森林分类器得到舌体初始轮廓；最后通过gvfsnake模型对初始轮廓进行演化，逼近舌体真实边缘。本发明最终达到对图像中的舌象进行分割的目的，即将舌体从复杂背景中提取出来。
颜色校正过程采用改进的灰度世界算法，针对舌象图像特有的三个颜色通道均值不同的特征调整相应的参数值，以符合舌象图像的特征；图像分割部分采用最大类间方差、色调阈值分割和rgb三色分量方差分割相结合的方法进行分割；确定初始轮廓部分通过对各连通域的凸包进行特征提取后，采用训练好的随机森林分类器得到舌体区域的轮廓；舌象图像分割则以上一部的舌体轮廓作为初始轮廓，采用gvfsnake模型对曲线进行演化，逼近舌体真实轮廓，最终达到舌象图像分割的目的。
相对于现有技术，本发明的有益效果为：
1、本发明根据舌象特征利用改进的灰度世界算法对图像进行颜色校正，减少了因采集环境颜色偏色带来的影响。
2、本发明采用最大类间方差、色调阈值分割和rgb三色分量方差分割相结合的方法进行图像分割；根据连通域凸包的空间位置特征、颜色特征、形状特征和轮廓特征，采用随机森林分类器进行提取舌体初始轮廓；采用gvfsnake模型进行舌象图像精细分割；这些方法都是针对开放环境下的舌象特点进行处理、分析，这些方法相辅相成，环环相扣，结合起来使用能够更加精确地分割出舌象，具有分割精度高、鲁棒性强的优点。
3、本发明使用的方法运算量小、程序实现简单，随机森林分类器不需要训练大量数据，特征提取算法复杂度低，处理速度快，对数据集的依赖程度不高。
4、本发明识别、分割的准确率高，综合考虑了舌象图像空间位置、颜色、形状和纹理等多种特征，综合考虑图像整体信息，相比其他方法采用单一特征而言，本发明具有更高的识别率和准确率。
本发明针对开放环境下舌象图像的特点分多个步骤进行处理、分析，达到舌象图像分割的目的。本发明的图像分割、特征提取、分类器识别、gvfsnake模型分割几个步骤在功能上相辅相成，图像分割主要利用阈值分割的方法进行粗筛选，算法复杂度低，特征提取利用图像的空间位置、颜色、形状和轮廓特征，算法实现简单，随机森林分类器识别舌体初始轮廓，鲁棒性强，舌象图像分割在初始轮廓的基础上利用主动轮廓模型进行分割，分割的精确度较高，运算速度快，取得了新的技术效果。相对于最接近的现有技术，具有明显的区别特征，且具有较明显的优点：本发明的方法不需要采用穷举的方法转换到霍夫空间投票，算法复杂度低，运算量小，也不需要进行有监督的训练大量数据，因此实现较为简单，对数据集的依赖程度不高，在准确率与处理速度上更优越，更具有实用性。
附图说明
图1为开放环境下采集的图像；
图2为本发明步骤5)得到的舌体初始轮廓；(a线为初始轮廓)
图3为本发明步骤6)舌象图像分割的最终结果
具体实施方式
下面结合附图和实施例对本发明作进一步详细描述。
本发明实施例包括如下步骤：
1)输入在开放环境下采集的图像a，如图1所示；
2)对采集的图像a进行颜色校正，得到校正图像b；
在步骤2)中，所述对采集的图像a进行颜色校正的具体方法如下：
(1)在标准光照环境下采集舌象图像s1，计算舌象图像s1的rgb三个颜色通道均值分别与舌象图像s1的整体均值ks的比值αr、αg、αb：
其中，舌象图像s1的整体均值ks＝(ravgs+gavgs+bavgs)/3；ravgs、gavgs、bavgs分别为标准光照环境下采集舌象图像s1的rgb三个颜色通道的均值；所述标准光照采集环境可为d65光源，色温为6500k；
本实施例中，各个参数的取值如下：αr＝1.09、αg＝0.956、αb＝0.94。
(2)按下式调整图像a的rgb三个颜色通道的均值，得到校正图像b；
k＝(ravg+gavg+bavg)/3
其中，k为图像a的整体均值，ravg、gavg、bavg分别为图像a的rgb三个颜色通道的均值；rd、gd、bd为校正图像b每个像素点的rgb三个颜色通道的值，rs、gs、bs为图像a每个像素点的rgb三个颜色通道的值；
3)对步骤2)得到的校正图像b进行图像分割，具体方法如下：
(1)将校正图像b转换到灰度空间图像fb1，按照最大类间方差法对灰度空间图像fb1进行阈值分割，得到分割图像b1’，并对分割图像b1’运用形态学运算平滑连通域，得到图像b1；
在步骤3)第(1)部分中，所述将校正图像b转换到灰度空间图像fb1，按照最大类间方差法对灰度空间图像fb1进行阈值分割，得到分割图像b1’，并对分割图像b1’运用形态学运算平滑连通域，得到图像b1的具体步骤如下：
a)将校正图像b转换到灰度空间图像fb1；
b)对于灰度空间图像fb1，按照最大类间方差法得到的阈值t对灰度空间图像fb1进行阈值分割，得分割图像b1’像素点的rgb取值fx,y(r,g,b)：
其中，fb1(x,y)表示灰度空间图像fb1像素点的取值；
c)对分割图像b1’运用形态学运算中的闭运算平滑连通域得到图像b1；图像b1中的像素值g1(x,y)为：
g1(x,y)＝dilate(erode(f1(x,y),element))
其中，f1(x,y)为分割图像b1’中的像素值，element定义为形态学运算中的结构元素；dilate定义为形态学运算中的膨胀操作；erode定义为形态学运算中的腐蚀操作。
本实施例中，element定义为[11×11]的椭圆结构
(2)将校正图像b转换到hsv颜色空间图像fb2，对图像fb2的h通道进行阈值分割，得到分割图像b2’，并对分割图像b2’运用形态学运算平滑连通域，得到图像b2；
在步骤3)第(2)部分中，所述将校正图像b转换到hsv颜色空间图像fb2，对图像fb2的h通道进行阈值分割，得到分割图像b2’，并对分割图像b2’运用形态学运算平滑连通域，得到图像b2的具体步骤如下：
a)将校正图像b转换到hsv颜色空间图像fb2；
b)采用下式对图像fb2进行色调阈值分割，得到分割图像b2’像素点的rgb取值fx,y(r,g,b)：
其中，hx,y表示图像fb2中h通道像素点取值，t1和t2表示设定的阈值；
本实施例中，t1＝7、t1＝29
c)对分割图像b2’运用形态学运算中的闭运算平滑连通域得到图像b2；图像b2中的像素值g2(x,y)为：
g2(x,y)＝dilate(erode(f2(x,y),element))
其中，f2(x,y)为分割图像b2’中的像素值，element定义为形态学运算中的结构元素；dilate定义为形态学运算中的膨胀操作；erode定义为形态学运算中的腐蚀操作。
本实施例中，element定义为[11×11]的椭圆结构
(3)采用rgb三色分量方差法对校正图像b进行阈值分割,得到分割图像b3’，并对分割图像b3’运用形态学运算平滑连通域，得到图像b3。
在步骤3)第(3)部分中，所述采用rgb三色分量方差法对校正图像b进行阈值分割,得到分割图像b3’，并对分割图像b3’运用形态学运算平滑连通域，得到图像b3的的具体步骤如下：
a)对于校正图像b，假设其大小为m×n，将校正图像b中每个像素点的rgb取值进行归一化操作，其取值范围为[0,1]，采用下式对校正图像b中的每个像素点计算rgb三色分量方差gate(m,n)，并对校正图像b进行分割，得到分割图像b3’像素点的rgb取值fm,n(r,g,b)：
gate(m,n)＝(rm,n-gm,n)+(bm,n-gm,n)×6+(rm,n+gm,n+bm,n)/3
其中，rm,n表示校正图像b中像素点(m,n)的r通道的取值，gm,n表示校正图像b中像素点(m,n)的g通道的取值，bm,n表示校正图像b中像素点(m,n)的b通道的取值；
b)对分割图像b3’运用形态学运算中的闭运算平滑连通域得到图像b3；图像b3中的像素值g3(x,y)为：
g3(x,y)＝dilate(erode(f3(x,y),element))
其中，f3(x,y)为分割图像b3’中的像素值，element定义为形态学运算中的结构元素；dilate定义为形态学运算中的膨胀操作；erode定义为形态学运算中的腐蚀操作。
本实施例中，element定义为[11×11]的椭圆结构
(4)对图像b1、图像b2、图像b3三个图像进行逻辑“与”运算，得到图像c。
4)对步骤3)得到的图像c中的各个连通域进行特征提取，具体方法如下：
(1)对于图像c中的每一个连通域，计算该连通域的凸包si；
(2)提取凸包区域的空间位置特征；
在步骤4)第(2)部分中，所述提取凸包区域的空间位置特征具体步骤如下：
a)计算凸包si与图像c的面积比
b)计算连通域质心ci离图像中心c0的欧式距离
c)计算连通域最小右边界矩形的长宽比rscale＝w/h，其中w表示最小外接矩形的宽，h表示最小外接矩形的长。
本实施例中，共获得3维空间位置特征
(3)提取凸包区域的颜色特征；
在步骤4)第(3)部分中，所述提取凸包区域的颜色特征具体步骤如下：
a)将校正图像b转换到hsv颜色空间图像c’；
b)凸包si对应于图像c’中的区域，计算区域的非均匀量化颜色特征。
所述非均匀量化颜色特征定义如下：
根据人体的感知特性对h，s，v三个颜色通道进行非均匀量化，将色调h分为8份，饱和度s分为2份，亮度v分为1份。其中色调h的取值范围为[0,360]，饱和度s的取值范围为[0,1]，亮度v的取值范围为[0,1]。
v＝{0v∈(0.15,1]
按照以上的量化等级，将h、s、v三个个颜色分量合成为一维特征矢量：
l＝hqsqv+sqv+v
其中qs和qv分别是分量s和v的量化级数，这里qs＝2，qv＝1，因此上式可表示为：
根据上式，l取值范围为[0,1,…,19]，计算l获得20柄的一维直方图，通过统计并归一化得到颜色特征。
本实施例中，共获得20维颜色特征
(4)提取凸包区域的形状特征；
在步骤4)第(4)部分中，所述提取凸包区域的形状特征具体为计算凸包si的7个hu不变矩mi(i∈[1,7])
所述hu不变矩定义如下：
对于大小为m×n的图像f(x,y)，f(x,y)的二维(p+q)阶矩定义为：
相应的(p+q)阶中心矩定义为：
由ηpq表示的归一化中心矩定义为：
hu用归一化中心矩的线性组合构造具有平移、伸缩、旋转不变的7个hu不变矩
本实施例中，共获得7维形状特征
(5)提取凸包区域的纹理特征；
在步骤4)第(5)部分中，所述提取凸包区域的纹理特征具体具体步骤如下：
a)将校正图像b转换到灰度空间图像c”；
b)凸包si对应于图像c”中的区域，分别计算区域的水平、垂直、对角偏移量三种灰度
共生矩阵的联合概率密度分布；
设灰度共生矩阵的联合概率密度分布记为[pmn]l×l，其中l为灰度取值范围，m＝[0,l-1],n＝[0,l-1]；根据联合概率密度分布计算得到6个纹理特征，所述6个纹理特征包括角二阶矩pcon、对比度pasm、熵pent、逆差矩pidm、中值pmean和灰度相关pcor。
本实施例中，共获得18维纹理特征
5)对步骤4)提取到的各个凸包的特征，用训练好的随机森林分类器计算舌体位置，如图2所示；
所述随机森林分类器定义为：
随机森林是一种以决策树作为基础分类器的分类算法。随机森林分类器由多棵决策树构成，这些决策树彼此独立，而且在训练样本的选择、特征子集的选择和树的生长过程中引入了随机性，即每棵决策树的构造都是从原数据集中随机抽取出一部分作为样本子集，然后再从样本子集中随机选取部分特征子集进行处理，并且所有决策树都自然生长，不进行剪枝。因此，多颗决策树可以并行生长，这样随机森林分类器相对于其他分类器具有更短的训练和识别时间，并且由于引入了随机性，分类器对噪声也不敏感。随机森林的最终决策是通过所有决策树投票的方法而得，因此，具有更高的检测率与鲁棒性。
随机森林采用bagging方法选取每个训练子集，即采用有放回的抽取方式，这样会导致原数据集中的部分样本可能重复出现在新的训练子集中，而另外一部分样本一次也不出现，这些没有被选择到的数据叫做“outofbag(oob)”数据，一般占所有数据的1/3。
对于模型中的每个节点，决策树的构造通过随机抽取特征子集来分裂生长，随机抽取的数量为特征总数的开方。对于分类问题，决策树生长过程中节点最佳分裂的度量选择基尼指数，对于一个含有m类样本的集合s，集合s的纯度用基尼值度量为：
其中pk为第k类样本占总样本数量的比例，k＝1,2,...,m。集合s的纯度与计算得到的基尼值成反比。对于节点的分裂，假定用属性v对集合s进行分裂，其中属性v为具有n个可能取值的离散属性{v1,v2,...,vn}，则分裂后的集合s会产生n个分支，其中每个分支对应于原始集合s中在属性v上取值相同的所有样本，记为sv。则根据公式可计算出sv的基尼纯度，再对每个分支赋予权重|sv|/|s|，其中|sv|为第v个分支的样本数量，|s|表示所有样本数量。于是，可计算出用属性v对集合s进行分裂所获得的基尼指数：
从公式中可以看出，基尼指数值越小，意味着所有分支节点的加权纯度越小，对集合s纯度的提升越明显。因此，在候选属性集合a中，选择使得分裂后gindex值最小的属性作为最佳分裂属性，即：
(1)综上所述，训练随机森林分类器的具体方法如下：
a)对于训练集中的所有图像，对每个图像按照步骤2)进行颜色校正然后按照步骤3)对校正图像进行图像分割，最后对分割图像的各个凸包按照步骤4)提取特征并人工标记样本；
b)在舌象训练样本中用bagging方法形成n个子样本；
c)对于每个子样本，随机选择p个属性作为节点分裂的候选属性，其中p取值为所有属性的开方；
d)在候选属性中计算基尼指数，选择吉尼指数最小的属性对决策树进行分裂；
f)重复步骤d)直到吉尼指数小于某一规定阈值；
g)重复步骤c)到f)，直到生成n棵决策树；
h)对于未知样本的分类，随机森林输出其决策树多数投票结果。
(2)对步骤4)提取到的各个凸包的特征，用步骤5)第(1)部分训练好的随机森林分类器输出最终舌体位置，其具体步骤如下：
a)对于各个凸包，随机森林分类器输出其为舌体的可能值pi；
b)所有pi值里最大值对应的凸包即为最终舌体区域d。
本实施例中，特征维度为48维，随机选择属性数p＝7，树的最大分裂深度为25，吉尼指数阈值为0.01，随机森林决策树个数为100。
6)以步骤5)得到的舌体区域d为基础，对图像b进行舌象图像分割，分割结果如图3所示。舌象图像分割的具体方法如下：
(1)将图像b转换到灰度空间图像b’；
(2)截取图像b’中的舌象图像分割区域，具体步骤如下：
a)计算舌体区域d对应的最小右边界矩形r；
b)将最小右边界矩形r向外扩散20个像素点，得到舌象图像分割区域r’。
(3)对舌象图像分割区域r’进行边缘增强，具体步骤如下：
a)计算舌象图像分割区域r’水平与垂直的梯度矢量场，分别为u和v；
b)将水平与垂直梯度矢量场叠加到舌象图像分割区域r’，得到边缘增强后的最终舌象图像分割区域r”，即区域r”的最终灰度取值为r’+u+v。
(4)以舌体区域d对应的凸包为初始轮廓，采用gvfsnake模型对最终舌象图像分割区域r”进行分割。
所述gvfsnake模型定义为：
在传统主动轮廓模型(snake)的基础上，采用梯度矢量流作为曲线演化的外部能量。传统snake模型是在图像区域内定义的一条轮廓曲线，其参数表达形式为x(s)＝[x(s),y(s)],s∈[0,1]，其中s是用傅里叶变换形式描述轮廓边界的变量，x(s)和y(s)表示曲线在图像区域内的坐标。该轮廓曲线的运动形式是通过最小化能量函数：
其中，α和β分别为曲线弹性和刚性的加权系数，相应的，x′(s)和x″(s)分别为曲线x(s)关于s的一阶导数和二阶导数。
要使能量函数e最小化，曲线x(s)必须满足欧拉(euler)方程：
为了求解该方程，引入时间变量t，将x(s)视为s和时间t的函数x(s,t)，则上述动态演化方程可表示为对t求偏微分：
当xt(s,t)稳定时，xt(s,t)＝0，满足欧拉方程，曲线演化完毕。
在上述公式的基础上，gvfsnake模型定义一个新的静态外部能量
xt(s,t)＝αx″(s,t)-βx″″(s,t)+v
对于新的静态外部能量gvf，相应的，将梯度矢量场定义为v(x,y)＝(u(x,y),v(x,y))，通过最小化下述能量方程得到：
其中，下标x和y表示的是图像在x和y方向上的偏导数。参数μ是控制能量函数在第一项和第二项之间权衡的正则化参数，根据图像中的噪声而定，当图像中的噪声较大时，μ取大值。
相应地，最小化能量方程可以转化为求解欧拉方程组：
其中，▽2为拉普拉斯算子。
进一步的，类似于传统snake模型的求解，欧拉方程组的求解可将分量u和v看作是时间t的函数：
当方程组趋于稳定的解即为欧拉方程组的期望解。由于方程组是解耦的，因此可以作为单独的标量偏导数求解方程。为了方便求解，改写方程组为：
其中，
本实施例中，α＝0.05、β＝0.1、μ＝0.2。
以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此。
- 还没有人留言评论。精彩留言会获得点赞！