CN110598830A - 基于标签蚁群的联合多细胞跟踪方法 - Google Patents

Abstract

本发明公开了一种基于标签蚁群的联合多细胞跟踪方法。首先，采用无标签的侦察蚁群在当前帧生成一组细胞候选蚁群，在此基础上构建一个二分图，由上帧估计出的细胞状态与当前产生的细胞候选蚁群组成；其次，利用标签蚁群并以航迹最优次模式分配为目标函数最优实现帧间的匹配，多细胞状态由进化的标签蚁群所近似的多贝努利参数得到，而细胞的谱系树由其确定的帧间匹配上的轨迹信息素场提取得到；最后，针对航迹断裂，提出了一种四步骤航迹恢复方法，以实现断裂航迹间的关联。本发明可以利用轨迹信息素来提取细胞的谱系数，用食物信息素来估计细胞的状态；与现有技术相比，细胞***精确率、召回率要明显提高。

Description

基于标签蚁群的联合多细胞跟踪方法

技术领域

本发明属于自动跟踪领域，更具体的涉及一种基于标签蚁群的联合多细胞跟 踪方法。

背景技术

对于细胞行为的定量分析需要可靠的自动跟踪技术，且能自动估计出与每个 细胞相关的状态、形态学参数及对应的谱系树。人工跟踪对实验员的专业水平要 求高，且耗时、结果难复制、并行跟踪难、易出错等等，所以一体化的自动跟踪 技术是细胞跟踪近年来主要的研究重点。

目前，多细胞自动跟踪方法有多种，如：分割关联方法，模型进化方法，滤 波方法。分割关联方法是目前常用的解决低信噪比多细胞***跟踪问题的有效途 径，但分割的精度、所采用的数据关联方法对于跟踪精度至关重要。模型跟踪方 法能够跟踪细胞的变形并自动产生分割，此方法容易受到近邻情形影响，如细胞 重叠、聚蚁群，粘连等。滤波方法通常归为贝叶斯框架，如多假设方法，粒子滤 波器，联合的概率数据关联算法和基于随机有限蚁群滤波方法，这些方法最明显 的特点是依赖严谨的数学理论基础，就细胞跟踪领域，对于细胞***、新进入细 胞，细胞身份等识别与估计的研究还不够深。

在细胞***、变形和测量不确定性情形下，数据关联一直是该领域挑战性难 题，与之紧密相关的是分割和***事件的检测。目前的分割方法主要有阈值方法， 水平蚁群方法，图形模型法和一些专门开发的细胞跟踪软件蚁群成器。最近，机 器学***。对于细胞跟踪，常用的数 据关联方法有整数规划法，最小成本流等算法，但这些算法的性能在很大程度上 取决于细胞分割的精度，且把数据关联与滤波算法独立处理。

本发明提出了一种基于标签蚁群的多细胞数据关联与状态估计的联合跟踪 方法。本发明首先采用无标签的侦察蚁群在当前帧生成一组细胞候选蚁群，在此 基础上构建一个二分图，它由上帧估计出的细胞状态与当前产生的细胞候选蚁群 组成，利用标签蚁群决策行为并以航迹最优次模式分配为目标函数最优实现帧间 的匹配，最后，多细胞状态由进化的标签蚁群所近似的多贝努利参数得到，而细 胞的谱系树由其确定的帧间匹配上的轨迹信息素场提取得到。针对航迹断裂，提 出了一种四步骤航迹恢复方法，以实现断裂航迹间的关联。实验结果说明，本发 明所提出的方法能很好地跟踪细胞***、变形且测量不确定性等跟踪难题，与目 前常用的多细胞跟踪方法相比，跟踪性能大幅提高。

蚁群：表示一个由多个蚂蚁所构成的蚁群合，其中每个蚂蚁的决策行为是随 机的，且行为简单，但蚂蚁个体间相互协作，可共同完成某个复杂任务。

标签：细胞标签用以标识细胞的身份和唯一性，用l＝(k₁,i₁)来表示，其中，k₁表示该细胞首次出现的时刻，i₁表示在该时刻出现的第i₁个细胞。当该细胞在时 刻k₂发生***时，其标签用(j＝1,2)来表示。为便于在图中标 记，通常将标记成

发明内容

1、发明目的。

本发明为了解决细胞***、变形和测量不确定性情形下等跟踪难题，提出了 一种基于标签蚁群的多细胞数据关联与状态估计的联合跟踪方法。首先，采用无 标签的侦察蚁群在当前帧生成一组细胞候选蚁群，在此基础上构建一个二分图， 它由上帧估计出的细胞状态与当前产生的细胞候选蚁群组成，利用标签蚁群并以 航迹最优次模式分配为目标函数最优实现帧间的匹配，最后，多细胞状态由进化 的标签蚁群所近似的多贝努利参数得到，而细胞的谱系树由其确定的帧间匹配上 的轨迹信息素场提取得到。最后，针对航迹断裂，提出了一种四步骤航迹恢复方 法，以实现断裂航迹间的关联。

2、本发明所采用的技术方案。

本发明公开了一种基于标签蚁群的联合多细胞跟踪方法，按照如下步骤进行：

多贝努利随机有限集定义及其近似：

多贝努利随机有限集：多贝努利随机有限集X是由一组固定数目M个的单 贝努利随机有限集X^(ζ)(ζ＝1,...,M)组成的并，而每个单贝努利随机有限集X^(ζ)可 以用其存在概率r^(ζ)∈(0,1)及其对应的概率密度函数π^(ζ)来表示，于是，多贝努利 随机有限集X可用贝努利参数集来表示；标签多贝努利随机有限集 可用{(r^(ζ),π^(ζ)):ζ∈Ψ}表示，其中Ψ表示索引集合，若单个贝努利分量返回一个非 空集，则赋予其标签；

标签多贝努利随机有限集近似：标签多贝努利随机有限集可用标签多贝努利 蚁群来近似，即其中，Ψ_a表示蚁群索引集合，表示标签蚁 群的存在概率，表示标签蚁群的离散概率密度分布，若大于某阈 值，则认为该蚁群会生存下来，并与细胞标签关联，其中定义为从蚁群索 引集合到细胞标签的身份映射；给定一个标签蚁群其存在概率可计算 为

其中，s为面积积分变量，g(s)为启发式似然函数，τ^(ζ)为标签蚁群对应的信 息素场，R^(ζ)为标签蚁群的影响区域；

