基于粒子滤波的脊柱MRI神经束和神经节分割方法

### 基于粒子滤波的脊柱MRI神经束和神经节分割方法 #### 1. 引言 神经通路的映射和定位对于脊柱病理诊断、治疗规划以及图像引导干预至关重要。高分辨率MRI技术的发展使得我们能够观察到硬脊膜内的神经束，它们穿过椎间孔并离开椎管。神经束与液体和骨骼有良好的对比度，但与骨髓和肌肉的强度往往相似。手动分割神经和神经节既具有挑战性又非常耗时。 目前，自动分割起源于硬脊膜囊并离开椎管的脊柱神经束对于诊断和手术规划具有重要意义。然而，高分辨率脊髓造影MR图像中神经的强度、对比度、形状和方向的变化使得分割成为一项具有挑战性的任务。 传统的管状结构分割方法，如区域生长、主动轮廓等，在神经分割任务中存在诸多问题。区域生长方法容易出现泄漏问题，即分割会“泄漏”到附近的结构中；主动轮廓方法需要良好的初始化，并且容易陷入局部极小值。而一些提取中心线的方法需要大量的用户交互，且在低对比度结构中假阳性率较高。 本文提出了一种基于粒子滤波的自动跟踪方法用于神经分割。粒子滤波具有多假设、灵活的动态和多样的似然模型等优点，非常适合处理神经分割任务中低对比度和神经方向多变的问题。 #### 2. 方法 ##### 2.1 粒子滤波基础 粒子滤波是一种顺序蒙特卡罗跟踪方法。假设在跟踪算法的第t步，粒子（状态）表示为$h_t$，并且状态遵循一阶马尔可夫链，即$p(h_t|h_{1:t - 1}) = p(h_t|h_{t - 1})$，其中$h_{1:t - 1}$表示第t步的状态历史。设$z_t$是第t步的基于图像的观察值，并且在给定状态的情况下，不同时间点的观察值是独立的，即$p(z_{1:t}|h_{1:t}) = p(z_t|h_t) \cdot p(z_{1:t - 1}|h_{1:t - 1})$。 跟踪算法的每一步估计后验分布$p(h_t|z_{1:t})$，通过一组K个加权样本$\{h_t^{(k)}, w_t^{(k)}\}_{k = 1}^K$非参数地表示。可以证明$p(h_t|z_{1:t}) \propto p(z_t|h_t) \cdot p(h_t|z_{1:t - 1})$，即粒子分布取决于似然函数$p(z_t|h_t)$和先验项$p(h_t|z_{1:t - 1})$。 粒子滤波器通过从第$t - 1$步生成的样本集$\{h_{t - 1}^{(k)}, w_{t - 1}^{(k)}\}_{k = 1}^K$生成第t步的样本集$\{h_t^{(k)}, w_t^{(k)}\}_{k = 1}^K$来维护后验分布。具体来说，根据权重$\{w_{t - 1}^{(k)}\}$从$\{h_{t - 1}^{(k)}, w_{t - 1}^{(k)}\}$中抽取一个样本$\langle h_{t - 1}, w_{t - 1} \rangle$，并将其传播为$\langle h_t^{(k)}, w_t^{(k)} \rangle$。状态向量$h_t^{(k)}$从$p(h_t|h_{t - 1})$中采样，其权重$w_t^{(k)}$通过似然$p(z_t|h_t^{(k)})$对$w_{t - 1}$进行重新缩放得到。每一步对权重进行归一化，使其总和为1。 ##### 2.2 神经轨迹的粒子表示 我们将每个粒子建模为围绕3D中心线曲线的管状结构。中心线设计为贝塞尔曲线，通过引入半径函数形成管状结构。 n次贝塞尔曲线由$n + 1$个控制点定义。我们选择使用三次曲线，即$c(\tau) = (1 - \tau)^3p_0 + 3(1 - \tau)^2\tau p_1 + 3(1 - \tau)\tau^2p_2 + \tau^3p_3$，其中$\tau \in [0, 1]$是参数化变量。 半径函数$r(\cdot)$沿线段二次变化，同样通过贝塞尔曲线定义：$r(\tau) = (1 - \tau)^2r_0 + 2(1 - \tau)\tau r_1 + \tau^2r_2$。此外，我们还维护线段内的平均图像强度$\mu$。状态向量$h = (p_0, p_1, p_2, p_3, r_0, r_1, r_2, \mu)$完全描述了相应的线段。这种构造可以处理具有可变方向性、厚度和对比度的管状结构，如神经束和神经节。 ##### 2.3 动态模型 动态模型描述了如何从算法前一步生成的状态向量$h_{t - 1}$构建状态向量$h_t$，即采样概率$p(h_t|h_{t - 1})$。 - 我们将$h_t$的第一个中心线控制点设置为$h_{t - 1}$的最后一个控制点：$p_{0,t} = p_{3,t - 1}$，以确保轨迹的连续性。 - 将$p_{1,t}$放置在$(p_{2,t - 1}p_{3,t - 1})$线上。$p_{0,t}$和$p_{1,t}$之间的距离$\ell$从$[0, L]$中均匀抽取，其中$L$是算法的一个参数。形式上，$p_{1,t} = p_{0,t} + \ell\hat{n}$，其中$\hat{n}$是$(p_{2,t - 1}, p_