专利名称：基于khm算法的视频显微图像细胞自动跟踪方法
技术领域：
本发明属于数字图像处理领域，涉及一种视频显微图像细胞自动跟踪的方法，该
方法可用于医学图像处理。
背景技术：
细胞动态学是众多生物医学处理中重要的一个方面，它直接引导着人类健康的走
向。细胞运动分析是细胞动态学研究中重要的组成部分之一。分析细胞运动，通常是对目
标细胞进行追踪，进而对细胞的运动速度、位移、轨迹、形状等属性进行定性或定量的分析，
它是进行细胞活性、细胞迁徙和细胞趋向性等细胞生物学和生物制药研究的有效方法和必
要手段，在生物学、药理学和病理学方面都具有十分重要的研究意义和实用价值。 目前，国内细胞跟踪的实现方法主要还是显微设备辅助下的人工定时检测和记
录，其自动化程度较低， 一方面要求工作人员付出辛苦的劳动，具有较高的劳动强度；另一
方面由于细胞形式的复杂多变性，如细胞分裂、结合、聚集、交叉和视频图像本身的噪声和
干扰等因素，使得人工观测过程更加困难，并且有较高的出错率。 近几年，国外研究机构在细胞自动跟踪方面取得了一些研究成果，借助图像分析 软件研制出一批细胞自动跟踪系统，如University of Virginia、University ofAberdeen 禾口 European Molecular Biology Laboratory等，但在商业化方面还不成熟，这些设备虽 然能够通过图像处理软件和分析软件计算出细胞在图像上的位置，但对于细胞的分裂与结 合、细胞聚合和多细胞运动轨迹的交叉等特殊情况难以处理。
发明内容
本发明的目的在于针对现有细胞跟踪方法的不足，提出一种基于KHM算法的视频 显微图像细胞自动跟踪方法，以实现对运动细胞进行自动跟踪，提高对细胞运动中出现特 殊情况时的处理能力。 实现本发明的技术方案是将细胞运动跟踪问题看作多目标跟踪问题，利用KHM
算法对预处理后的细胞二值图像进行聚类，获得所有细胞的位置，以最近邻方法为准则进
行目标细胞关联，当关联过程中出现一对多，多对一的特殊情况时，计算不同关联结果的差
异代价函数，以差异代价函数最小的关联结果作为最终关联结果，完成运动细胞的自动跟
踪，提高对细胞运动中出现特殊情况时的处理能力。具体实现步骤如下 1)对细胞运动图像进行预处理，利用分水岭方法完成目标细胞与背景以及目标细
胞之间的二值分割，并从获得的目标二值图像中提取目标细胞； 2)根据目标细胞运动特点，为每个目标细胞建立匀速运动学动力模型，并对模型 中的状态向量及误差协方差矩阵进行初始化； 3)利用卡尔曼滤波对当前时刻所有目标细胞在下一时刻的位置进行预测，以这些 预测值作为初始点，利用KHM算法对下一时刻二值图像进行聚类，获得图像中所有细胞中 心位置；
4)以最近邻方法为准则进行目标细胞关联，当关联过程中出现一对多，多对一的 特殊情况时，计算不同关联结果的差异代价函数，以差异代价函数最小的关联结果作为最 终关联结果，确定每个目标细胞在下一时刻的位置，完成运动细胞当前时刻到下一时刻的 足艮足宗； 5)更新步骤2)中的目标细胞模型的状态向量及误差协方差矩阵。
6)循环步骤3)到步骤5)直到跟踪到最后一帧图像。
本发明与现有的技术相比具有以下优点 1.由于KHM算法具有对初始值不敏感的特性和聚类性能良好的优点，本发明利用 它对细胞二值图像进行聚类，能够较为准确的得到细胞中心位置； 2.本发明由于引入差异代价函数，可以在出现细胞的分裂与结合、细胞聚合和多
细胞运动轨迹的密集交叉等特殊情况时，能够准确的完成目标细胞的跟踪； 3.仿真结果表明，本发明的目标细胞跟踪位置误差在很小的范围内，有较好的跟
踪效果。
图1是本发明的流程示意图； 图2是本发明使用的细胞视频序列原图； 图3是本发明实现的图像预处理后得到的细胞二值图像； 图4是本发明实现的部分细胞跟踪效果示意图； 图5是本发明实现的部分细胞跟踪位置误差示意图。
具体实施例方式参照图l，本发明的具体实现过程如下 步骤l，对细胞运动图像进行预处理，利用分水岭方法完成细胞目标与背景以及细 胞目标之间的二值分割，并从获得的二值图像中提取目标细胞。 la)采用0tsu方法计算图像阈值，对输入图像进行阈值处理，得到初始的目标与 背景二值分割结果； lb)将检测到的目标轮廓作为分水岭的内部控制标记。获取内部控制标记的方法 是，对于阈值化方法检测到的初始目标轮廓，将轮廓所在像素点灰度值赋为l，将轮廓以外 像素点灰度值赋为O，得到的二值图像，将该二值图像作为分水岭变换的内部控制标记；
lc)将内部控制标记依次进行距离变换和分水岭变换，并将得到的分水线结果作 为分水岭的外部控制标记；利用内部控制标记和外部控制标记，对输入图像的梯度图像进 行修正； ld)将修正后的梯度图像进行分水岭变换，检测图像中的目标轮廓，提取目标细胞。 步骤2，根据目标细胞运动特点，为每个目标细胞建立匀速运动学动力模型，并对 模型中的状态向量及其误差协方差矩阵进行初始化； 本模型是建立在一个二维平面上，以X轴表示水平方向，以Y轴表示垂直方向，以 图像左下角为原点来建立坐标系。由于相邻两幅图像的时间间隔较短，近似认为目标在相
5'卩两幅图像之间是匀速运动，状态向量定义为X = (x y Xvy》T，量测向量定义为Z = (x y)
T， x， y分别为目标中心横坐标和纵坐标，xv， yv分别为目标在x轴方向和y轴方向的速度， 状态转移矩阵和量测矩阵分别定义为<formula>formula see original document page 6</formula>
坐标，xv 协方差矩阵初始化为
其中，T为帧间时间间隔。 2a)建立目标运动模型为 Xjk+1) = FiXi(k) + r，i(k)
其中，Xi(k)是第i个目标的状态向量，&和rt分别是第i个目标的状态转移矩 阵和噪声矩阵W(k)是与第i个目标关联的系统噪声；
2b)建立测量模型为
Zi (k) = HA (k) +Vi (k) 其中，Zi(k)是第i个目标的量测向量，Hi是第i个目标的量测矩阵，Vi(k)是与第 i个目标关联的量测噪声。 2c)状态向量及其误差协方差矩阵进行初始化；
状态向量X二 (x y xv y》T中，x， y初始化为目标细胞在第一帧中的横坐标和纵 yv初始化为0 ;
「0.22 0 0 0 0 0.22 0 0 0 0 0.22 0 0 0 0 0.22 步骤3，利用卡尔曼滤波对当前时刻所有目标细胞在下一时刻的位置进行预测，以 这些预测值作为初始点，利用KHM算法对下一时刻二值图像进行聚类，获得图像中所有细 胞中心位置； 3a)利用卡尔曼滤波器进行状态估计和预测， 一步预测值J^H和预测误差的协方 差阵^w分别为
4—》维—Z+rgr 其中，,表示在k-1时刻对目标的状态估计，i表示k-1时刻误差协方差矩 阵，F为状态转移矩阵，r为噪声矩阵，Q为系统噪声W的方差； 3b)以预测位置为初始聚类中心，利用KHM方法对细胞二值图像进行聚类，获得细 胞中心位置； KHM方法是1999年由Zhang提出的一种基于中心迭代的聚类算法，用所有数据点 到每个聚类中心的调和平均值的和作为目标函数，具有对初始值不敏感的特性和聚类性能 良好的优点，其要点是
假设是一组数据元素X二 {XpX^…X山其中Xi二 {xn，xi2，... ，xim}表示具有m 个属性的数据对象，聚类个数为k， KHM算法的目标函数表示为
欧式
<formula>formula see original document page 7</formula>
这里Q二 {Cll， c12， . . . ， ClJ是聚类中心，d(Xi，C》是距离测度，简写为du，采用
距离，聚类中心的更新公式如下<formula>formula see original document page 7</formula> <formula>formula see original document page 7</formula> 通过对聚类中心的不断更新，使得目标函数不断减小直到稳定，获得最终的聚类结果。 相对于其它基于中心点迭代的聚类算法，KHM算法用数据点与所有聚类中心的距
离的调和平均替代了数据点与聚类中心的最小距离，从而引入了聚类中心对数据点的条件
概率和每次迭代过程中数据点的动态加权，实际上起到了将"硬"聚类"软化"的作用。 本发明将KHM算法应用在对细胞二值图像的聚类中，具体步骤为 3bl)将卡尔曼滤波器所得的一步量测预测值作为初始聚类中心，将细胞二值图像
中非零像素点作为待聚类的数据点； 3b2)计算每个聚类中心对二值图像中非零像素点的条件概率PKHM，以及迭代权函
<formula>formula see original document page 7</formula> <formula>formula see original document page 7</formula><formula>formula see original document page 7</formula> 该权函数W表示每次迭代过程中，数据点对聚类中心的影响程度，该条件概率PKHM
表示数据点归属于聚类中心所代表的类的可能性，其中Cj是聚类中心，Xi是二值图像中非
零像素点； 3b3)利用条件概率PKM，迭代权函数w，计算新的聚类中心<formula>formula see original document page 7</formula> 3b4)重复步骤3b2)和3b3)，直到收敛或达到最大迭代次数，停止迭代，得到聚类 中心，即细胞中心位置。
步骤4，进行目标细胞的关联，确定每个目标细胞在下一时刻的位置，完成运动细 胞当前时刻到下一时刻的跟踪； 4a)计算目标细胞一步量测预测值与聚类中心之间的马氏距离，利用最近邻关联 准则，用离目标细胞一步量测预测位置最近的聚类中心作为目标细胞的更新位置，完成目 标细胞的关联； 4b)若前后两帧图像出现细胞增加，细胞消失，细胞分裂等特殊情况时，前后两帧
图像中细胞数目会发生变化，出现一对多，多对一的特殊情况，通过计算不同关联结果的差
异代价函数，以差异代价函数最小的关联结果作为最终关联结果，确定目标细胞在当前时
刻的位置，差异代价函数可以表示为 cosTi, j = a D (i ， j) + P P (i ， j) + Y S (i ， j) 其中，D(i， j)为关联细胞与目标细胞之间马氏距离，P(i， j)为关联细胞与目标细 胞之间的周长差，S(i，j)为关联细胞与目标细胞之间面积差，差异代价函数越小，则这两个 细胞关联的可能性就越大，a为距离权值，13为周长差权值，Y为面积差权值，三个权值分 别取值为a = 0. 5， P = 0. 2， Y = 0. 3 ; 步骤5，更新步骤2中的目标细胞模型的状态向量及误差协方差矩阵。
利用更新公式更新状态向量》^和误差协方差矩阵声^ ，分别为
之『J^—i+《"-Zi^-,)
》举- - AC (《# + W《 其中，K为滤波增益，《=,T^(/^维—,/^ 1 ，式^为一步预测值，,为预测 误差的协方差阵，zk为新获得的量测值，即细胞中心位置，H为量测矩阵，R为量测噪声V的 协方差矩阵。 步骤6，循环步骤3到步骤5，不断的进行目标细胞从当前时刻到下一时刻的跟踪， 直到跟踪目标细胞到最后一帧图像，停止跟踪。
本发明的效果可以通过以下仿真进一步说明
1.仿真内容 采用蒙特卡罗仿真实验50次，对本发明进行运动细胞跟踪实验，其中图像序列来
自于显微细胞视频序列，如图2所示。图2中的图像序列信息如下 图像每帧大小650像素X515像素 图像序列数目74帧 图像帧间时间间隔T = 1/24s硬件平台为Intel Core2 Duo CPU E6550 @ 2. 33GHZ、2GB RAM
软件平台为MATLAB 7. 0
2.仿真结果及分析 利用分水岭方法对图2中显微细胞图像进行二值分割，得到的二值图像，如图3所 示。 利用本发明方法对目标细胞跟踪的效果，如图4所示，。 利用50次蒙特卡罗仿真实验的仿真结果求得目标细胞的位置跟踪误差，如图5所
从图4中可以看出，本发明方法的跟踪轨迹非常接近目标的真实运动轨迹，从图5中可以看出本发明方法的位置跟踪误差在很小的范围内。仿真结果表明，本发明能够对视频显微图像中的运动细胞进行自动跟踪，并取得较好的跟踪效果。
权利要求
一种基于KHM算法的视频显微图像细胞自动跟踪方法，包括如下步骤1)对细胞运动图像进行预处理，利用分水岭方法完成目标细胞与背景以及目标细胞之间的二值分割，并从获得的目标二值图像中提取目标细胞；2)根据目标细胞运动特点，为每个目标细胞建立匀速运动学动力模型，并对模型中的状态向量及误差协方差矩阵进行初始化；3)利用卡尔曼滤波对当前时刻所有目标细胞在下一时刻的位置进行预测，以这些预测值作为初始点，利用KHM算法对下一时刻二值图像进行聚类，获得图像中所有细胞中心位置；4)以最近邻方法为准则进行目标细胞关联，当关联过程中出现一对多，多对一的特殊情况时，计算不同关联结果的差异代价函数，以差异代价函数最小的关联结果作为最终关联结果，确定每个目标细胞在下一时刻的位置，完成运动细胞当前时刻到下一时刻的跟踪；5)更新步骤2)中的目标细胞模型的状态向量及误差协方差矩阵。6)循环步骤3)到步骤5)直到跟踪到最后一帧图像。
2. 根据权利要求1所述的细胞自动跟踪方法，其中步骤3)所述的用KHM算法对下一时 刻二值图像进行聚类，按如下步骤进行2a)将卡尔曼滤波器所得的一步量测预测值作为初始聚类中心，将细胞二值图像中非 零像素点作为待聚类的数据点；2b)计算每个聚类中心对二值图像中非零像素点的条件概率PKHM，以及迭代权函数w :(K)2&柳) = 7该权函数w表示每次迭代过程中，数据点对聚类中心的影响程度，该条件概率P 表示 数据点归属于聚类中心所代表的类的可能性，其中Cj是聚类中心，Xi是二值图像中非零像素点；2C)利用条件概率P^，迭代权函数W，计算新的聚类中心<formula>formula see original document page 2</formula>2d)重复步骤2b)和2c)，直到收敛或达到最大迭代次数，停止迭代，得到聚类中心，即 细胞中心位置。
3.根据权利要求1所述的细胞自动跟踪方法，其中步骤4)所述的以最近邻方法为准则 进行目标细胞关联，按如下步骤进行3a)计算目标细胞一步量测预测值与聚类中心之间的马氏距离，利用最近邻关联准则， 用离目标细胞一步量测预测位置最近的聚类中心作为目标细胞的更新位置，完成目标细胞的关联；3b)若前后两帧图像出现细胞增加，细胞消失，细胞分裂等情况时，前后两帧图像中细 胞数目会发生变化，出现一对多，多对一的特殊情况，通过计算不同关联结果的差异代价函 数，以差异代价函数最小的关联结果作为最终关联结果，确定目标细胞在当前时刻的位置， 差异代价函数可以表示为cosTi,j = aD(i， j) + PP(i， j) + YS(i， j)其中，D(i， j)为关联细胞与目标细胞之间马氏距离，P(i， j)为关联细胞与目标细胞之 间的周长差，S(i，j)为关联细胞与目标细胞之间面积差，差异代价函数越小，则这两个细胞 关联的可能性就越大，a为距离权值，13为周长差权值，Y为面积差权值，三个权值分别取 值为a = 0. 5， P = 0. 2， Y = 0. 3。
全文摘要
本发明公开了一种基于KHM算法的视频显微图像细胞自动跟踪方法，主要解决现有跟踪方法自动化程度较低以及出现特殊情况时无法准确跟踪的不足。其步骤包括1)利用分水岭方法完成细胞图像的二值分割，提取目标细胞；2)为每个目标细胞建立匀速运动学动力模型，并初始化；3)利用KHM算法对下一时刻二值图像进行聚类，获得所有细胞中心位置；4)通过目标细胞的关联，确定目标细胞在下一时刻的位置，完成运动细胞的跟踪；5)更新步骤2)中目标细胞模型的状态向量及其误差协方差矩阵；6)循环步骤3)到步骤5)直到跟踪到最后一帧图像。本发明具有确定细胞中心位置准确；在特殊情况下跟踪效果好的优点，可用于医学图像处理。
文档编号G06T7/20GK101719278SQ20091025442
公开日2010年6月2日 申请日期2009年12月21日 优先权日2009年12月21日
发明者刘芳, 吴建设, 尚荣华, 李阳阳, 王爽, 郑堃, 韩红, 马文萍 申请人:西安电子科技大学