细胞帧间匹配问题定义：

表示在第k帧估计得到的m个细胞位置，表示在第 k+1帧得到的n个细胞候选集，于是构建一个二分图G＝(U,V,E)，其中， E表示从U到V上的边匹配假设的集合，帧间匹配的问题就转化成

且其中，D^(C)(·)表示对于匹配 假设C前提下的成本函数，f(·)表示从U到V上的匹配假设映射，f^-1(·)表示f的 逆映射，|·|表示某集合的势；

成本函数D采用的是航迹最优次模式分配(OSPA-T)，即：给定两组带标签 的集合X和Y，元素个数分别为m和n，X＝{(x₁,l₁),...,(x_m,l_m)}表示上一帧已估计的 带标签l的细胞状态x的集合，Y＝{(y₁,s₁),...,(y_n,s_n)}表示在当前帧与之对应的匹配 假设带标签s的候选量测y的子集，若m≤n,则OSPA-T定义为

其中，d^c((x,l),(y,s))＝min(c,d((x,l),(y,s)))表示截距距离，参数c和p分别表示对 两集合元素个数差异和带标签状态变量差异的敏感度正系数(均大于零)，ψ_n表 示在集合Y上的可行性扰动集合，ψ(i)表示ψ_n中的一个扰动，基距离 d((x,l),(y,s))＝(d(x,y)^p′+d(l,s)^p′)^1/p′，p′表示基距离阶参数，d(x,y)表示两个状态x,y 之间的p′范数距离，d(l,s)表示标签误差，若l＝s，则d(l,s)＝0，否则 d(l,s)＝λ,λ∈[0,c]；若m＞n，有D_p,c(X,Y)＝D_p,c(Y,X)；

基于无标签侦察蚁群生成候选集

在第k+1帧，给定侦察蚁群初始随机分布在当前的图像区域，为对侦察蚂蚁 行为进行建模，引入了一种新的蚂蚁动力学***，即从混沌行为逐步演化至确定 性行为，从而发现食物源可能出现的区域；侦察蚂蚁的动力学行为定义如下

其中，表示侦察蚂蚁在第t次迭代时混沌移动变量，h_a表示蚂蚁自组织 因子，取值范围在0和1之间，表示侦察蚂蚁在第t次迭代时的状态，参数ψ_d表示蚂蚁在搜索空间的搜索范围，参数V_a表示控制搜索空间的不同区域，a、b为 常数，p_a(t-1)表示的是该侦察蚂蚁及其近邻蚂蚁在前t-1次迭代里得到的最优解；

表示侦察蚂蚁在第t次迭代时的近邻蚂蚁集合，假定侦察蚂蚁及其近 邻蚂蚁在当前的状态、前t-1次迭代最优状态已知，则侦察蚂蚁在前t次迭代里 得到的最优状态为

其中，符号∨表示取大运算，g^-1(·)表示对似然函数g求逆运算，似然函数g的 定义为

其中，I^l表示归一化的灰度图像在像素l的强度值，B_q表示以像素q为中心 的方形区域所包括的像素集合，N_(q)表示像素q的近邻像素集合，符号|·|表示集 合里元素的数目，符号∧表示求小运算；在第t次迭代上，若侦察蚂蚁在像素q 上，则有

对于先后两次迭代中蚁群集合的蚂蚁状态视为两状态集合和则将OSPA-T简化为最优次模式分配(OSPA)距离于是， 当m_a≤n_a时，有

若m_a＞n_a，定义

步骤1：初始化：给定若干组侦察蚁群，侦察蚂蚁的动力学行为中的混沌移 动变量y_a(0)，自组织因子h_a，侦察蚂蚁初始状态搜索范围参数ψ_d，搜索 空间的不同区域参数V_a，常数a、b等已知；

步骤2：对于每一组侦察蚁群，每个侦察蚂蚁按动力学行为模型移动，按上 述似然函数定义来感知侦察蚂蚁所在相应像素的似然函数值，并与其近邻侦察蚂 蚁进行信息交换得到截止目前迭代次得到的最优状态；

步骤3：若侦察蚂蚁所在的像素位置是已探知可能细胞所在位置区域，则将 其对应的似然函数设为一很小值；

步骤4：计算本次迭代侦察蚁群状态集合与上一次迭代侦察蚁群状态集合的 最优次模式分配距离的差值，若差值小于某一个阈值，则终止该侦察蚁群的 侦察工作，进而通过聚类重构，生成若干个侦察子蚁群；

步骤5：重复步骤2-4，直到没有新的侦察子蚁群生成为止；将已生成的每 个侦察子蚁群的中心作为这样，共生成n个候选量测

更进一步，基于标签蚁群的匹配构建与参数估计联合跟踪具体方法如下：

将标签多贝努利随机有限集参数集{(r^(ζ),π^(ζ)):ζ∈Ψ}与标签多贝努利蚁群参数集存在一一对应关系，且标签贝努利蚁群与细胞有一一对应 关系，即

假定标签为ζ蚁群以一个生存概率会从上一帧位置移到下一帧所在位置，除 了其标签的唯一性不变以外，生存概率定义为：

其中，H，W表示当前帧图像的高与宽，表示在第k帧标签蚁群 中心离图像边界最短的高度、宽度方向距离，表示该标签蚁群在第k帧的存在 概率；

同样，标签蚁群在下一帧k+1发生***的概率定义为

其中，exp(·)表示指数函数，分别表示细胞在第k帧和初始 帧估计所得到的轮廓宽高比，分别表示细胞在第k帧和初始帧 的平均强度，υ为0和1之间的调节系数，表示对应细胞的局部窗口面积大小，表示在第k帧所估计出的细胞面积大小；

基于上述定义，基于标签蚁群的帧间匹配与参数估计一体化方法描述如下：

步骤1：初始化：假定在上一帧k共有m个细胞估计，则有 为 细胞i的估计位置，其所对应的标签多贝努利蚁群参数集可表示为对于每个贝努利蚁群，蚁群规模均为由无标签侦察蚁群生成n个候选量测集 表示候选量测j位置；

步骤2：对于任一标签贝努利蚁群任选其中一标签蚂蚁(与所在的蚁 群标签一致)，随机生成的两随机数κ_S∈(0,1)和κ_M∈(0,1)，若满足：

