CN107194937B - 一种开放环境下中医舌象图像分割方法 - Google Patents

一种开放环境下中医舌象图像分割方法 Download PDF

Info

Publication number
CN107194937B
CN107194937B CN201710392075.2A CN201710392075A CN107194937B CN 107194937 B CN107194937 B CN 107194937B CN 201710392075 A CN201710392075 A CN 201710392075A CN 107194937 B CN107194937 B CN 107194937B
Authority
CN
China
Prior art keywords
image
tongue
segmentation
convex hull
region
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201710392075.2A
Other languages
English (en)
Other versions
CN107194937A (zh
Inventor
黄晓阳
王彦晖
郑丰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xiamen University
Original Assignee
Xiamen University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xiamen University filed Critical Xiamen University
Priority to CN201710392075.2A priority Critical patent/CN107194937B/zh
Publication of CN107194937A publication Critical patent/CN107194937A/zh
Application granted granted Critical
Publication of CN107194937B publication Critical patent/CN107194937B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/11Region-based segmentation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/24Classification techniques
    • G06F18/243Classification techniques relating to the number of classes
    • G06F18/24323Tree-organised classifiers
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/12Edge-based segmentation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/136Segmentation; Edge detection involving thresholding
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/149Segmentation; Edge detection involving deformable models, e.g. active contour models
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/187Segmentation; Edge detection involving region growing; involving region merging; involving connected component labelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/40Analysis of texture
    • G06T7/41Analysis of texture based on statistical description of texture
    • G06T7/45Analysis of texture based on statistical description of texture using co-occurrence matrix computation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/90Determination of colour characteristics
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20076Probabilistic image processing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20081Training; Learning
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20112Image segmentation details
    • G06T2207/20116Active contour; Active surface; Snakes
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20172Image enhancement details
    • G06T2207/20192Edge enhancement; Edge preservation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30004Biomedical image processing

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Data Mining & Analysis (AREA)
  • Software Systems (AREA)
  • Probability & Statistics with Applications (AREA)
  • Mathematical Physics (AREA)
  • Computing Systems (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Artificial Intelligence (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Image Analysis (AREA)
  • Other Investigation Or Analysis Of Materials By Electrical Means (AREA)

Abstract

一种开放环境下中医舌象图像分割方法,涉及医学图像处理。提供可快速并准确地分割出可供后续分析的舌象,以便于进行后续舌象分析的一种开放环境下中医舌象图像分割方法。1)输入在开放环境下采集的图像A;2)对采集的图像A进行颜色校正,得到校正图像B;3)对步骤2)得到的校正图像B进行图像初步分割;4)对步骤3)得到的图像C中的各个连通域进行特征提取;5)对步骤4)提取到的各个凸包的特征,用随机森林分类器计算舌体位置;6)以步骤5)得到的舌体区域D为基础,对图像B进行开放环境下中医舌象图像分割方法。

Description

