Distribution of Zero-effort Miss Distance Estimation Errors in Continuous-time Controlled System With Mode Mismatch
-
摘要: 针对大机动目标拦截中零控脱靶量（Zero-effort miss distance，ZEM）估计误差分布的求解问题，在线性估计器与独立模式辨识器的配置下，提出了一种解析求解方法.在存在固定的模式辨识延迟下，ZEM的估计误差服从有偏的高斯分布.最后，通过与蒙特卡洛仿真的比较，验证了本文方法的正确性.Abstract: For the problem of solving the error distribution of zero-effort miss distance (ZEM) in highly maneuvering target interception, an analytical approach is proposed for the configuration of a linear estimator and a separate mode decision-maker. In the case of a fixed mode decision delay, the error of ZEM follows a biased Gaussian distribution. Finally, the correctness of the proposed method is verified by comparison with Monte Carlo simulation.
-
乳腺癌是导致女性死亡的主要癌症之一 [1], 它严重影响着女性的身心健康. 早期诊断和早期综合治疗是防止乳腺癌的最有效手段 [2].由于成本低廉, 性价比高、无辐射以及无创伤等原因, 超声成像技术已经成为检测乳腺癌的重要手段.然而, 乳腺超声(Breast ultrasound, BUS)图像往往需要有临床经验的医生来进行判读, 医生判读的准确性又会受到各种内部或外部因素的影响.为了克服这种缺陷, 帮助医生提高诊断的准确率和客观性, 计算机辅助诊断(Computer aided diagnosis, CAD)系统越来越多的被应用于临床实践中 [3-5].
分割是乳腺超声图像CAD系统中重要的一个环节, 通过分割可以定位肿瘤区域的位置, 为后续肿瘤特征提取以及良恶性分类提供必要的信息.由于超声图像具有高噪声、低对比度、边缘模糊不清等特点, 超声图像的分割成为图像处理领域中一个难度较高、亟待解决的问题.近年来, 许多基于不同模型的乳腺肿瘤分割方法已经被提了出来, 例如:区域增长方法(Region growing)[6]、 马尔科夫随机场方法(Markov random field, MRF)[7]、 主动轮廓方法(Active contour)[8]、 神经网络方法(Neural network)[9]、 水平集方法(Level set)[10]等. 国内的学者也提出了许多不同的方法 [11-13].虽然现有的一些分割方法通过结合多种先验知识, 已经取得了较好的效果, 但是这些方法仍然存在一些问题:
1) 乳腺超声图像对比度较低, 而且在图像的脂肪层以及肌肉层存在较多的与肿瘤强度接近的低回声区域, 现有的依赖空间域特征的方法不能有效地将肿瘤与低回声区域区分开来, 特别是当低回声区域与肿瘤区域的边界较为接近时, 现有的方法往往会出现过分割.
2) 现有的乳腺超声图像分割方法都将单帧静态超声图像作为研究对象, 然而单帧静态图像只能反映肿瘤某一侧面的信息.实际临床中, 医生在扫查和诊断中用的是超声视频序列, 这个视频序列与单帧静态图像相比, 可以提供更完整、更全面的肿瘤信息.
为了解决上述问题, 本文提出了一种基于多域先验的协同分割模型来对乳腺超声序列进行分割. 协同分割(Co-segmentation)方法最早由Rother等 [14]提出, 是一种无监督学习的分割方法, 目的是将多幅具有相同或相似目标的图像分割为前景和背景.目前越来越多的人开始对协同分割方法进行研究. 现有的大多数协同分割方法是基于马尔可夫随机场(MRF)的优化问题 [15], 其基本思想是建立起包括图像内部能量项和图像间全局能量项的能量方程, 然后将能量方程转化到MRF问题上进行优化. 这类方法通常是在全局能量项上进行改进, 使全局能量项更合理或者使整个模型的优化更便利. 除了基于MRF的方法外, 还有一些其他的协同分割方法被提出. Joulin等 [16]从聚类的角度出发, 通过将图像内部聚类与全局图像聚类结合起来, 提出了基于聚类的协同分割方法. Kim等 [17]从各向异性热扩散的原理出发, 将协同分割问题视作温度最大化问题来解决. Meng等 [18]通过将多幅图像构建成图, 然后使用最短路径的方法来进行协同分割. 现有的协同分割方法大都是在自然图像上进行分割的, 方法使用的特征也是从自然图像的角度考虑的. 与自然图像相比, 乳腺超声图像有较低的分辨率和对比度, 图像中各个组织之间的边界比较模糊, 图像内部有较多的噪声, 所以已有的协同分割方法在乳腺超声图像中并不适用.
本文提出的多域协同分割模型充分结合了乳腺超声图像的特点.一方面, 通过考虑乳腺超声图像在空间域的强度分布、 位置和姿态信息, 与频率域的边缘信息相结合, 建立起超声图像内部的能量关系; 另一方面, 本文结合了乳腺肿瘤图像的特征来对协同分割框架中的能量项进行了设计, 建立起超声序列之间的全局能量模型.该模型可以有效地利用序列图像特点进行乳腺肿瘤分割.
1. 基于多域先验的协同分割模型
本文提出的多域协同分割模型分割步骤如图 1所示.模型处理的是乳腺超声序列.在使用多域协同分割模型分割之前, 需要对输入序列进行预分割.肿瘤在乳腺超声图像中表现为"颜色相对较暗、团块状"区域.然而, 图像中的肌肉和脂肪组织也有类似的外观.预分割可以给出肿瘤所在区域的位置, 使肿瘤区域成为前景区域, 肌肉、脂肪和其他区域成为背景区域.这是个粗定位的过程, 无法获得肿瘤准确的边界.文献[19]提出了一种基于单帧显著性的乳腺肿瘤检测方法, 文章首先根据医学先验对乳腺肿瘤进行定位, 然后根据乳腺的背景特征以及强度特征构建起解剖学线索和对比度线索, 最终得到乳腺肿瘤的显著性映射结果. 因为文献[20]中的方法具有全自动、 定位准确等特点, 本文用其进行预分割.预分割算法的复杂度为 O$(K^3)$, 其中, K为图像中超像素块的数量.
1.1 分割模型概述
本文提出的多域协同分割模型包括两个部分:内部能量项和全局能量项.为了便于后续表示, 使用$I = \{ {I_k}\} _{k = 1}^M$ 表示含有M帧图像的乳腺超声序列, 使用 ${{x}_{k}} = \{ x_i^k\} _{i = 1}^N$ 表示第${I_k}$帧图像对应的二值标签集合, 其中图像${I_k}$含有的像素个数为N, $x_i^k$表示图像${I_k}$的第i个像素对应的二值标签, $x_i^k$取值为1时其为前景, 取值为0时其为背景.本文的目标是通过多域协同分割模型来实现更优分割, 也就是对标签集合 ${{X}}(\{ {{{x}}_{{k}}}\} _{k = 1}^M)$实现最优分配, 这个问题可以通过最小化式(1)中的能量方程来实现:
$E({{X}}) = {E_{\rm In}} + {E_{\rm Global}} $
(1) 式中, ${E_{\rm In}}$表示单帧图像内部的能量项, 用来约束单帧图像内前景像素与背景像素的关系, ${E_{\rm Global}}$ 表示序列之间的全局能量项, 用来构建序列之间的能量关系.其中${E_{\rm In}}$又分为空域能量项和频域能量项两部分, 所以式(1)可以表示为
$E({{X}}) = \lambda {E_{SC}} + {E_{FC}} + {E_{\rm Global}} $
(2) 式中, ${E_{SC}}$表示空域能量项, 用来对肿瘤的强度、 位置和姿态进行描述. ${E_{FC}}$表示频域能量项, 用来对肿瘤的边缘信息进行描述, ${E_{SC}}$和${E_{FC}}$共同建立起超声图像内部的能量关系. $\lambda$ 是空域能量项的权值, 用于调节空域能量项与频域能量项在模型中所占的比重.以下分别对模型中的各项进行详细说明.
1.2 空域能量项
由于肿瘤区域与低回声区域在强度上具有较大的相似性, 仅使用强度信息很难将之区分, 并容易产生过分割.文献[4]从单帧肿瘤图像的特点出发, 提取了肿瘤姿态、 位置和强度分布等空间域信息来区分肿瘤与正常组织区域, 本文利用文献[4]中的思想, 并将其扩展到序列图像中, 得到了如下的空域能量项表示:
${{E}_{SC}}(X)=\sum\limits_{k=1}^{M}{\left[ \sum\limits_{i}{-\ln (P_{i}^{k}H_{i}^{k})+\sum\limits_{(i, j)\in {{N}_{i}}}{H_{(i, j)}^{k}}} \right]}$
(3) 式中, $P_i^k$是对图像${I_k}$在像素i处的肿瘤姿态与位置描述. $H_i^k$是对图像${I_k}$在像素i处的强度描述, $H_{(i, j)}^k$中像素i处的邻域像素不连续性程度的描述, ${N_i}$表示图像${I_k}$中像素i的邻域集合.
1.2.1 姿态与位置描述
肿瘤的姿态和位置是区分良恶性肿瘤的重要特征, 能不随图像亮度和对比度的变化而变化 [4]. 由于本文使用的所有乳腺超声图像中均只含有一个肿瘤, 所以这里使用2维椭圆高斯方程来表示肿瘤的姿态与位置:
$P_{i}^{k}(x_{i}^{k})=G_{i}^{k}x_{i}^{k}+(1-G_{i}^{k})\bar{x}_{i}^{k}$
(4) $\begin{align} & G_{i}^{k}(i)=\exp \{-[a{{({{i}_{x}}-{{i}_{{{x}_{0}}}})}^{2}}+2b({{i}_{x}}-{{i}_{{{x}_{0}}}})({{i}_{y}}-{{i}_{{{y}_{0}}}})+c{{({{i}_{y}}-{{i}_{{{y}_{0}}}})}^{2}}] \\ & +c{{({{i}_{y}}-{{i}_{{{y}_{0}}}})}^{2}}]\} \\ \end{align}$
(5) 式中, 使用的椭圆由预分割结果中的前景拟合得到, 拟合时, 将前景区域的质心做为椭圆的中心, 将与前景区域具有相同标准二阶中心矩的椭圆的长轴作为椭圆长轴, 将与前景区域具有相同标准二阶中心矩的椭圆的短轴作为椭圆短轴, 将与前景区域具有相同标准二阶中心距的椭圆的长轴与水平轴的夹角作为椭圆的倾角. $({i_{{x_0}}}, {i_{{y_0}}})$和$({i_x}, {i_y})$分别表示图像${I_k}$中椭圆的中心坐标和像素i 处的坐标. $a, b, c$用来控制椭圆的姿态, 参数$a, b, c$由式(6)给出:
$a=\frac{{{\cos }^{2}}\theta }{2\sigma _{x}^{2}}+\frac{{{\sin }^{2}}\theta }{2\sigma _{y}^{2}}, b=-\frac{\sin 2\theta }{4\sigma _{x}^{2}}+\frac{\sin 2\theta }{4\sigma _{y}^{2}}, c=-\frac{{{\sin }^{2}}\theta }{2\sigma _{x}^{2}}+\frac{{{\cos }^{2}}\theta }{2\sigma _{y}^{2}}$
(6) 式中, $\theta$表示椭圆的倾角, ${\sigma _x}$和${\sigma _y}$分别表示1/2长度的椭圆长轴与短轴.椭圆高斯方程能够将分割结果有效地约束在自身的范围内, 排除远离肿瘤区域中与肿瘤区域相似的低回声区域.
1.2.2 强度约束
一些方法在使用强度分布对乳腺超声图像进行分割时, 会先假设图像的强度满足一定的先验分布, 例如: 瑞利(Rayleigh)分布、 伽玛(Gamma)分布或者指数(Exponential)分布等, 由于乳腺超声图像在采集时可能使用不同的设备, 或者使用不同的参数来进行采集, 这些预先定义的强度分布与图像的实际情况会有一定的偏差.因此, 本节在对强度进行定义的时候, 并没有假设图像的强度分布, 而是使用图像的前景与背景直方图来进行表示, 这种表示方法比预先假设强度分布具有更好的鲁棒性. 强度表示的一阶定义如式(7):
$H_i^k(x_i^k) = \Pr (y_i^k|obj)x_i^k + \Pr (y_i^k|bkg)\overline {x}_i^k $
(7) 式中, $y_i^k$表示图像$I_i^k$在像素i处的强度值, $\Pr (y_i^k|obj)$和$\Pr (y_i^k|bkg)$分别表示图像${I_k}$中的像素i属于前景和背景的类条件概率密度, 条件概率密度的值由前景直方图和背景直方图来确定. 图像强度的二阶定义由相邻像素不连续性的程度来定义, 其表示形式如式(8):
$H_{(i, j)}^k = [x_i^k \ne x_j^k]\exp \left(\frac{ -{\left\| {y_i^k -y_j^k} \right\|^2}}{2{\sigma ^2}}\right) $
(8) 式中, $[\phi]$为指示函数, 当 $\phi$为真时, 取值为1; 当$\phi$为假时, 取值为0, $y_i^k$和$y_j^k$分别表示图像${I_k}$ 在像素i和j处的强度值. $H_{(i, j)}^k$表示图像${I_k}$在像素$i $处的不连续性, 当i处比较平滑时, $H_{(i, j)}^k$值相对较小, 当i处不平滑时, $H_{(i, j)}^k$ 值相对越大. $\sigma$ 用来控制不连续性的大小, 实验表明, 当 $\sigma$ 取15时能达到最好效果.
1.3 频域能量项
因为乳腺超声图像在采集时使用的参数会有变动, 其亮度与对比度也会随之变化. 传统的基于空间域的分割方法在检测边缘时, 往往会受到这些变化的影响, 进而不能得到满意的结果. 相位一致性 [20]与零交叉 [21]是两种重要的频率域边缘检测方法, 能够不随图像亮度与对比度的变化而变化, 目前已经被广泛地应用于图像处理领域 [22].本部分将相位一致性与零交叉方法相结合, 构建起了模型的频域能量项:
$\begin{align} &{{E}_{FC}}\left( X \right)=\sum\limits_{k=1}^{M}{\sum\limits_{\left( i, j \right)\in {{N}_{i}}}{ED_{\left( i, j \right)}^{k}=}} \\ &\sum\limits_{k=1}^{M}{\sum\limits_{\left( i, j \right)\in {{N}_{i}}}{\left[ x_{i}^{k}\ne x_{j}^{k} \right]\times ZC_{\left( i, j \right)}^{k}\left( 1-\max \left( PCM_{i}^{k}, PCM_{j}^{k} \right) \right)}} \\ \end{align}$
(9) 式中, $ZC$和$PCM$分别表示零交叉检测与最大能量响应对应的相位一致性值. 本节中相位一致性的计算采用文献[4]中的方法.该方法将图像与一系列不同尺度不同方向的 Log-Gabor滤波器做卷积, 然后进行相位一致性的求解. 滤波器最大响应处对应的尺度和方向能够反映图像局部区域的平滑程度以及方向. 如果一个点在较大的尺度上有较大的响应, 则此点附近的区域较为平滑; 相反, 如果一个点在较小的尺度上有较大的响应, 则此点附近的区域不平滑 [4]. 在方向方面, 如果一个点在方向角$\theta$ 上有较大的响应, 则此点处对应的方向角接近于$\theta$. 因此, 肿瘤边缘的尺度和方向与滤波器的尺度和方向相同时, 滤波器的局部响应最大. 本文中, 选择一个点处Log-Gabor滤波器最大实部${e_{ns, no}}$对应的尺度$(ns)$ 和方向$(no)$作为该点的边缘尺度和方向.
$(s_i^k, n_i^k) = \arg \mathop {\max }\limits_{ns, no} \{ {e_{ns, no}}(i)\}, i \in {I_k} $
(10) 零交叉检测ZC的形式如式(11):
$ZC_{(i, j)}^{k}=\text{ }\frac{1}{2}[1-\text{sgn}({{e}_{(s_{i}^{k}, n_{i}^{k})}}(i))\cdot \text{sgn}({{e}_{(s_{j}^{k}, n_{j}^{k})}}(j))], j\in {{N}_{i}}, i\in {{I}_{k}}$
(11) $ZC_{(i, j)}^k$能够决定点i和j之间是否存在一条边, 如果符号函数${e_{s_i^k, n_i^k}}$与${e_{s_j^k, n_j^k}}$ 相同, 则$ZC_{(i, j)}^k$为0, 说明两点之间不存在边, 如果不相同, 则$ZC_{(i, j)}^k$ 为1, 说明两点之间存在边.
最大能量响应对应的相位一致性值PCM的形式如式(12):
$PCM_{i}^{k}=\cos [{{\varphi }_{(s_{i}^{k}, n_{i}^{k})}}(i)]\text{ }PC_{i}^{k}, i\in {{I}_{k}}$
(12) 式中, $PCM_i^k$能够检测到相位0或$\pi $处的边缘, 由于其构成具有软间隔特性, 所以可以检测到模糊边缘(相位接近0或$\pi $). ${\varphi _{(s_i^k, n_i^k)}}$表示图像${I_k}$在像素i处的相位, $PC_i^k$表示图像${I_k}$在像素i处的相位一致性值. $PC_i^k$的值可以通过文献[23]中的方法得到.
1.4 全局能量项
乳腺超声序列中, 相邻的乳腺超声图像在前景区域具有相似性. 协同分割方法正是从序列图像前景区域的相似性出发构建起全局能量项的, 这种方法通过对多帧图像的前景和背景建模, 得到更准确的前景与背景分布, 进而实现比单帧分割更准确的序列分割.
现在的一些协同分割方法 [24]在构建全局能量项时, 通常是对两幅图像使用强度直方图来进行约束. 由于乳腺超声序列中图像数量多于两幅, 这种只能处理成对图像的方法不能够充分地利用序列信息.因此, 本部分使用了文献[15]中的基于高斯混合模型(Gaussian mixture model, GMM)的方法来定义全局能量项. GMM有以下优点: 1) GMM是一种概率分布模型, 具有较高的准确率; 2) GMM能够详细地描述一个目标不同区域的分布情况.本文中全局能量项的定义如式(13):
${{E}_{\text{Global}}}(X)=\sum\limits_{k=1}^{M}{\sum\limits_{i}{-\ln [\Pr (z_{i}^{k}|{{\theta }_{coobj}})}}x_{i}^{k}+\Pr (z_{i}^{k}|\theta _{bkg}^{k})\overline{x_{i}^{k}}]$
(13) 式中, $\Pr ( \cdot )$为GMM的概率分布, $z_i^k$表示图像${I_k}$在像素i处的特征(强度特征与PCM特征), ${\theta _{coobj}}$ 为图像序列中所有图像前景的GMM参数, $\theta _{bkg}^k$为图像${I_k}$的背景像素的GMM参数, GMM参数的计算过程如下: 首先, 使用K-means算法分别将序列图像的所有前景区域分为$Fc$个类别, 将序列中每帧图像的背景区域分为$Bc$ 个类别, 类别个数即为GMM中高斯模型的个数.然后, 将每个类别中所有像素特征的均值、 协方差矩阵以及本类别中像素占总像素的比例分别作为本类别对应高斯的均值、 协方差矩阵和高斯权重.这样就可以得到前景和背景的GMM参数.实验结果表明, 当$Fc = 3$, $Bc = 5$时有较好的效果.考虑到乳腺超声序列中前景相似的特性, 对${E_{\rm Global}}$ 的最小化能够保证图像序列前景区域的一致性与连续性.
1.5 模型优化
本文的目标是从乳腺超声序列中分割出肿瘤, 这是一个二分割问题. Boykov等提出的交互式图割 [25]为本文模型的优化提供了一种有效的手段. 在图割方法中, 标准的能量方程形式如式(14)所示 [26]:
$E\left( L \right) = \sum\limits_{p \in P} {{D_p}\left( {{L_p}} \right) + \sum\limits_{(p, q) \in N} {{V_{p, q}}\left( {{L_p}, {L_q}} \right)} } $
(14) 式中, $L = \left\{ {{L_p}|p \in P} \right\}$为图像P的标签集合, ${D_p}\left( \cdot \right)$为数据惩罚项, ${V_{p, q}}$为平滑项, N为所有相邻像素对的集合.本文提出的多域协同分割模型中, 式(3)空域能量项中的$P_i^kH_i^k$对应数据惩罚项, $H_{\left( {i, j} \right)}^k$对应平滑项, 式(9)频域能量项对应平滑项, 式(13)全局能量项对应数据惩罚项.这样, 本文提出的模型可以映射到标准能量方程的形式上去. 又因为本文模型满足文献[27]中的图表示形式, 即:
$\begin{align} &H_{\left( i, j \right)}^{k}\left( 0, 0 \right)+H_{\left( i, j \right)}^{k}\left( 1, 1 \right)+ED_{\left( i, j \right)}^{k}\left( 0, 0 \right)+ \\ &ED_{\left( i, j \right)}^{k}\left( 1, 1 \right)\le H_{\left( i, j \right)}^{k}\left( 0, 1 \right)+H_{\left( i, j \right)}^{k}\left( 1, 0 \right)+ED_{\left( i, j \right)}^{k}\left( 0, 1 \right)+ED_{\left( i, j \right)}^{k}\left( 1, 0 \right) \\ \end{align}$
所以本文提出的模型能够使用最大流最小割算法1来进行全局最优化, 模型复杂度为O$\left( {{n^2}m\left| C \right|} \right)$, 其中n为顶点数, m为边数, $\left| C \right|$为最小割的数量.
1.6 算法实现
因为在开始求解模型时, 乳腺超声序列的前景和背景是通过预分割得到的, 预分割结果比较粗糙, 可以通过迭代的方式来对乳腺超声序列进行分割.迭代终止的条件可以为: GMM的参数满足收敛条件, 或者达到最大迭代次数(通常使用5 $\sim$ 10次迭代).完整的分割算法流程如下:
步骤 1. 使用文献[19]中的方法对乳腺图像序列进行预分割;
步骤 2. 使用图割算法对式(1)中的${E_{\rm In}}$进行最小化, 产生初步分割结果;
步骤 3. 初始化GMM的参数${\theta _{coobj}}$, $\theta _{bkg}^k$;
步骤 4. 使用图割算法对式(2)进行最小化;
步骤 5. 根据式(2)的优化结果更新参数$\theta _{coobj}$, $\theta _{bkg}^k$;
步骤 6. 循环执行步骤3和步骤5直到满足收敛条件或者达到最大迭代次数为止.
2. 实验结果
2.1 数据集
为了验证本文提出方法在乳腺超声序列图像上的分割效果, 本文对50组乳腺超声序列进行了分割. 50组序列中包括25组恶性序列和25组良性序列, 每组序列包括5 $\sim$ 12帧图像, 总计400帧图像, 图像的平均大小为600像素 × 480像素. 本文采用的所有乳腺超声图像序列由哈尔滨医科大学附属第二医院超声科提供, 并由超声专家对肿瘤区域进行手工标注, 作为金标准验证分割的效果.本文实验在一台4核2.13 GHz CPU、12 GB内存的计算机上运行, 使用Matlab来进行算法实现, 平均处理一帧图像的时间为10.02 s.
2.2 评价标准
为了对分割效果进行评价, 本文使用了面积和边界两种误差评价标准.其中面积评价标准包括: 真阳性比(True positive ratio, TPR)、假阳性比(False positive ratio, FPR)和总体相似度(Similarity ratio, SIR):
$\text{TPR}=\frac{\left| {{A}_{m}}\cap {{A}_{r}} \right|}{\left| {{A}_{m}} \right|}, \text{FPR}=\frac{\left| {{A}_{m}}\cap {{A}_{r}}-{{A}_{m}} \right|}{\left| {{A}_{m}} \right|}\text{, SIR}=\frac{\left| {{A}_{m}}\cap {{A}_{r}} \right|}{\left| {{A}_{m}}\cup {{A}_{r}} \right|}$
(15) 其中, ${A_m}$表示医生手工标定的肿瘤区域, ${A_r}$为算法分割出的肿瘤区域. TPR指标越高, 则分割结果覆盖手工标定区域的程度越高; FPR指标越低, 则覆盖的错误区域越少; SIR指标越高, 则分割结果与手工标定的区域越接近 [28-29].
豪斯多夫距离(Hausdorff distance, HD)和平均绝对距离(Mean absolute distance, MD)用来评价两个边界的差异程度 [28]. 这里将超声专家标注的边界记为$Q = \left\{ {{q_1}, {q_2}, \cdots , {q_\eta }} \right\}$将分割得到乳腺肿瘤边界记为$P = \{ {p_1}, {p_2}, \cdots, p_j, {p_\alpha }\} $, Q和P对应于同一个肿瘤边界. 对于P上的任意一点${p_j}$, 定义:
$d({p_j}, Q) = \mathop {\min }\limits_w \left\| {{p_j} -{q_w}} \right\|, w = 1, \cdots, \eta $
(16) 其中, $\left\| \cdot \right\|$表示2维欧氏距离, $\alpha $和$\eta $ 分别表示边界P和Q上点的个数. HD和MD的定义如下:
$\text{HD}=\underset{j}{\mathop{\max }}\, d({{p}_{j}}, Q), j=1, \cdots , \alpha $
(17) $\text{MD}=\frac{\sum\limits_{j=1}^{\alpha }{d({{p}_{j}}, Q)}}{\alpha }$
(18) HD表示两个边界的最大不一致程度, MD表示两个轮廓的平均不一致程度. 本文中使用所有肿瘤图像的平均豪斯多夫距离(Average Hausdorff distance, AHD)以及平均绝对误差距离(Average mean absolute distance, AMD)来对分割结果的边界进行评估, AHD与AMD值越大, 分割得到的边界与手工标定的边界差别越大.
2.3 结果
为了全面展示本文提出的模型, 本节首先展示模型的部分中间结果以及模型中参数的变化对最终分割结果的影响. 然后, 将本文提出的模型与三种最新的全自动分割方法 [4, 30-31]进行对比分析, 从模型对与肿瘤相似区域的分割能力, 模型处理低对比度、边界模糊的超声图像的能力以及模型处理单帧分割不容易区分的能力等方面对模型进行评价(由于篇幅有限, 每个序列取3帧图像表示). 实验中, 所有的实验数据均是在相应方法各自的最优参数条件下获得.
图 2显示了本文提出算法的部分中间结果, 图 2(a) 是原始图像序列, 图 2(b) 是原始图像序列的预分割结果, 预分割能够定位肿瘤的位置, 并进行简单分割. 图 2(c) 是乳腺肿瘤序列在空间域中得到的2D高斯椭圆, 用来对肿瘤的姿态与位置进行约束, 从图 2中可以看出, 高斯椭圆能够准确地限制肿瘤区域的位置, 排除背景噪声的干扰.图 2(d) 是乳腺肿瘤序列对应的强度分布序列, 强度图能够有效地利用空间域中的强度信息来表示肿瘤, 图 2(e) 是乳腺肿瘤序列在频率域对应的边缘检测结果, 从图像中可以看出, 最大能量响应处对应的相位一致性值能够明显地表示出肿瘤图像中的边缘信息.图 2(c) 和(d)为空间域中数据项对应的中间结果, 图 2(e) 为频率域的部分中间结果.图 2(f) 是本文提出的多域协同分割模型最后的分割结果.
参数$\lambda $用来控制空间域信息和频率域信息对分割结果的影响.当$\lambda $ 较小时, 空域能量项对分割结果的影响较小, 频域能量项对结果的影响较大, 分割结果对边缘信息更为敏感, 容易将疑似边缘区域作为肿瘤边缘. 随着$\lambda $增大, 空域能量项对分割结果的影响逐渐增大.当$\lambda $ 取较大值时, 空域能量项会对分割结果产生较大影响, 分割结果对强度信息更为敏感, 容易将与肿瘤区域近似的组织区域误判断为肿瘤, 产生过分割. 图 3(b) $\sim$ 图 3(d) 分别是采用$\lambda= 0.05$, $\lambda = 10$和$\lambda = 0.3$ 对图 3(a)的序列分割的结果. 采用过小的$\lambda $值($\lambda = 0.05$)时, 模型过多考虑了频率域信息, 而采用过大的$\lambda $值($\lambda = 10$)时, 分割产生了过多细小的区域和曲折的边界曲线, 即过分割. 实验表明, 当$\lambda = 0.3$时, 模型分割出的区域与实际肿瘤区域最相近.
图 4所示的序列中, 肿瘤图像内部含有较多与肿瘤区域相似的正常组织. 图 4(b) $\sim$ (e) 分别是文献[28]提出的基于人工神经网络的模型, 文献[29]提出的基于目标识别的分割模型, 文献[4]提出的基于图分割的模型与本文提出的多域协同分割模型对图 4 (a)所示序列的分割结果. 通过对比图 4 (f)医生标定的肿瘤区域可以看出, 文献[28]中的分割结果中包含较多的正常组织, 文献[4, 29]则正确的区分了肿瘤区域与正常组织. 本文提出的模型最接近医生手工标定的肿瘤区域.
图 4 4 含较多与肿瘤区域相似的正常组织的超声序列分割结果Fig. 4 Segmentation results of normal breast ultrasound sequence similar to tumor regions图 5展示了4种分割方法在低对比度、边界模糊的乳腺超声序列上的分割结果. 通过对比医生标定的肿瘤区域可以看出, 本文提出的协同分割模型对于低对比度、 边界模糊的超声图像仍可以准确地分割(图 5 (e)), 文献[4]中的方法也得到了相似的结果, 但是其上边缘部分没有完整地收敛到肿瘤的边界. 文献[29]中的方法只能识别部分肿瘤区域, 文献[28]中的方法虽然将肿瘤区域正确识别, 但其结果中包含了较多的正常组织.
图 5 低对比度、边界模糊的超声序列的分割结果Fig. 5 Segmentation results of breast ultrasound sequence with low contrast and blurry boundaries图 6和图 7分别是4种分割方法在单帧不容易区分的乳腺超声序列上的分割结果. 观察这些序列中的图像可以看出, 这些序列中部分肿瘤图像对比度较低, 肿瘤内部不均匀, 有部分低回声区域同时也有部分较亮的区域(图 6的2帧图像, 图 7的第2帧图像).在使用单帧分割时, 往往会将肿瘤内部较亮与较暗区域的分界处作为边界, 造成误分割(图 6 (b) $\sim$ (d), 7 (b) $\sim$ (d)). 本文提出的多域协同分割模型能够通过结合序列图像中所有图像的前景信息来进行分割, 利用图像序列的全局信息对分割结果进行纠正, 从而得到正确的结果(图 6 (e)和图 7(e) ).
图 6 单帧不易分割的超声序列的分割结果Fig. 6 Segmentation results of breast ultrasound sequence hard to segment for single frame
图 7 单帧不易分割的超声序列的分割结果Fig. 7 Segmentation results of breast ultrasound sequence hard to segment for single frame表 1列出了文献[4, 28-29]中的方法与本文提出的方法在同一数据集上的整体表现. 从结果中可以看出, 文献 [4]中的方法有稍高的TPR (82.61 %)值, 这表明其模型分割得到的肿瘤区域与实际的肿瘤区域有较高的重合部分. 不过, 本文提出的模型有较低的FPR (17.54 %)值和较高的SIR(73.63 %)值, 这表明本文提出的模型能够有效地区分与肿瘤区域相似的正常组织区域.从边界评价标准上来看, 本文的方法也远低于文献[4, 28-29]中的结果, 这表明本文模型得到的分割结果更接近于实际肿瘤边界. 4类方法的受试者工作特征(Receiver operating characteristic curve, ROC)曲线见图 8, 表 2列出了4类方法的AUC (Area under roc curve)值, 从表 2中可以看出, 本文模型的AUC面积最大, 这表明本文模型有最好的分割效果.
表 1 4类方法在所有乳腺超声图像上的整体分割结果Table 1 Segmentation results of four methods on all breast ultrasound images表 2 4 类分割方法的AUC 值Table 2 Results of AUC of four segmentation methods3. 结论
本文提出了一种多域协同分割模型来对乳腺超声图像序列进行分割, 本文的主要贡献有: 1)通过将空域与频域信息结合, 既可以利用空域中肿瘤的姿态、位置与强度特征, 又可以通过频域来得到肿瘤的边缘特征, 使结果的准确性有较大的提升; 2)结合了乳腺肿瘤图像的特征来对现有协同分割框架中的能量项进行设计, 建立起乳腺超声图像序列模型, 有效地利用序列图像特点进行乳腺肿瘤分割.实验结果表明, 本文的方法可以准确地对超声图像序列进行分割.在未来的工作中, 如何更好地利用图像序列前景之间的相似性, 例如肿瘤在形状上的关系, 以及分割多肿瘤图像, 是今后进一步的研究方向.
-
表 1 仿真参数
Table 1 Simulation parameters
参数类型 参数名称 单位 值(范围) 弹目参数 ${V_p}$ m/s 2 300 ${V_e}$ m/s 2 700 $a_p^{\max}$ g 30 $a_e^{\max}$ g 12, 15 ${\tau _p}$ s 0.2 ${\tau _e}$ s 0.2 观测参数 T s 0.01 ${\sigma _\theta }$ mrad 5 ${\sigma _a }$ $\rm m/{\rm s}^2$ 1 场景参数 $r_0$ m 15 000 ${\phi_p}(0)$ rad 均匀分布($-\pi /18, \pi /18$) ${\phi_e}(0)$ rad $ > \pi /2$且满足碰撞三角形 ${u _p}(0)$ g 0 ${a_p}(0)$ g 0 ${a_e}(0)$ g $a_e^{\max}$ 目标机动方式 - 随机乒乓 估计器参数 $s_w$ $\rm g^2\rm{/Hz}$ 1 $t_{sw}$ s 2 $\Delta t$ s 0.1 初始状态 - ${{\hat {\pmb x}}_0} = {[0, 0, 0, 0]^{\rm T}}$ 初始估计误差 - ${{\tilde {\pmb x}}_0} = {[0, 0, a_e^{\max}, 0]^{\rm T}}$ 初始估计协方差 - ${{{\tilde P}}_0} = \left[{\begin{array}{*{20}{c}} 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0\\ 0 & 0 &{{{(a_e^{\max})}^2}}& 0\\ 0 & 0 & 0 & 0 \end{array}} \right]$ 
下载: 导出CSV
符号说明 $P$ 导弹 $E$ 目标 $\tau_p$, $\tau_e$ 导弹和目标控制系统的时间常数 $a_p^{\max}, a_e^{\max}$ 导弹和目标最大横向加速度 ${V_p}, {V_e}$ 导弹和目标的飞行速度 ${u_p}, {u_e}$ 导弹和目标的横向加速度指令 $r$ 弹目相对距离 ${t_{sw}}$ 目标模式切换时刻 $t$ 仿真时间 ${t_f}$ 终止时刻 g 重力加速度, $9.8\rm m/{\rm{s}^2}$ $m$ 目标的运动模式 ${m_1}, {m_2}$ 目标在模式切换时刻前后的运动模式 $\Delta m$ 目标运动模式改变量, $\Delta m = {m_2} - {m_1}$ $T$ 离散采样时间间隔 ${\sigma _\theta }$ 测角精度 ${\sigma _a}$ 导弹加速度测量精度 ${s_w}$ 目标指令加速度误差的功率谱密度 $\Delta t$ 目标运动模式辨识延迟 ${\tilde {\pmb x}}$ 状态估计误差 ${\pmb \xi }, {\Sigma}$ 状态估计误差的均值和方差 $\eta (t)$ ZEM估计误差 $\mu, {\sigma ^2}$ ZEM估计误差的均值和方差
下载: 导出CSV
-
[1] Shinar J, Turetsky V. Three-dimensional validation of an integrated estimation/guidance algorithm against randomly maneuvering targets. Journal of Guidance, Control, and Dynamics, 2014, 32(3):1034-1039 doi: 10.2514/1.40542 [2] Shinar J, Shima T. Nonorthodox guidance law development approach for intercepting maneuvering targets. Journal of Guidance, Control, and Dynamics, 2002, 25(4):658-666 doi: 10.2514/2.4960 [3] Shinar J, Turetsky V. Meeting the challenges of modern interceptor guidance by non-conventional approaches. In: Proceedings of the 17th Mediterranean Conference on Control and Automation. Thessaloniki, Greece: IEEE, 2009. 1563-1568 http://ieeexplore.ieee.org/document/5164770/ [4] Shinar J, Oshman Y, Turetsky V. On the need for integrated estimation/guidance design for hit-to-kill accuracy. In: Proceedings of the 2003 American Control Conference. Denver, CO, USA: IEEE, 2003, 1: 402-407 [5] Witsenhausen H S. Separation of estimation and control for discrete time system. Proceedings of the IEEE, 1971, 59(11):1557-1566 doi: 10.1109/PROC.1971.8488 [6] Alexandre P. Multiple Model Estimation and Detection for Adaptive Guidance of Hybrid Systems [Master dissertation], McGill University, Canada, 2004. [7] Glizer V Y, Turetsky V, Shinar J. Terminal cost distribution in discrete-time controlled system with disturbance and noise-corrupted state information. International Journal of Applied Mathematics, 2012, 42(1):52-59 http://d.old.wanfangdata.com.cn/OAPaper/oai_doaj-articles_9728187a2dae195516565edae35386af [8] Shinar J, Glizer V Y, Turetsky V. Distribution of terminal cost functional in continuous-time controlled system with noise-corrupted state information. In: Proceedings of the 27th Convention of Electrical & Electronics Engineers in Israel. Eilat, Israel: IEEE, 2011. 1-5 http://ieeexplore.ieee.org/xpls/icp.jsp?arnumber=6377008 [9] Shinar J, Glizer V Y, Turetsky V. Terminal state distribution of continuous-time system with random disturbance and noise-corrupted information. International Journal of Applied Mathematics, 2015, 45(2):77-84 [10] Glizer V Y, Turetsky V, Shinar J. Distribution of terminal cost functional in discrete-time controlled system with noise-corrupted state information. In: Proceedings of the 2012 IEEE 27th Convention of Electrical & Electronics Engineers in Israel (IEEEI). London, UK: IEEE, 2012. [11] Moldavskaya E, Shinar J. Distribution of the zero-effort miss distance estimation error in interception problems. Applied Mathematics and Computation, 2015, 269:217-231 doi: 10.1016/j.amc.2015.07.073 [12] Shinar J, Turetsky V, Oshman Y. New logic-based estimation/guidance algorithm for improved homing against randomly maneuvering targets. In: Proceedings of the 2014 AIAA Guidance, Navigation and Control Conference. Providence, Rhode Island: AIAA, 2004. 1-17 doi: 10.2514/6.2004-4886 [13] Dionne D, Michalska H, Shinar J, et al. Decision-directed adaptive estimation and guidance for an interception endgame. Journal of Guidance, Control, and Dynamics, 2015, 29(4):970-980 doi: 10.2514/1.16050 [14] Shinar J, Turetsky V, Oshman Y. Integrated estimation/guidance design approach for improved homing against randomly maneuvering targets. Journal of Guidance, Control, and Dynamics, 2007, 30(1):154-161 doi: 10.2514/1.22916 [15] 范红旗.主动寻的制导中机动目标运动模式辨识技术博士学位论文, 国防科学技术大学, 中国, 2008.Fan Hong-Qi. Technology on Maneuvering Target Motion Mode Identification in Active Homing Guidance [Ph. D. dissertation], National University of Defense Technology, China, 2008. [16] Author(s) Zhu Y L, Fan H Q, Fan J P, Lu Z Q, Fu Q. Target turning maneuver detection using high resolution doppler profile. IEEE Transactions on Aerospace and Electronic Systems, 2012, 48(1): 762-779 http://ieeexplore.ieee.org/xpl/articleDetails.jsp?arnumber=6129669 [17] Bizup D F, Brown D E. Maneuver detection using the radar range rate measurement. IEEE Transactions on Aerospace and Electronic Systems, 2004, 40(1):330-336 doi: 10.1109/TAES.2004.1292169 [18] Hughes E J, Leyland M. Target manoeuvre detection using radar glint. Electronics Letters, 1998, 34(17):1695-1696 doi: 10.1049/el:19981188 [19] Fan H Q, Zhu Y L, Fu Q. Impact of mode decision delay on estimation error for maneuvering target interception. IEEE Transactions on Aerospace and Electronic Systems, 2011, 47(1):702-711 doi: 10.1109/TAES.2011.5705700 [20] Turetsky V, Shinar J. Missile guidance laws based on pursuit-evasion game formulations. Automatica, 2003, 39(4):607-618 doi: 10.1016/S0005-1098(02)00273-X [21] 李运迁, 齐乃明.基于零控脱靶量的大气层内拦截弹制导律.宇航学报, 2010, 31(7):1768-1774 doi: 10.3873/j.issn.1000-1328.2010.07.011Li Yun-Qian, Qi Nai-Ming. A zero-effort miss distance-based guidance law for endoatmoshperic interceptor. Journal of Astronautics, 2010, 31(7):1768-1774 doi: 10.3873/j.issn.1000-1328.2010.07.011 [22] Blair W D, Watson G A, Kirubarajan T, Bar-Shalom Y. Benchmark for radar allocation and tracking in ECM. IEEE Transactions on Aerospace and Electronic Systems, 1998, 34(4):1097-1114 doi: 10.1109/7.722694 [23] Li X R, He C. Model-set design, choice, and comparison for multiple-model estimation. In: Proceedings of SPIE 3809 International Symposium on Optical Science, Engineering, and Instrumentation. Denver, CO, USA: SPIE, 1999. 501-513 http://spie.org/x648.xml?product_id=364047 [24] Baram Y, Sandell N. An information theoretic approach to dynamical systems modeling and identification. IEEE Transactions on Automatic Control, 1978, 23(1):61-66 doi: 10.1109/TAC.1978.1101690 期刊类型引用(2)
1. 张华扬，王铁超. 基于观测器的切换模糊系统的事件触发控制. 模糊系统与数学. 2019(01): 92-101 . 
百度学术2. 朱芳来，蔡明，郭胜辉. 离散切换系统观测器存在性讨论及降维观测器设计. 自动化学报. 2017(12): 2091-2099 . 本站查看其他类型引用(4)
-
下载:
计量
- 文章访问数: 2168
- HTML全文浏览量: 251
- PDF下载量: 836
- 被引次数: 6