(1)：则该标签蚂蚁能从上一帧生存下来，并从其所 对应的细胞估计位置出发，选择候选量测的概率为

其中，1_Y(X)表示包含函数，若则等于1，否则为0；表示第 个细胞估计位置所对应的候选量测子集，表示在轨迹上的轨迹信息 素强度，θ(·)表示轨迹信息素的身份特征函数，它与细胞相关；表示启 发式函数，并有d(·,·)表示距离函数，α₁,β₁表示相对重要性系数；

(2)：说明该蚂蚁不仅能从上一帧生存下来，而且还 要应对细胞发生***可能，于是，该标签蚂蚁复制成两只一样的标签蚂蚁，先后 从其所对应的细胞估计位置出发，选择不同候选量测的概率为

(3)：或该蚂蚁不动，说明在下一帧它所对应的细胞 可能消失；

当标签蚂蚁按照(1)或(2)决策后，则在经过的轨迹上释放一定量的 信息素：

其中，ρ表示信息素残留系数，表示在当前帧k+1对应标签贝努利蚁群 在上一次迭代更新后所得到的存在概率，表示细胞的影响区域；

步骤3：当m个标签贝努利蚁群按照步骤2完成相应的决策后，就构成了一 个从U到V之间的一个匹配假设C，并计算所对应的OSPA-T值；

步骤4：按步骤2-3，进行次循环，在次循环中，得到本次迭代最佳匹 配假设C，其对应的候选集表示为对任按 步骤2所得到的标签贝努利子蚁群开始细胞状态估计，即觅食；若觅食蚂蚁 在像素p上，它移动到像素q的概率定义为

其中，α₂,β₂表示相对重要性系数，N_(p)表示像素p近邻像素集合，表 示在当前帧k+1像素q上的食物信息素强度；伴随着上述的决策，标签蚂蚁同时 在像素q上释放一定量的信息素

其中，ε₁为阈值系数；

当标签贝努利子蚁群中所有蚂蚁完成上述觅食决策后，对应的食物信息 素场进行更新

其中，表示在像素q上总的信息素输入，表示在像素q上总 的信息素扩散输入，定义为

其中，D表示扩散常系数，取值范围在(0,1)之间，κ为控制常系数，表 示像素q和之间的距离，R_q表示像素q的影响区域；

利用多次迭代后，由所得到的食物信息素场和启发式函数，得到标签贝努利 子蚁群的贝努利参数近似

其中，表示当前帧k+1觅食标签蚁群的概率密度离散分布， 表示蚂蚁的状态所对应的权重，表示标签贝努利子蚁群中的 蚂蚁数；对于食物信息素场，进行二值化、形态学和填充三步简单处理，即可得 到候选集对应细胞的轮廓及其位置

步骤5：对于截止目前迭代所得到的最佳匹配假设C，对其对应的候选集位 置，轨迹信息素场，食物信息素场进行全局更新

其中，Q₁，Q₂为信息素加强常数，表示截止目前迭代所得到的最佳匹配假 设C对应的位置更新后的候选集集合；

步骤6：重复步骤2-5，直至收敛条件满足；同时，对所得到的最佳匹配迭 代所得到的食物信息素场，进行二值化、形态学和填充三步简单处理，若得到的 标签子蚁群存在概率可得到对应细胞的轮廓及其位置估计；

对于所得到的帧间最佳匹配，将其轨迹信息素场提取，并保留轨迹信息素的 身份特征，通过帧帧顺序匹配，形成细胞谱系树；

对于未匹配的候选量测，视为当前帧可能出现的新生项，对应生成相互独立 的标签贝努利子蚁群，按步骤4对其对应的食物信息素场构建，若新生的标签子 蚁群存在概率并不低于0.5，则提取出相应的细胞的轮廓及其位置估计；

对于多个标签贝努利子蚁群所得到的位置估计近邻且相互距离小于某阈值， 则将存在概率大的标签贝努利子蚁群保留，其他满足条件的标签贝努利子蚁群直 接去除并不再考虑；同时，对于标签更新，将最先出现的标签赋予存在概率大的 标签贝努利子蚁群。

更进一步，基于四阶段的航迹恢复具体方法如下：

基于四阶段的航迹恢复的相似度标准如下：

相似度标准1：χ⁺，χ^-分别表示航迹χ在诞生和结束时刻所对应的细胞状 态，对于两航迹χ₁和χ₂，若结束时刻先于则相似度标准1定义为

其中，R^-表示满足帧间最大移动距离约束条件的航迹集合，且其任一航迹χ^-结束时刻先于 表示对应的航迹χ₁在诞生时刻细胞的面积估计，表 示对应的航迹χ₂在结束时刻细胞的面积估计，∩表示两细胞面积交集，∪表示 两细胞面积求和；

相似度标准2：对于航迹χ_i，其在诞生时刻k所对应的细胞状态为对 于航迹在时刻k和k-1所对应的细胞状态分别为和则相似度标准 2定义为

其中，C_h(·)和C_w(·)分别表示细胞轮廓的主轴和短轴的长度，ΔC_h(·,·)和ΔC_w(·,·)分别 表示细胞轮廓的主轴和短轴的长度差绝对值，ΔA(·,·)表示细胞面积差的绝对值；

具体步骤如下：

阶段1：给定一长度为3个采样间隔的时间窗口，从最后一帧开始反向滑动， 当其右侧边与航迹的起始位置相交时，对应的帧若为第k帧，此相交的航迹 将与其它在这时间窗口内(从第k-3帧至第k-1帧)结束的航迹进行关联；若第 k-1帧内有相应的航迹候选(即该航迹在第k-1帧结束)，且满足帧间最大移动距 离约束条件，则按相似度标准1选择最大一个与之关联；否则，进入阶段2；

阶段2：排除阶段1已关联的若非空，且第k-2，k-3帧内有相应 的航迹候选(即该航迹在第k-2，k-3帧结束)，并满足帧间最大移动距离约束 条件，则按相似度标准1选择最大一个与之关联；

阶段3：排除阶段1，2已关联的若仍非空，与近邻(满足帧间最 大移动距离约束条件)的完整航迹在第k，k-1帧进行相似度标准2计算，若 相似度值大于某给定阈值，则航迹χ_i并入航迹对应细胞***；