一种开放环境下中医舌象图像分割方法
技术领域
本发明涉及医学图像处理,尤其是涉及一种开放环境下中医舌象图像分割方法。
背景技术
传统的中医舌诊,主要依靠医生的肉眼观察进行分析判断,诊断结果往往与医生的知识水平和经验积累有关,近年来,信息技术的发展推动了中医舌象分析客观化、数字化和自动化的进程。国内外学者对此进行了许多有益的探索,并开发了一些舌象分析系统,取得了较好的效果(余兴龙,竺子民,金国潘,等.中医舌诊自动识别系统[J].仪器仪表学报,1994,15(1):67-71;蒋依吾.电脑化中医舌诊系统[J].中国中西医结合杂志,2000,20(2):145-147;Cai Y.A novel imaging system for tongue inspection[C].IEEE Imtc,Alaska.IEEE,2002:159–163;Zhang H Z,Wang K Q,Zhang D,Pang B,Huang B.Computeraided tongue diagnosis system[J].Conf Proc IEEE Eng Med Bio Soc,2005,7:6754-6757),但这些数字化舌诊分析系统的拍摄环境通常是固定的,即在密闭、光照稳定的环境下拍摄采集,这样能够得到高质量的舌象图像,但存在一定局限性,表现在采集设备的特殊性与不易携带。在舌象图像分割方面,香港理工大学David.Zhang等人在这一领域取得了较新的研究成果,进行了许多有益的探索(Pang B,Zhang D,Wang K.The bi-ellipticaldeformable contour and its application to automated tongue segmentation inChinese medicine[J].IEEE Transactions on Medical Imaging,2005,24(8):946;NingJ,Zhang D,Wu C,et al.Automatic tongue image segmentation based on gradientvector flow and region merging[J].Neural Computing and Applications,2012,21(8):1-8;Cui Z,Zuo W,Zhang H,et al.Automated Tongue Segmentation Based on 2DGabor Filters and Fast Marching[M]//Intelligence Science and Big DataEngineering.2013:328-335;WU K,ZHANG D.Robust tongue segmentation by fusingregion-based and edge-based approaches[J].Expert Systems with Applications,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
Figure BDA0001307749070000031
Figure BDA0001307749070000032
中,舌象图像S1的整体均值KS=(Ravgs+Gavgs+Bavgs)/3;Ravgs、Gavgs、Bavgs分别为标准光照环境下采集舌象图像S1的RGB三个颜色通道的均值;所述标准光照采集环境可为D65光源,色温为6500K;(2)按下式调整图像A的RGB三个颜色通道的均值,得到校正图像B;
K=(Ravg+Gavg+Bavg)/3
Figure BDA0001307749070000033
Figure BDA0001307749070000034
Figure BDA0001307749070000035
其中,K为图像A的整体均值,Ravg、Gavg、Bavg分别为图像A的RGB三个颜色通道的均值;Rd、Gd、Bd为校正图像B每个像素点的RGB三个颜色通道的值,Rs、Gs、Bs为图像A每个像素点的RGB三个颜色通道的值;
在步骤3)中,所述对步骤2)得到的校正图像B进行图像初步分割的具体方法可为:
(1)将校正图像B转换到灰度空间图像fB1,按照最大类间方差法(参见文献:OHTSUN.A threshold selection method from gray-level histograms[J].System Man&Cybernetics IEEE Transactions on,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的面积比
Figure BDA0001307749070000041
删除面积比小于0.02的区域;
b)计算连通域质心Ci离图像中心C0的欧式距离
Figure BDA0001307749070000042
其中Cix、Ciy分别表示凸包Si质心的横坐标与纵坐标,C0x、C0y分别表示图像C中心的横坐标与纵坐标;
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]。
Figure BDA0001307749070000051
Figure BDA0001307749070000052
V={0v∈(0.15,1]
按照以上的量化等级,将H、S、V三个个颜色分量合成为一维特征矢量:
l=HQsQv+SQv+V
其中Qs和Qv分别是分量S和V的量化级数,这里Qs=2,Qv=1,因此上式可表示为:
Figure BDA0001307749070000053
根据上式,l取值范围为[0,1,…,19],计算l获得20柄的一维直方图,通过统计并归一化我们能够得到20个颜色特征。
(4)提取凸包区域的形状特征;
在步骤4)第(4)部分中,所述提取凸包区域的形状特征具体为计算凸包Si的7个Hu不变矩mi(i∈[1,7])(参见文献:Hu M.Visual Pattern Recognition by MomentInvariants[J].Information Theory Ire Transactions on,1962,8(2):179-187.)。
所述Hu不变矩定义如下:
对于大小为M×N的图像f(x,y),f(x,y)的二维(p+q)阶矩定义为:
Figure BDA0001307749070000054
相应的(p+q)阶中心矩定义为:
Figure BDA0001307749070000061
其中
Figure BDA0001307749070000065
由ηpq表示的归一化中心矩定义为:
Figure BDA0001307749070000062
其中γ=(p+q)/2+1
Hu用归一化中心矩的线性组合构造具有平移、伸缩、旋转不变的7个Hu不变矩
Figure BDA0001307749070000063
Figure BDA0001307749070000064
(5)提取凸包区域的纹理特征;
在步骤4)第(5)部分中,所述提取凸包区域的纹理特征具体具体步骤如下:
a)将校正图像B转换到灰度空间图像C”;
b)凸包Si对应于图像C”中的区域,分别计算区域的水平、垂直、对角偏移量三种灰度
共生矩阵的联合概率密度分布;
设灰度共生矩阵的联合概率密度分布记为[Pmn]L×L,其中L为灰度取值范围,m=[0,L-1],n=[0,L-1];根据联合概率密度分布计算得到6个纹理特征(参见文献:HARALICKR M,SHANMUGAM K,DINSTEIN I.Textural features for image classification[J].Systems Man&Cybernetics IEEE Transactions on,2010,smc-3(6):610–621.),所述6个纹理特征包括角二阶矩PCON、对比度PASM、熵PENT、逆差矩PIDM、中值PMEAN和灰度相关PCOR
Figure BDA0001307749070000071
在步骤5)中,所述随机森林分类器的定义为:
随机森林(参见文献:Surhone L M,Tennoe M T,Henssonow S F,et al.RandomForest[J].Machine Learning,2001,45(1):5-32)是一种以决策树作为基础分类器的分类算法。随机森林分类器由多棵决策树构成,这些决策树彼此独立,而且在训练样本的选择、特征子集的选择和树的生长过程中引入了随机性,即每棵决策树的构造都是从原数据集中随机抽取出一部分作为样本子集,然后再从样本子集中随机选取部分特征子集进行处理,并且所有决策树都自然生长,不进行剪枝。因此,多颗决策树可以并行生长,这样随机森林分类器相对于其他分类器具有更短的训练和识别时间,并且由于引入了随机性,分类器对噪声也不敏感。随机森林的最终决策是通过所有决策树投票的方法而得,因此,具有更高的检测率与鲁棒性。
所述随机森林采用Bagging方法选取每个训练子集,即采用有放回的抽取方式,这样会导致原数据集中的部分样本可能重复出现在新的训练子集中,而另外一部分样本一次也不出现,这些没有被选择到的数据叫做“out of bag(OOB)”数据,一般占所有数据的1/3。
对于模型中的每个节点,决策树的构造通过随机抽取特征子集来分裂生长,随机抽取的数量为特征总数的开方。对于分类问题,决策树生长过程中节点最佳分裂的度量选择基尼指数,对于一个含有m类样本的集合S,集合S的纯度用基尼值度量为:
Figure BDA0001307749070000072
其中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进行分裂所获得的基尼指数:
Figure BDA0001307749070000081
从上式中可以看出,基尼指数值越小,意味着所有分支节点的加权纯度越小,对集合S纯度的提升越明显。因此,在候选属性集合A中,选择使得分裂后Gindex值最小的属性作为最佳分裂属性,即:
Figure BDA0001307749070000082
训练随机森林分类器的具体方法可为:
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对应的凸包为初始轮廓,采用GVF Snake模型(参见文献:Xu C,Prince J L.Snakes,shapes,and gradient vector flow[J].IEEE Transactions onImage Processing A Publication of the IEEE Signal Processing Society,1998,7(3):359-69.)对最终舌象图像分割区域R”进行分割,得到舌象分割图像。
所述GVF Snake模型定义为:
在传统主动轮廓模型(Snake)的基础上,采用梯度矢量流作为曲线演化的外部能量。传统Snake模型是在图像区域内定义的一条轮廓曲线,其参数表达形式为X(s)=[x(s),y(s)],s∈[0,1],其中s是用傅里叶变换形式描述轮廓边界的变量,x(s)和y(s)表示曲线在图像区域内的坐标。该轮廓曲线的运动形式是通过最小化能量函数:
Figure BDA0001307749070000091
其中,α和β分别为曲线弹性和刚性的加权系数,相应的,X′(s)和X″(s)分别为曲线X(s)关于s的一阶导数和二阶导数。
要使能量函数E最小化,曲线X(s)必须满足欧拉(Euler)方程:
Figure BDA0001307749070000092
为了求解该方程,引入时间变量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,满足欧拉方程,曲线演化完毕。
在上述公式的基础上,GVF Snake模型定义一个新的静态外部能量
Figure BDA0001307749070000093
即梯度矢量流。取代演化方程中的
Figure BDA0001307749070000094
新的动态演化方程可表示为:
Xt(S,t)=αX″(s,t)-βX″″(s,t)+V
对于新的静态外部能量GVF,相应的,将梯度矢量场定义为V(x,y)=(u(x,y),v(x,y)),通过最小化下述能量方程得到:
Figure BDA0001307749070000101
其中,下标x和y表示的是图像在x和y方向上的偏导数。参数μ是控制能量函数在第一项和第二项之间权衡的正则化参数,根据图像中的噪声而定,当图像中的噪声较大时,μ取大值。
相应地,最小化能量方程可以转化为求解欧拉方程组:
Figure BDA0001307749070000102
其中,▽2为拉普拉斯算子。
进一步的,类似于传统Snake模型的求解,欧拉方程组的求解可将分量u和v看作是时间t的函数:
Figure BDA0001307749070000103
当方程组趋于稳定的解即为欧拉方程组的期望解。由于方程组是解耦的,因此可以作为单独的标量偏导数求解方程。为了方便求解,改写方程组为:
Figure BDA0001307749070000104
其中,
Figure BDA0001307749070000105
c1(x,y)=b(x,y)+fx(x,y),c2(x,y)=b(x,y)+fy(x,y)。fx(x,y)和fy(x,y)可由任意的图像梯度算子求得,进而求得b(x,y)、c1(x,y)和c2(x,y),作为迭代过程中固定系数。
本发明针对开放环境下中医舌象图像的特点对图像进行舌象图像分割,首先对图像进行颜色校正的预处理,减少因外界光源色温带来的影响;然后为了确定初始轮廓,先对图像进行分割,得到多个连通区域;并对各连通区域进行特征提取,最终通过训练好的随机森林分类器得到舌体初始轮廓;最后通过GVF Snake模型对初始轮廓进行演化,逼近舌体真实边缘。本发明最终达到对图像中的舌象进行分割的目的,即将舌体从复杂背景中提取出来。
颜色校正过程采用改进的灰度世界算法,针对舌象图像特有的三个颜色通道均值不同的特征调整相应的参数值,以符合舌象图像的特征;图像分割部分采用最大类间方差、色调阈值分割和RGB三色分量方差分割相结合的方法进行分割;确定初始轮廓部分通过对各连通域的凸包进行特征提取后,采用训练好的随机森林分类器得到舌体区域的轮廓;舌象图像分割则以上一部的舌体轮廓作为初始轮廓,采用GVF Snake模型对曲线进行演化,逼近舌体真实轮廓,最终达到舌象图像分割的目的。
相对于现有技术,本发明的有益效果为:
1、本发明根据舌象特征利用改进的灰度世界算法对图像进行颜色校正,减少了因采集环境颜色偏色带来的影响。
2、本发明采用最大类间方差、色调阈值分割和RGB三色分量方差分割相结合的方法进行图像分割;根据连通域凸包的空间位置特征、颜色特征、形状特征和轮廓特征,采用随机森林分类器进行提取舌体初始轮廓;采用GVF Snake模型进行舌象图像精细分割;这些方法都是针对开放环境下的舌象特点进行处理、分析,这些方法相辅相成,环环相扣,结合起来使用能够更加精确地分割出舌象,具有分割精度高、鲁棒性强的优点。
3、本发明使用的方法运算量小、程序实现简单,随机森林分类器不需要训练大量数据,特征提取算法复杂度低,处理速度快,对数据集的依赖程度不高。
4、本发明识别、分割的准确率高,综合考虑了舌象图像空间位置、颜色、形状和纹理等多种特征,综合考虑图像整体信息,相比其他方法采用单一特征而言,本发明具有更高的识别率和准确率。
本发明针对开放环境下舌象图像的特点分多个步骤进行处理、分析,达到舌象图像分割的目的。本发明的图像分割、特征提取、分类器识别、GVF Snake模型分割几个步骤在功能上相辅相成,图像分割主要利用阈值分割的方法进行粗筛选,算法复杂度低,特征提取利用图像的空间位置、颜色、形状和轮廓特征,算法实现简单,随机森林分类器识别舌体初始轮廓,鲁棒性强,舌象图像分割在初始轮廓的基础上利用主动轮廓模型进行分割,分割的精确度较高,运算速度快,取得了新的技术效果。相对于最接近的现有技术,具有明显的区别特征,且具有较明显的优点:本发明的方法不需要采用穷举的方法转换到霍夫空间投票,算法复杂度低,运算量小,也不需要进行有监督的训练大量数据,因此实现较为简单,对数据集的依赖程度不高,在准确率与处理速度上更优越,更具有实用性。
附图说明
图1为开放环境下采集的图像;
图2为本发明步骤5)得到的舌体初始轮廓;(a线为初始轮廓)
图3为本发明步骤6)舌象图像分割的最终结果
具体实施方式
下面结合附图和实施例对本发明作进一步详细描述。
本发明实施例包括如下步骤:
1)输入在开放环境下采集的图像A,如图1所示;
2)对采集的图像A进行颜色校正,得到校正图像B;
在步骤2)中,所述对采集的图像A进行颜色校正的具体方法如下:
(1)在标准光照环境下采集舌象图像S1,计算舌象图像S1的RGB三个颜色通道均值分别与舌象图像S1的整体均值KS的比值αr、αg、αb
Figure BDA0001307749070000121
Figure BDA0001307749070000122
其中,舌象图像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
Figure BDA0001307749070000123
Figure BDA0001307749070000124
Figure BDA0001307749070000125
其中,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):
Figure BDA0001307749070000131
其中,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):
Figure BDA0001307749070000132
其中,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
Figure BDA0001307749070000141
其中,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的面积比
Figure BDA0001307749070000142
删除面积比小于0.02的区域;
b)计算连通域质心Ci离图像中心C0的欧式距离
Figure BDA0001307749070000151
其中Cix、Ciy分别表示凸包Si质心的横坐标与纵坐标,C0x、C0y分别表示图像C中心的横坐标与纵坐标;
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]。
Figure BDA0001307749070000152
Figure BDA0001307749070000153
V={0v∈(0.15,1]
按照以上的量化等级,将H、S、V三个个颜色分量合成为一维特征矢量:
l=HQsQv+SQv+V
其中Qs和Qv分别是分量S和V的量化级数,这里Qs=2,Qv=1,因此上式可表示为:
Figure BDA0001307749070000154
根据上式,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)阶矩定义为:
Figure BDA0001307749070000161
相应的(p+q)阶中心矩定义为:
Figure BDA0001307749070000162
其中
Figure BDA0001307749070000163
由ηpq表示的归一化中心矩定义为:
Figure BDA0001307749070000164
其中γ=(p+q)/2+1
Hu用归一化中心矩的线性组合构造具有平移、伸缩、旋转不变的7个Hu不变矩
Figure BDA0001307749070000165
Figure BDA0001307749070000166
本实施例中,共获得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
Figure BDA0001307749070000171
本实施例中,共获得18维纹理特征
5)对步骤4)提取到的各个凸包的特征,用训练好的随机森林分类器计算舌体位置,如图2所示;
所述随机森林分类器定义为:
随机森林是一种以决策树作为基础分类器的分类算法。随机森林分类器由多棵决策树构成,这些决策树彼此独立,而且在训练样本的选择、特征子集的选择和树的生长过程中引入了随机性,即每棵决策树的构造都是从原数据集中随机抽取出一部分作为样本子集,然后再从样本子集中随机选取部分特征子集进行处理,并且所有决策树都自然生长,不进行剪枝。因此,多颗决策树可以并行生长,这样随机森林分类器相对于其他分类器具有更短的训练和识别时间,并且由于引入了随机性,分类器对噪声也不敏感。随机森林的最终决策是通过所有决策树投票的方法而得,因此,具有更高的检测率与鲁棒性。
随机森林采用Bagging方法选取每个训练子集,即采用有放回的抽取方式,这样会导致原数据集中的部分样本可能重复出现在新的训练子集中,而另外一部分样本一次也不出现,这些没有被选择到的数据叫做“out of bag(OOB)”数据,一般占所有数据的1/3。
对于模型中的每个节点,决策树的构造通过随机抽取特征子集来分裂生长,随机抽取的数量为特征总数的开方。对于分类问题,决策树生长过程中节点最佳分裂的度量选择基尼指数,对于一个含有m类样本的集合S,集合S的纯度用基尼值度量为:
Figure BDA0001307749070000172
其中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进行分裂所获得的基尼指数:
Figure BDA0001307749070000181
从公式中可以看出,基尼指数值越小,意味着所有分支节点的加权纯度越小,对集合S纯度的提升越明显。因此,在候选属性集合A中,选择使得分裂后Gindex值最小的属性作为最佳分裂属性,即:
Figure BDA0001307749070000182
(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对应的凸包为初始轮廓,采用GVF Snake模型对最终舌象图像分割区域R”进行分割。
所述GVF Snake模型定义为:
在传统主动轮廓模型(Snake)的基础上,采用梯度矢量流作为曲线演化的外部能量。传统Snake模型是在图像区域内定义的一条轮廓曲线,其参数表达形式为X(s)=[x(s),y(s)],s∈[0,1],其中s是用傅里叶变换形式描述轮廓边界的变量,x(s)和y(s)表示曲线在图像区域内的坐标。该轮廓曲线的运动形式是通过最小化能量函数:
Figure BDA0001307749070000191
其中,α和β分别为曲线弹性和刚性的加权系数,相应的,X′(s)和X″(s)分别为曲线X(s)关于s的一阶导数和二阶导数。
要使能量函数E最小化,曲线X(s)必须满足欧拉(Euler)方程:
Figure BDA0001307749070000192
为了求解该方程,引入时间变量t,将X(s)视为s和时间t的函数X(s,t),则上述动态演化方程可表示为对t求偏微分:
Figure BDA0001307749070000193
当Xt(S,t)稳定时,Xt(S,t)=0,满足欧拉方程,曲线演化完毕。
在上述公式的基础上,GVF Snake模型定义一个新的静态外部能量
Figure BDA0001307749070000194
即梯度矢量流。取代演化方程中的-▽Eext,新的动态演化方程可表示为:
Xt(S,t)=αX″(s,t)-βX″″(s,t)+V
对于新的静态外部能量GVF,相应的,将梯度矢量场定义为V(x,y)=(u(x,y),v(x,y)),通过最小化下述能量方程得到:
Figure BDA0001307749070000201
其中,下标x和y表示的是图像在x和y方向上的偏导数。参数μ是控制能量函数在第一项和第二项之间权衡的正则化参数,根据图像中的噪声而定,当图像中的噪声较大时,μ取大值。
相应地,最小化能量方程可以转化为求解欧拉方程组:
Figure BDA0001307749070000202
其中,▽2为拉普拉斯算子。
进一步的,类似于传统Snake模型的求解,欧拉方程组的求解可将分量u和v看作是时间t的函数:
Figure BDA0001307749070000203
当方程组趋于稳定的解即为欧拉方程组的期望解。由于方程组是解耦的,因此可以作为单独的标量偏导数求解方程。为了方便求解,改写方程组为:
Figure BDA0001307749070000204
其中,
Figure BDA0001307749070000205
c1(x,y)=b(x,y)+fx(x,y),c2(x,y)=b(x,y)+fy(x,y)。fx(x,y)和fy(x,y)可由任意的图像梯度算子求得,进而求得b(x,y)、c1(x,y)和c2(x,y),作为迭代过程中固定系数。
本实施例中,α=0.05、β=0.1、μ=0.2。
以上所述,仅为本发明的具体实施方式,但本发明的保护范围并不局限于此。

Claims (7)

1.一种开放环境下中医舌象图像分割方法,其特征在于包括以下步骤:
1)输入在开放环境下采集的图像A;
2)对采集的图像A进行颜色校正,得到校正图像B;
3)对步骤2)得到的校正图像B进行图像初步分割;
4)对步骤3)得到的图像C中的各个连通域进行特征提取;
5)对步骤4)提取到的各个凸包的特征,用随机森林分类器计算舌体位置;
所述随机森林分类器的定义为:
随机森林以决策树作为基础分类器的分类算法,随机森林分类器由多棵决策树构成,在训练样本的选择、特征子集的选择和树的生长过程中引入随机性,即每棵决策树的构造都是从原数据集中随机抽取出一部分作为样本子集,再从样本子集中随机选取部分特征子集进行处理,所有决策树都自然生长,多颗决策树并行生长,随机森林的最终决策是通过所有决策树投票的方法而得,具有更高的检测率与鲁棒性;
所述随机森林采用Bagging方法选取每个训练子集,即采用有放回的抽取方式;
对于模型中的每个节点,决策树的构造通过随机抽取特征子集来分裂生长,随机抽取的数量为特征总数的开方,对于分类问题,决策树生长过程中节点最佳分裂的度量选择基尼指数,对于一个含有m类样本的集合S,集合S的纯度用基尼值度量为:
Figure FDA0002177746280000011
其中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进行分裂所获得的基尼指数:
Figure FDA0002177746280000012
从上式中看出,基尼指数值越小,意味着所有分支节点的加权纯度越小,对集合S纯度的提升越明显,在候选属性集合A中,选择使得分裂后Gindex值最小的属性作为最佳分裂属性,即:
Figure FDA0002177746280000021
所述训练随机森林分类器的具体方法为:
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对应的凸包为初始轮廓,采用GVF Snake模型对最终舌象图像分割区域R”进行分割,得到舌象分割图像;
所述GVF Snake模型定义为:
在传统主动轮廓模型的基础上,采用梯度矢量流作为曲线演化的外部能量;传统Snake模型是在图像区域内定义的一条轮廓曲线,其参数表达形式为X(s)=[x(s),y(s)],s∈[0,1],其中s是用傅里叶变换形式描述轮廓边界的变量,x(s)和y(s)表示曲线在图像区域内的坐标,该轮廓曲线的运动形式是通过最小化能量函数:
Figure FDA0002177746280000031
其中,α和β分别为曲线弹性和刚性的加权系数,相应的,X′(s)和X″(s)分别为曲线X(s)关于s的一阶导数和二阶导数;
要使能量函数E最小化,曲线X(s)满足欧拉方程:
Figure FDA0002177746280000032
为了求解该方程,引入时间变量t,将X(s)视为s和时间t的函数X(s,t),则上述动态演化方程表示为对t求偏微分:
Figure FDA0002177746280000033
当Xt(S,t)稳定时,Xt(S,t)=0,满足欧拉方程,曲线演化完毕;
在上述公式的基础上,GVF Snake模型定义一个新的静态外部能量
Figure FDA0002177746280000034
即梯度矢量流,取代演化方程中的
Figure FDA0002177746280000035
新的动态演化方程表示为:
Figure FDA0002177746280000036
对于新的静态外部能量
Figure FDA0002177746280000037
相应将梯度矢量场定义为V(x,y)=(u(x,y),v(x,y)),通过最小化下述能量方程得到:
Figure FDA0002177746280000038
其中,下标x和y表示的是图像在x和y方向上的偏导数,参数μ是控制能量函数在第一项和第二项之间权衡的正则化参数,根据图像中的噪声而定,当图像中的噪声较大时,μ取大值;
相应地,最小化能量方程转化为求解欧拉方程组:
Figure FDA0002177746280000041
其中,
Figure FDA0002177746280000042
为拉普拉斯算子;
进一步的,类似于传统Snake模型的求解,欧拉方程组的求解将分量u和v看作是时间t的函数:
Figure FDA0002177746280000043
当方程组趋于稳定的解即为欧拉方程组的期望解,由于方程组是解耦的,因此作为单独的标量偏导数求解方程,改写方程组为:
Figure FDA0002177746280000044
其中,
Figure FDA0002177746280000045
c1(x,y)=b(x,y)*fx(x,y),c2(x,y)=b(x,y)*fy(x,y);fx(x,y)和fy(x,y)由任意的图像梯度算子求得,进而求得b(x,y)、c1(x,y)和c2(x,y),作为迭代过程中固定系数。
2.如权利要求1所述一种开放环境下中医舌象图像分割方法,其特征在于在步骤2)中,所述对采集的图像A进行颜色校正的具体方法如下:
(1)在标准光照环境下采集舌象图像S1,计算舌象图像S1的RGB三个颜色通道均值分别与舌象图像S1的整体均值KS的比值αr、αg、αb
Figure FDA0002177746280000046
Figure FDA0002177746280000047
其中,舌象图像S1的整体均值KS=(Ravgs+Gavgs+Bavgs)/3;
Ravgs、Gavgs、Bavgs分别为标准光照环境下采集舌象图像S1的RGB三个颜色通道的均值;
(2)按下式调整图像A的RGB三个颜色通道的均值,得到校正图像B;
K=(Ravg+Gavg+Bavg)/3
Figure FDA0002177746280000048
Figure FDA0002177746280000051
Figure FDA0002177746280000052
其中,K为图像A的整体均值,Ravg、Gavg、Bavg分别为图像A的RGB三个颜色通道的均值;Rd、Gd、Bd为校正图像B每个像素点的RGB三个颜色通道的值,Rs、Gs、Bs为图像A每个像素点的RGB三个颜色通道的值。
3.如权利要求1所述一种开放环境下中医舌象图像分割方法,其特征在于在步骤3)中,所述对步骤2)得到的校正图像B进行图像初步分割的具体方法为:
(1)将校正图像B转换到灰度空间图像fB1,按照最大类间方差法对灰度空间图像fB1进行阈值分割,得到分割图像B1’,并对分割图像B1’运用形态学运算平滑连通域,得到图像B1;
(2)将校正图像B转换到HSV颜色空间图像fB2,对图像fB2的H通道进行阈值分割,得到分割图像B2’,并对分割图像B2’运用形态学运算平滑连通域,得到图像B2;
(3)采用RGB三色分量方差法对校正图像B进行阈值分割,得到分割图像B3’,并对分割图像B3’运用形态学运算平滑连通域,得到图像B3;
(4)对图像B1、图像B2、图像B3三个图像进行逻辑“与”运算,得到图像C。
4.如权利要求3所述一种开放环境下中医舌象图像分割方法,其特征在于在步骤3)第(1)、(2)、(3)部分中,所述将分割结果图像运用形态学运算平滑连通域,其具体步骤为对分割结果图像B(x,y)进行先腐蚀再膨胀:
B(x,y)=dilate(erode(B(x,y),element))
其中,element定义为形态学运算中的结构元素;dilate定义为形态学运算中的膨胀操作;erode定义为形态学运算中的腐蚀操作。
5.如权利要求1所述一种开放环境下中医舌象图像分割方法,其特征在于在步骤4)中,所述对步骤3)得到的图像C中的各个连通域进行特征提取的具体方法为:
(1)对于图像C中的每一个连通域,计算该连通域的凸包Si
(2)提取凸包区域的空间位置特征;
(3)提取凸包区域的颜色特征;
(4)提取凸包区域的形状特征;
(5)提取凸包区域的纹理特征。
6.如权利要求5所述一种开放环境下中医舌象图像分割方法,其特征在于在步骤4)第(2)部分中,所述提取凸包区域的空间位置特征具体步骤如下:
a)计算凸包Si与图像C的面积比
Figure FDA0002177746280000061
删除面积比小于0.02的区域;
b)计算连通域质心Ci离图像中心C0的欧式距离
Figure FDA0002177746280000062
其中Cix、Ciy分别表示凸包Si质心的横坐标与纵坐标,C0x、C0y分别表示图像C中心的横坐标与纵坐标;
c)计算连通域最小右边界矩形的长宽比Rscale=W/H,其中W表示最小外接矩形的宽,H表示最小外接矩形的长。
7.如权利要求5所述一种开放环境下中医舌象图像分割方法,其特征在于在步骤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];
Figure FDA0002177746280000063
Figure FDA0002177746280000064
V={0 v∈(0.15,1]}
按照以上的量化等级,将H、S、V三个个颜色分量合成为一维特征矢量:
l=HQsQv+SQv+V
其中Qs和Qv分别是分量S和V的量化级数,Qs=2,Qv=1,因此上式表示为:
Figure FDA0002177746280000071
根据上式,l取值范围为[0,1,…,19],计算l获得20柄的一维直方图,通过统计并归一化能够得到20个颜色特征;
在步骤4)第(4)部分中,所述提取凸包区域的形状特征具体为计算凸包Si的7个Hu不变矩mi,其中i∈[1,7];
所述Hu不变矩定义如下:
对于大小为M×N的图像f(x,y),f(x,y)的二维(p+q)阶矩定义为:
Figure FDA0002177746280000072
相应的(p+q)阶中心矩定义为:
Figure FDA0002177746280000073
其中
Figure FDA0002177746280000074
由ηpq表示的归一化中心矩定义为:
Figure FDA0002177746280000075
其中γ=(p+q)/2+1
Hu用归一化中心矩的线性组合构造具有平移、伸缩、旋转不变的7个Hu不变矩
Figure FDA0002177746280000076
其中i∈[1,7];
Figure FDA0002177746280000077
在步骤4)第(5)部分中,所述提取凸包区域的纹理特征具体具体步骤如下:
a)将校正图像B转换到灰度空间图像C”;
b)凸包Si对应于图像C”中的区域,分别计算区域的水平、垂直、对角偏移量三种灰度共生矩阵的联合概率密度分布;
c)设灰度共生矩阵的联合概率密度分布记为[Pmn]L×L,其中L为灰度取值范围,m=[0,L-1],n=[0,L-1];根据联合概率密度分布计算得到6个纹理特征,所述6个纹理特征包括角二阶矩PCON、对比度PASM、熵PENT、逆差矩PIDM、中值PMEAN和灰度相关PCOR
Figure FDA0002177746280000081
CN201710392075.2A 2017-05-27 2017-05-27 一种开放环境下中医舌象图像分割方法 Active CN107194937B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710392075.2A CN107194937B (zh) 2017-05-27 2017-05-27 一种开放环境下中医舌象图像分割方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710392075.2A CN107194937B (zh) 2017-05-27 2017-05-27 一种开放环境下中医舌象图像分割方法