阶段4：排除阶段1，2，3已关联的若仍非空，则认为该起始的航 迹(且航迹寿命至少为3个采样间隔)为正式航迹，否则，去除该航迹。

3、本发明所产生的技术效果。

(1)本发明提出了一种基于标签蚁群的多细胞数据关联与状态估计的联合跟 踪方法，数据关联与状态估计并行且交互式进行；解决了基于标签蚁群的匹配构 建问题，不同于一般的1-1匹配，本发明能处理富余的候选蚁群量测，并通过蚁 群觅食行为不断更新最佳匹配所对应的候选蚁群量测位置，即所匹配的候选蚁群 量测是动态变化的。

(2)本发明提出了一种新的匹配目标函数，即OSPA-T，综合考虑势，位置， 标签匹配误差。

(3)本发明可以利用轨迹信息素来提取细胞的谱系数，用食物信息素来估计 细胞的状态；与现有技术相比，细胞***精确率、召回率要明显提高。

附图说明

图1为多细胞跟踪序列(第一行为序列1，第二行为序列2)。

图2为序列1跟踪结果(细胞谱系树)，(a)三维谱系树视角1，(b)三维谱 系树视角2，(c)三维谱系树视角3，(d)三维谱系树视角4。

图3为序列2跟踪结果(细胞谱系树)，(a)三维谱系树视角1，(b)三维谱 系树视角2，(c)三维谱系树视角3，(d)三维谱系树视角4。

具体实施方式

实施例1

1.多贝努利随机有限蚁群定义及其近似：

多贝努利随机有限蚁群：多贝努利随机有限蚁群X是由一组固定数目(如M 个)的单贝努利随机有限蚁群X^(ζ)(ζ＝1,...,M)组成的并，而每个单贝努利随机有 限蚁群X^(ζ)可以用其存在概率r^(ζ)∈(0,1)及其对应的概率密度函数π^(ζ)来表示，于 是，多贝努利随机有限蚁群X可用贝努利参数蚁群来表示。标签多 贝努利随机有限蚁群可用{(r^(ζ),π^(ζ)):ζ∈Ψ}表示，其中Ψ表示索引蚁群合，若单个 贝努利分量返回一个非空蚁群，则赋予其标签。

标签多贝努利随机有限蚁群近似：标签多贝努利随机有限蚁群可用标签多贝 努利蚁群来近似，即其中，Ψ_a表示蚁群索引蚁群合，表示 标签蚁群的存在概率，表示标签蚁群的离散概率密度分布，若大 于某阈值，则认为该蚁群会生存下来，并与细胞标签关联，其中定义为从 蚁群索引蚁群合到细胞标签的身份映射。给定一个标签蚁群其存在概率可计算为

其中，s为面积积分变量，g(s)为启发式似然函数，τ^(ζ)为标签蚁群对应的信 息素场，R^(ζ)为标签蚁群的影响区域。

2.细胞帧间匹配问题定义

表示在第k帧估计得到的m个细胞位置，表示在第 k+1帧得到的n个细胞候选蚁群，于是构建一个二分图G＝(U,V,E)，其中， E表示从U到V上的边(匹配假设)的蚁群合，帧间匹配的问题就 转化成且 其中，D^(C)(·)表示对于匹配假设C前提下的成本函数，f(·) 表示从U到V上的匹配假设映射，f^-1(·)表示f的逆映射，|·|表示某蚁群合的势。

成本函数D采用的是航迹最优次模式分配(OSPA-T)，即：给定两组带标签 的蚁群合X和Y，元素个数分别为m和n，X＝{(x₁,l₁),...,(x_m,l_m)}表示上一帧已估计 的带标签l的细胞状态x的蚁群合，Y＝{(y₁,s₁),...,(y_n,s_n)}表示在当前帧与之对应的 匹配假设带标签s的候选量测y的子蚁群，若m≤n,则OSPA-T定义为

其中，d^c((x,l),(y,s))＝min(c,d((x,l),(y,s)))表示截距距离，参数c和p分别表示对 两蚁群合元素个数差异和带标签状态变量差异的敏感度正系数(均大于零)，ψ_n表示在蚁群合Y上的可行性扰动蚁群合，ψ(i)表示ψ_n中的一个扰动，基距离 d((x,l),(y,s))＝(d(x,y)^p′+d(l,s)^p′)^1/p′，p′表示基距离阶参数，d(x,y)表示两个状态x,y 之间的p′范数距离，d(l,s)表示标签误差，若l＝s，则d(l,s)＝0，否则 d(l,s)＝λ,λ∈[0,c]。若m＞n，有D_p,c(X,Y)＝D_p,c(Y,X)。

3.基于无标签侦察蚁群生成候选蚁群

在第k+1帧，给定侦察蚁群初始随机分布在当前的图像区域，为对侦察蚂蚁 行为进行建模，引入了一种新的蚂蚁动力学***，即从混沌行为逐步演化至确定 性行为，从而发现食物源(细胞)可能出现的区域。侦察蚂蚁的动力学行为定义 如下

其中，表示侦察蚂蚁在第t次迭代时混沌移动变量，h_a表示蚂蚁自组织 因子，取值范围在0和1之间，表示侦察蚂蚁在第t次迭代时的状态，参数ψ_d表示蚂蚁在搜索空间的搜索范围，参数V_a表示控制搜索空间的不同区域，a、b为 常数，p_a(t-1)表示的是该侦察蚂蚁及其近邻蚂蚁在前t-1次迭代里得到的最优解。

表示侦察蚂蚁在第t次迭代时的近邻蚂蚁蚁群合，假定侦察蚂蚁及 其近邻蚂蚁在当前的状态、前t-1次迭代最优状态已知，则侦察蚂蚁在前t次迭 代里得到的最优状态为

其中，符号∨表示取大运算，g^-1(·)表示对似然函数g求逆运算，似然函数g 的定义为

其中，I^l表示归一化的灰度图像在像素l的强度值，B_q表示以像素q为中心 的方形区域所包括的像素蚁群合，N_(q)表示像素q的近邻像素蚁群合，符号|·|表 示蚁群合里元素的数目，符号∧表示求小运算。在第t次迭代上，若侦察蚂蚁在 像素q上，则有

对于先后两次迭代中蚁群蚁群合的蚂蚁状态视为两状态蚁群合 和则将OSPA-T简化为最优次模式分配(OSPA) 距离于是，当m_a≤n_a时，有

若m_a＞n_a，定义

基于上述定义，对于第k+1帧细胞图像，多贝努利随机有限侦察蚁群的生成 方法如下：

步骤1：初始化：给定若干组侦察蚁群，侦察蚂蚁的动力学行为中的混沌移 动变量y_a(0)，自组织因子h_a，侦察蚂蚁初始状态搜索范围参数ψ_d，搜索 空间的不同区域参数V_a，常数a、b等已知。

步骤2：对于每一组侦察蚁群，每个侦察蚂蚁按动力学行为模型移动，按上 述似然函数定义来感知侦察蚂蚁所在相应像素的似然函数值，并与其近邻侦察蚂 蚁进行信息交换得到截止目前迭代次得到的最优状态。

步骤3：若侦察蚂蚁所在的像素位置是已探知可能细胞所在位置区域，则将 其对应的似然函数设为一很小值。

步骤4：计算本次迭代侦察蚁群状态蚁群合与上一次迭代侦察蚁群状态蚁群 合的最优次模式分配距离的差值，若差值小于某一个阈值，则终止该侦察蚁 群的侦察工作，进而通过聚类重构，生成若干个侦察子蚁群。

步骤5：重复步骤2-4，直到没有新的侦察子蚁群生成为止。将已生成的每 个侦察子蚁群的中心作为这样，共生成n个候选量测

4.基于标签蚁群的匹配构建与参数估计联合跟踪

本发明将标签多贝努利随机有限蚁群参数蚁群{(r^(ζ),π^(ζ)):ζ∈Ψ}与标签多贝努利蚁群参数蚁群存在一一对应关系，且标签贝努利蚁群与细 胞有一一对应关系，即

假定标签为ζ蚁群以一个生存概率会从上一帧位置移到下一帧所在位置，除 了其标签的唯一性不变以外，生存概率定义为：

其中，H，W表示当前帧图像的高与宽，表示在第k帧标签蚁群 中心离图像边界最短的高度、宽度方向距离，表示该标签蚁群在第k帧的存在 概率。

同样，标签蚁群在下一帧k+1发生***的概率定义为

其中，exp(·)表示指数函数，分别表示细胞在第k帧和初始 帧估计所得到的轮廓宽高比，分别表示细胞在第k帧和初始帧 的平均强度，υ为0和1之间的调节系数，表示对应细胞的局部窗口面积大小，表示在第k帧所估计出的细胞面积大小。

基于上述定义，基于标签蚁群的帧间匹配与参数估计一体化方法描述如下：

步骤1：初始化：假定在上一帧k共有m个细胞估计，则有 为 细胞i的估计位置，其所对应的标签多贝努利蚁群参数蚁群可表示为对于每个贝努利蚁群，蚁群规模均为由无标签侦察蚁群生成n个候选量测蚁 群 表示候选量测j位置。

步骤2：对于任一标签贝努利蚁群任选其中一标签蚂蚁(与所在的蚁 群标签一致)，随机生成的两随机数κ_S∈(0,1)和κ_M∈(0,1)，若满足：

(1)：则该标签蚂蚁能从上一帧生存下来，并从其所 对应的细胞估计位置出发，选择候选量测的概率为

其中，1_Y(X)表示包含函数，若则等于1，否则为0；表示第 个细胞估计位置所对应的候选量测子蚁群，表示在轨迹上的轨迹信 息素强度，θ(·)表示轨迹信息素的身份特征函数，它与细胞相关；表示 启发式函数，并有d(·,·)表示距离函数，α₁,β₁表示相对重要性系 数。

(2)：说明该蚂蚁不仅能从上一帧生存下来，而且还 要应对细胞发生***可能，于是，该标签蚂蚁复制成两只一样的标签蚂蚁，先后 从其所对应的细胞估计位置出发，选择不同候选量测的概率为

(3)：或该蚂蚁不动，说明在下一帧它所对应的细胞 可能消失。

当标签蚂蚁按照(1)或(2)决策后，则在经过的轨迹上释放一定量的 信息素：

其中，ρ表示信息素残留系数，表示在当前帧k+1对应标签贝努利蚁群 在上一次迭代更新后所得到的存在概率，表示细胞的影响区域。

步骤3：当m个标签贝努利蚁群按照步骤2完成相应的决策后，就构成了一 个从U到V之间的一个匹配假设C，并计算所对应的OSPA-T值。

步骤4：按步骤2-3，进行次循环，在次循环中，得到本次迭代最佳匹 配假设C，其对应的候选蚁群表示为对任按步骤2所得到的标签贝努利子蚁群开始细胞状态估计，即觅食。若觅食蚂 蚁在像素p上，它移动到像素q的概率定义为

其中，α₂,β₂表示相对重要性系数，N_(p)表示像素p近邻像素蚁群合，表示在当前帧k+1像素q上的食物信息素强度。伴随着上述的决策，标签蚂蚁同 时在像素q上释放一定量的信息素

其中，ε₁为阈值系数。

当标签贝努利子蚁群中所有蚂蚁完成上述觅食决策后，对应的食物信息 素场进行更新

其中，表示在像素q上总的信息素输入，表示在像素q上总 的信息素扩散输入，定义为

其中，D表示扩散常系数，取值范围在(0,1)之间，κ为控制常系数，表 示像素q和之间的距离，R_q表示像素q的影响区域。

利用多次迭代后，由所得到的食物信息素场和启发式函数，得到标签贝努利 子蚁群的贝努利参数近似

其中，表示当前帧k+1觅食标签蚁群的概率密度离散分布， 表示蚂蚁的状态所对应的权重，表示标签贝努利子蚁群中的 蚂蚁数。对于食物信息素场，进行二值化、形态学和填充三步简单处理，即可得 到候选蚁群对应细胞的轮廓及其位置

步骤5：对于截止目前迭代所得到的最佳匹配假设C，对其对应的候选蚁群 位置，轨迹信息素场，食物信息素场进行全局更新

其中，Q₁，Q₂为信息素加强常数，表示截止目前迭代所得到的最佳匹配假 设C对应的位置更新后的候选蚁群蚁群合。

步骤6：重复步骤2-5，直至收敛条件满足。同时，对所得到的最佳匹配迭 代所得到的食物信息素场，进行二值化、形态学和填充三步简单处理，若得到的 标签子蚁群存在概率可得到对应细胞的轮廓及其位置估计。

对于所得到的帧间最佳匹配，将其轨迹信息素场提取，并保留轨迹信息素的 身份特征，通过帧帧顺序匹配，形成细胞谱系树。

对于未匹配的候选量测，视为当前帧可能出现的新生项，对应生成相互独立 的标签贝努利子蚁群，按步骤4对其对应的食物信息素场构建，若新生的标签子 蚁群存在概率并不低于0.5，则提取出相应的细胞的轮廓及其位置估计。

对于多个标签贝努利子蚁群所得到的位置估计近邻且相互距离小于某阈值， 则将存在概率大的标签贝努利子蚁群保留，其他满足条件的标签贝努利子蚁群直 接去除并不再考虑。同时，对于标签更新，将最先出现的标签赋予存在概率大的 标签贝努利子蚁群。

5.基于四步骤的航迹恢复

针对上述方法生成的航迹，可能出现航迹的断裂现象，为此，提出一种反向 滑动时间窗口航迹恢复四步骤法。所采用的关联标准，除了满足帧间最大移动距 离约束条件外，定义如下相似度标准：

相似度标准1：χ⁺，χ^-分别表示航迹χ在诞生和结束时刻所对应的细胞状 态，对于两航迹χ₁和χ₂，若结束时刻先于则相似度标准1定义为

其中，R^-表示满足帧间最大移动距离约束条件的航迹蚁群合，且其任一航 迹χ^-结束时刻先于 表示对应的航迹χ₁在诞生时刻细胞的面积估计， 表示对应的航迹χ₂在结束时刻细胞的面积估计，∩表示两细胞面积交蚁群， ∪表示两细胞面积求和。

相似度标准2：对于航迹χ_i，其在诞生时刻k所对应的细胞状态为对 于航迹在时刻k和k-1所对应的细胞状态分别为和则相似度标准 2定义为

其中，C_h(·)和C_w(·)分别表示细胞轮廓的主轴和短轴的长度，ΔC_h(·,·)和ΔC_w(·,·)分别 表示细胞轮廓的主轴和短轴的长度差绝对值，ΔA(·,·)表示细胞面积差的绝对值。

具体步骤如下：

步骤1：给定一长度为3个采样间隔的时间窗口，从最后一帧开始反向滑动， 当其右侧边与航迹的起始位置相交时，对应的帧若为第k帧，此相交的航迹 将与其它在这时间窗口内(从第k-3帧至第k-1帧)结束的航迹进行关联。若第k-1帧内有相应的航迹候选(即该航迹在第k-1帧结束)，且满足帧间最大移动距 离约束条件，则按相似度标准1选择最大一个与之关联。否则，进入步骤2。

步骤2：排除步骤1已关联的若非空，且第k-2，k-3帧内有相应 的航迹候选(即该航迹在第k-2，k-3帧结束)，并满足帧间最大移动距离约束 条件，则按相似度标准1选择最大一个与之关联。

步骤3：排除步骤1，2已关联的若仍非空，与近邻(满足帧间最 大移动距离约束条件)的完整航迹在第k，k-1帧进行相似度标准2计算，若 相似度值大于某给定阈值，则航迹χ_i并入航迹对应细胞***。

步骤4：排除步骤1，2，3已关联的若仍非空，则认为该起始的航 迹(且航迹寿命至少为3个采样间隔)为正式航迹，否则，去除该航迹。

采用本发明提供的技术，具有以下特点：

1)提出了一种基于标签蚁群的多细胞数据关联与状态估计的联合跟踪方法， 数据关联与状态估计并行且交互式进行。

2)本发明所涉及的基于标签蚁群的匹配构建难题，不同于一般的1-1匹配问 题，本发明能处理富余的候选蚁群量测，并通过蚁群觅食行为不断更新最佳匹配 所对应的候选蚁群量测位置，即所匹配的候选蚁群量测是动态变化的。

3)本发明提出了一种新的匹配目标函数，即OSPA-T，综合考虑势，位置， 标签匹配误差。

4)本发明可以利用轨迹信息素来提取细胞的谱系数，用食物信息素来估计细 胞的状态。

5)对比专利1ZL201610457030.4公开了一种基于带标签的多贝努利滤波器 的多细胞自动跟踪方法；对比专利2ZL201810585820.X公开了多贝努利随机有 限蚁群多细胞跟踪方法；对比专利3ZL201810815283.3公开了双层多贝努利随 机有限蚁群多细胞跟踪方法；本发明所设计的是基于蚁群启发的数据关联与状态 估计一体化的多细胞跟踪方法。本发明与所比对的专利1，2，3不同的是：将数 据关联与参数估计进行联合处理，所采用的理论、技术、方法完全不同。由下表 1所知，本发明与比对的专利1，2，3相比，细胞***精确率、召回率要明显提 高。

表1***事件性能比较(基于10次仿真结果统计)

评价指标：

其中，TP表示检测到正确细胞***事件的数目，FP表示误检测到细胞*** 事件的数目(虚警)，FN表示未检测到细胞***事件的数目。

上述实施例为本发明较佳的实施方式，但本发明的实施方式并不受上述实施 例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替 代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (3)

1.基于标签蚁群的联合多细胞跟踪方法，其特征在于：

多贝努利随机有限集定义及其近似：

多贝努利随机有限集：多贝努利随机有限集X是由一组固定数目M个的单贝努利随机有限集X^(ζ)(ζ＝1,...,M)组成的并，而每个单贝努利随机有限集X^(ζ)可以用其存在概率r^(ζ)∈(0,1)及其对应的概率密度函数π^(ζ)来表示，于是，多贝努利随机有限集X可用贝努利参数集来表示；标签多贝努利随机有限集可用{(r^(ζ),π^(ζ)):ζ∈Ψ}表示，其中Ψ表示索引集合，若单个贝努利分量返回一个非空集，则赋予其标签；

标签多贝努利随机有限集近似：标签多贝努利随机有限集可用标签多贝努利蚁群来近似，即其中，Ψ_a表示蚁群索引集合，表示标签蚁群的存在概率，表示标签蚁群的离散概率密度分布，若大于某阈值，则认为该蚁群会生存下来，并与细胞标签关联，其中定义为从蚁群索引集合到细胞标签的身份映射；给定一个标签蚁群其存在概率可计算为

其中，s为面积积分变量，g(s)为启发式似然函数，τ^(ζ)为标签蚁群对应的信息素场，R^(ζ)为标签蚁群的影响区域；

细胞帧间匹配问题定义：

表示在第k帧估计得到的m个细胞位置，表示在第k+1帧得到的n个细胞候选集，于是构建一个二分图G＝(U,V,E)，其中， E表示从U到V上的边匹配假设的集合，帧间匹配的问题就转化成

且其中，D^(C)(·)表示对于匹配假设C前提下的成本函数，f(·)表示从U到V上的匹配假设映射，f^-1(·)表示f的射，|·|表示某集合的势；

成本函数D采用的是航迹最优次模式分配(OSPA-T)，即：给定两组带标签的集合X和Y，元素个数分别为m和n，X＝{(x₁,l₁),...,(x_m,l_m)}表示上一帧已估计的带标签l的细胞状态x的集合，Y＝{(y₁,s₁),...,(y_n,s_n)}表示在当前帧与之对应的匹配假设带标签s的候选量测y的子集，若m≤n,则OSPA-T定义为

其中，d^c((x,l),(y,s))＝min(c,d((x,l),(y,s)))表示截距距离，参数c和p分别表示对两集合元素个数差异和带标签状态变量差异的敏感度正系数(均大于零)，ψ_n表示在集合Y上的可行性扰动集合，ψ(i)表示ψ_n中的一个扰动，基距离d((x,l),(y,s))＝(d(x,y)^p′+d(l,s)^p′)^1/p′，p′表示基距离阶参数，d(x,y)表示两个状态x,y之间的p′范数距离，d(l,s)表示标签误差，若l＝s，则d(l,s)＝0，否则d(l,s)＝λ,λ∈[0,c]；若m＞n，有D_p,c(X,Y)＝D_p,c(Y,X)；

基于无标签侦察蚁群生成候选集

在第k+1帧，给定侦察蚁群初始随机分布在当前的图像区域，为对侦察蚂蚁行为进行建模，引入了一种新的蚂蚁动力学***，即从混沌行为逐步演化至确定性行为，从而发现食物源可能出现的区域；侦察蚂蚁的动力学行为定义如下

其中，表示侦察蚂蚁在第t次迭代时混沌移动变量，h_a表示蚂蚁自组织因子，取值范围在0和1之间，表示侦察蚂蚁在第t次迭代时的状态，参数ψ_d表示蚂蚁在搜索空间的搜索范围，参数V_a表示控制搜索空间的不同区域，a、b为常数，p_a(t-1)表示的是该侦察蚂蚁及其近邻蚂蚁在前t-1次迭代里得到的最优解；

表示侦察蚂蚁在第t次迭代时的近邻蚂蚁集合，假定侦察蚂蚁及其近邻蚂蚁在当前的状态、前t-1次迭代最优状态已知，则侦察蚂蚁在前t次迭代里得到的最优状态为

其中，符号∨表示取大运算，表示对似然函数g求逆运算，似然函数g的定义为

其中，I^l表示归一化的灰度图像在像素l的强度值，B_q表示以像素q为中心的方形区域所包括的像素集合，N_(q)表示像素q的近邻像素集合，符号|·|表示集合里元素的数目，符号∧表示求小运算；在第t次迭代上，若侦察蚂蚁在像素q上，则有

对于先后两次迭代中蚁群集合的蚂蚁状态视为两状态集合和则将OSPA-T简化为最优次模式分配(OSPA)距离于是，当m_a≤n_a时，有

若m_a＞n_a，定义

步骤1：初始化：给定若干组侦察蚁群，侦察蚂蚁的动力学行为中的混沌移动变量y_a(0)，自组织因子h_a，侦察蚂蚁初始状态搜索范围参数ψ_d，搜索空间的不同区域参数V_a，常数a、b等已知；

步骤2：对于每一组侦察蚁群，每个侦察蚂蚁按动力学行为模型移动，按上述似然函数定义来感知侦察蚂蚁所在相应像素的似然函数值，并与其近邻侦察蚂蚁进行信息交换得到截止目前迭代次得到的最优状态；

步骤3：若侦察蚂蚁所在的像素位置是已探知可能细胞所在位置区域，则将其对应的似然函数设为一很小值；

步骤4：计算本次迭代侦察蚁群状态集合与上一次迭代侦察蚁群状态集合的最优次模式分配距离的差值，若差值小于某一个阈值，则终止该侦察蚁群的侦察工作，进而通过聚类重构，生成若干个侦察子蚁群；

步骤5：重复步骤2-4，直到没有新的侦察子蚁群生成为止；将已生成的每个侦察子蚁群的中心作为这样，共生成n个候选量测

2.根据权利要求1所述的基于标签蚁群的联合多细胞跟踪方法，其特征在于还包括：基于标签蚁群的匹配构建与参数估计联合跟踪

将标签多贝努利随机有限集参数集{(r^(ζ),π^(ζ)):ζ∈Ψ}与标签多贝努利蚁群参数集存在一一对应关系，且标签贝努利蚁群与细胞有一一对应关系，即

假定标签为ζ蚁群以一个生存概率会从上一帧位置移到下一帧所在位置，除了其标签的唯一性不变以外，生存概率定义为：

其中，H，W表示当前帧图像的高与宽，表示在第k帧标签蚁群中心离图像边界最短的高度、宽度方向距离，表示该标签蚁群在第k帧的存在概率；

同样，标签蚁群在下一帧k+1发生***的概率定义为

其中，exp(·)表示指数函数，分别表示细胞在第k帧和初始帧估计所得到的轮廓宽高比，分别表示细胞在第k帧和初始帧的平均强度，υ为0和1之间的调节系数，表示对应细胞的局部窗口面积大小，表示在第k帧所估计出的细胞面积大小；

基于上述定义，基于标签蚁群的帧间匹配与参数估计一体化方法描述如下：

步骤1：初始化：假定在上一帧k共有m个细胞估计，则有 为细胞i的估计位置，其所对应的标签多贝努利蚁群参数集可表示为对于每个贝努利蚁群，蚁群规模均为由无标签侦察蚁群生成n个候选量测集 表示候选量测j位置；

步骤2：对于任一标签贝努利蚁群任选其中一标签蚂蚁(与所在的蚁群标签一致)，随机生成的两随机数κ_S∈(0,1)和κ_M∈(0,1)，若满足：

(1)：则该标签蚂蚁能从上一帧生存下来，并从其所对应的细胞估计位置出发，选择候选量测的概率为

其中，1_Y(X)表示包含函数，若则等于1，否则为0；表示第个细胞估计位置所对应的候选量测子集，表示在轨迹上的轨迹信息素强度，θ(·)表示轨迹信息素的身份特征函数，它与细胞相关；表示启发式函数，并有d(·,·)表示距离函数，α₁,β₁表示相对重要性系数；

(2)：说明该蚂蚁不仅能从上一帧生存下来，而且还要应对细胞发生***可能，于是，该标签蚂蚁复制成两只一样的标签蚂蚁，先后从其所对应的细胞估计位置出发，选择不同候选量测的概率为

(3)：或该蚂蚁不动，说明在下一帧它所对应的细胞可能消失；

当标签蚂蚁按照(1)或(2)决策后，则在经过的轨迹上释放一定量的信息素：

其中，ρ表示信息素残留系数，表示在当前帧k+1对应标签贝努利蚁群在上一次迭代更新后所得到的存在概率，表示细胞的影响区域；

步骤3：当m个标签贝努利蚁群按照步骤2完成相应的决策后，就构成了一个从U到V之间的一个匹配假设C，并计算所对应的OSPA-T值；

步骤4：按步骤2-3，进行次循环，在次循环中，得到本次迭代最佳匹配假设C，其对应的候选集表示为对任按步骤2所得到的标签贝努利子蚁群开始细胞状态估计，即觅食；若觅食蚂蚁在像素p上，它移动到像素q的概率定义为

其中，α₂,β₂表示相对重要性系数，N_(p)表示像素p近邻像素集合，表示在当前帧k+1像素q上的食物信息素强度；伴随着上述的决策，标签蚂蚁同时在像素q上释放一定量的信息素

其中，ε₁为阈值系数；

当标签贝努利子蚁群中所有蚂蚁完成上述觅食决策后，对应的食物信息素场进行更新

其中，表示在像素q上总的信息素输入，表示在像素q上总的信息素扩散输入，定义为

其中，D表示扩散常系数，取值范围在(0,1)之间，κ为控制常系数，表示像素q和之间的距离，R_q表示像素q的影响区域；

利用多次迭代后，由所得到的食物信息素场和启发式函数，得到标签贝努利子蚁群的贝努利参数近似

其中，表示当前帧k+1觅食标签蚁群的概率密度离散分布，表示蚂蚁的状态所对应的权重，表示标签贝努利子蚁群中的蚂蚁数；对于食物信息素场，进行二值化、形态学和填充三步简单处理，即可得到候选集对应细胞的轮廓及其位置

步骤5：对于截止目前迭代所得到的最佳匹配假设C，对其对应的候选集位置，轨迹信息素场，食物信息素场进行全局更新

其中，Q₁，Q₂为信息素加强常数，表示截止目前迭代所得到的最佳匹配假设C对应的位置更新后的候选集集合；

步骤6：重复步骤2-5，直至收敛条件满足；同时，对所得到的最佳匹配迭代所得到的食物信息素场，进行二值化、形态学和填充三步简单处理，若得到的标签子蚁群存在概率可得到对应细胞的轮廓及其位置估计；

对于所得到的帧间最佳匹配，将其轨迹信息素场提取，并保留轨迹信息素的身份特征，通过帧帧顺序匹配，形成细胞谱系树；

对于未匹配的候选量测，视为当前帧可能出现的新生项，对应生成相互独立的标签贝努利子蚁群，按步骤4对其对应的食物信息素场构建，若新生的标签子蚁群存在概率并不低于0.5，则提取出相应的细胞的轮廓及其位置估计；

对于多个标签贝努利子蚁群所得到的位置估计近邻且相互距离小于某阈值，则将存在概率大的标签贝努利子蚁群保留，其他满足条件的标签贝努利子蚁群直接去除并不再考虑；同时，对于标签更新，将最先出现的标签赋予存在概率大的标签贝努利子蚁群。

3.根据权利要求1所述的基于标签蚁群的联合多细胞跟踪方法，其特征在于还包括基于四阶段的航迹恢复：

基于四阶段的航迹恢复的相似度标准如下：

相似度标准1： 分别表示航迹在诞生和结束时刻所对应的细胞状态，对于两航迹和若结束时刻先于则相似度标准1定义为

其中，R^-表示满足帧间最大移动距离约束条件的航迹集合，且其任一航迹结束时刻先于 表示对应的航迹在诞生时刻细胞的面积估计，表示对应的航迹在结束时刻细胞的面积估计，∩表示两细胞面积交集，∪表示两细胞面积求和；

相似度标准2：对于航迹其在诞生时刻k所对应的细胞状态为对于航迹在时刻k和k-1所对应的细胞状态分别为和则相似度标准2定义为

其中，C_h(·)和C_w(·)分别表示细胞轮廓的主轴和短轴的长度，ΔC_h(·,·)和ΔC_w(·,·)分别表示细胞轮廓的主轴和短轴的长度差绝对值，ΔA(·,·)表示细胞面积差的绝对值；

具体步骤如下：

阶段1：给定一长度为3个采样间隔的时间窗口，从最后一帧开始反向滑动，当其右侧边与航迹的起始位置相交时，对应的帧若为第k帧，此相交的航迹将与其它在这时间窗口内(从第k-3帧至第k-1帧)结束的航迹进行关联；若第k-1帧内有相应的航迹候选(即该航迹在第k-1帧结束)，且满足帧间最大移动距离约束条件，则按相似度标准1选择最大一个与之关联；否则，进入阶段2；

阶段2：排除阶段1已关联的若非空，且第k-2，k-3帧内有相应的航迹候选(即该航迹在第k-2，k-3帧结束)，并满足帧间最大移动距离约束条件，则按相似度标准1选择最大一个与之关联；

阶段3：排除阶段1，2已关联的若仍非空，与近邻(满足帧间最大移动距离约束条件)的完整航迹在第k，k-1帧进行相似度标准2计算，若相似度值大于某给定阈值，则航迹并入航迹对应细胞***；

阶段4：排除阶段1，2，3已关联的若仍非空，则认为该起始的航迹(且航迹寿命至少为3个采样间隔)为正式航迹，否则，去除该航迹。