Publications (2)

Publication Number Publication Date
CN107194937A CN107194937A (zh) 2017-09-22
CN107194937B true CN107194937B (zh) 2020-04-24

Family

ID=59875989

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710392075.2A Active CN107194937B (zh) 2017-05-27 2017-05-27 一种开放环境下中医舌象图像分割方法

Country Status (1)

Country Link
CN (1) CN107194937B (zh)

Families Citing this family (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108090507A (zh) * 2017-10-19 2018-05-29 电子科技大学 一种基于集成方法的医疗影像纹理特征处理方法
CN108492305B (zh) * 2018-03-19 2020-12-22 深圳牙领科技有限公司 一种嘴唇内侧轮廓线的分割方法、系统和介质
CN109460781B (zh) * 2018-10-20 2022-08-19 上海瑞轩食品有限公司 一种基于决策树归纳学习的牛排等级划分方法
CN111091371A (zh) * 2018-10-24 2020-05-01 北京意锐新创科技有限公司 快速支付方法和装置
CN109636864A (zh) * 2018-12-19 2019-04-16 新绎健康科技有限公司 一种基于颜色校正与深度卷积神经网络的舌分割方法及系统
CN109785345A (zh) * 2019-01-25 2019-05-21 中电健康云科技有限公司 图像分割方法及装置
CN109978895B (zh) * 2019-03-29 2021-08-13 北京峰云视觉技术有限公司 一种舌体图像分割方法和装置
CN110175998A (zh) * 2019-05-30 2019-08-27 沈闯 基于多尺度深度学习的乳腺癌图像识别方法、装置及介质
CN111340773B (zh) * 2020-02-24 2022-11-22 齐鲁工业大学 一种视网膜图像血管分割方法
CN111667500A (zh) * 2020-06-04 2020-09-15 天津市天中依脉科技开发有限公司 一种基于图像块先验的舌象分割算法
CN113840135B (zh) * 2021-09-03 2023-10-20 大连中科创达软件有限公司 色偏检测方法、装置、设备及存储介质
CN114037011B (zh) * 2021-11-08 2024-05-28 北京工业大学 一种中医舌色噪声标注样本的自动识别与清洗方法
CN116503392B (zh) * 2023-06-26 2023-08-25 细胞生态海河实验室 一种用于卵巢组织分析的卵泡区域分割方法
CN116542983B (zh) * 2023-07-07 2023-09-26 季华实验室 平顶光束特征检测方法、装置、设备和计算机存储介质

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101799870A (zh) * 2009-02-06 2010-08-11 天津市医学堂科技有限公司 舌象轮廓提取方法及其应用
CN102426583A (zh) * 2011-10-10 2012-04-25 北京工业大学 基于图像内容分析的中医舌象检索方法
CN103735253A (zh) * 2014-01-17 2014-04-23 厦门强本科技有限公司 一种基于移动终端的中医舌象分析系统及方法
CN106023151A (zh) * 2016-05-09 2016-10-12 厦门大学 一种开放环境下中医舌象目标检测方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101799870A (zh) * 2009-02-06 2010-08-11 天津市医学堂科技有限公司 舌象轮廓提取方法及其应用
CN102426583A (zh) * 2011-10-10 2012-04-25 北京工业大学 基于图像内容分析的中医舌象检索方法
CN103735253A (zh) * 2014-01-17 2014-04-23 厦门强本科技有限公司 一种基于移动终端的中医舌象分析系统及方法
CN106023151A (zh) * 2016-05-09 2016-10-12 厦门大学 一种开放环境下中医舌象目标检测方法

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
An Assessment Method of Tongue Image Quality Based on Random Forest in Traditional Chinese Medicine;Xinfeng Zhang et al;《ICIC 2015: Advanced Intelligent Computing Theories and Applications》;20150813;第730-737页 *
Application of Image Segmentation Technique in Tongue Diagnosis;ZHAI Xueming et al;《2009 International Forum on Information Technology and Applications》;20090517;第768-771页 *
Random Forests;Leo Breiman et al;《Machine Learning》;20011031;第45卷(第1期);第5-32页 *
Snakes, shapes, and gradient vector flow;Chenyang Xu et al;《IEEE Transactions on Image Processing》;19980331;第7卷(第3期);第359-369页 *
一种2型糖尿病中医证型的舌图像识别方法;阚红星等;《中国生物医学工程学报》;20161231;第35卷(第6期);摘要及正文第1.7节 *

Also Published As

Publication number Publication date
CN107194937A (zh) 2017-09-22

Similar Documents

Publication Publication Date Title
CN107194937B (zh) 一种开放环境下中医舌象图像分割方法
CN109154978B (zh) 用于检测植物疾病的系统和方法
CN108364288B (zh) 用于乳腺癌病理图像的分割方法和装置
US10839510B2 (en) Methods and systems for human tissue analysis using shearlet transforms
WO2022001571A1 (zh) 一种基于超像素图像相似度的计算方法
CN109978848B (zh) 基于多光源颜色恒常模型检测眼底图像中硬性渗出的方法
CN108537751B (zh) 一种基于径向基神经网络的甲状腺超声图像自动分割方法
CN111882561A (zh) 一种癌细胞识别诊断系统
Prajapati et al. Brain tumor detection by various image segmentation techniques with introduction to non negative matrix factorization
CN107862680B (zh) 一种基于相关滤波器的目标跟踪优化方法
Casanova et al. Texture analysis using fractal descriptors estimated by the mutual interference of color channels
CN107622280B (zh) 基于场景分类的模块化处方式图像显著性检测方法
Akshaya et al. Kidney stone detection using neural networks
Jenifa et al. Classification of cotton leaf disease using multi-support vector machine
Altaei et al. Brain tumor detection and classification using SIFT in MRI images
Anandgaonkar et al. Brain tumor detection and identification from T1 post contrast MR images using cluster based segmentation
CN114373079A (zh) 一种快速准确的探地雷达目标检测方法
CN113850792A (zh) 一种基于计算机视觉的细胞分类计数方法及系统
CN113989799A (zh) 一种宫颈异常细胞识别方法、装置和电子设备
CN109299295B (zh) 蓝印花布图像数据库搜索方法
Pankaja et al. Leaf recognition and classification using GLCM and hierarchical centroid based technique
CN115018820A (zh) 基于纹理加强的乳腺癌多分类方法
Ma et al. Deep attention network for melanoma detection improved by color constancy
Varjão et al. Citrus fruit quality classification using support vector machines
CN103871084A (zh) 蓝印花布图案识别方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant