一种肝脏CTA序列图像的血管分割方法

摘要

本发明公开一种肝脏CTA序列图像的血管分割方法，首先对输入三维肝脏序列图像进行对比度增强和平滑噪声预处理，然后，采用OOF和OFA算法增强肝脏血管及其边界，细化血管中心；根据血管的几何结构，自动搜索血管中心线的种子点，并提取肝脏血管的中心线、构建肝脏血管树；最后，结合快速行进法初步分割肝脏血管并计算相应的血管和背景灰度直方图，采用图割算法实现肝脏血管的精确分割。本发明针对低对比度、强噪声和模糊边界的CTA序列图像，充分利用血管的几何形状和灰度信息，能有效地准确分割肝脏血管。本发明的肝脏CTA序列图像的血管分割方法能够推广到其他的三维血管分割中。

主权项

一种肝脏CTA序列图像的血管分割方法，其特征在于包括以下步骤：步骤1：对原始三维肝脏图像进行预处理，提高血管的对比度，平滑CTA图像噪声首先，分别调节3D肝脏图像的窗宽/窗位到[400～800]/[200～500]范围，并将灰度值归一化到[0,255]，得到图像I₀(x)，其中x为三维空间位置坐标，然后采用各向异性滤波算法去除噪声：$\left\{\begin{array}{c}\frac{\partial I(x,t)}{\partial t}=div(c\left(x,t\right)·▿I)\\ I(x,0)=I_0\left(x\right)\end{array}\right.$其中，t为扩散时间，I(x,t)为t时刻的图像，div(·)为散度操作，为图像梯度，c(x,t)为非负的单调递减传导系数：$c(x,t)=e^{-{\left(\right|▿I|/k)}^2}$其中，k为扩散门限，k为[50～80]的自然数；步骤2：采用多尺度OOF和OFA算法增强血管及其相应的边界，细化血管中心OOF滤波器通过最小化球面通量$f(x,r,\hat{\rho })=\frac{1}{4{\mathrm{\pi r}}^2}{\int }_{\partial B_r}\left(\right(▿\left(G_{{\sigma }_0}*I\left(x+h\right)\right)·\hat{\rho }\left)\hat{\rho }\right)·\hat{n}dA={\hat{\rho }}^T{Q}_{r,x}\hat{\rho }$寻找最佳的投影方向，构建血管增强滤波器，其中，*为卷积运算，B_r表示以x球心，r为半径的球表面，为球表面的外向单位法向量，为球表面上的位置坐标，dA为球表面无穷小的面积分量，是方差为σ₀的高斯函数，I(x)为经各向异性滤波后得到的图像，为方向向量，Q_r,x为3×3通量矩阵，其第i行、第j列元素为$q_{r,x}^{i,j}=\frac{1}{4{\mathrm{\pi r}}^2}{\partial }_{i,j}\left(G_{{\sigma }_0}\right(x)*b_r(x)*I(x\left)\right)$其中，i,j∈{1,2,3}，b_r(x)为以x为球心、r为半径的球体阶跃函数，即在该球体内，b_r(x)取1，否则，b_r(x)取0；通过对通量矩阵Q_r,_x的特征值λ₁(x,r)、λ₂(x,r)、λ₃(x,r)和特征向量e₁(x,r)、e₂(x,r)、e₃(x,r)分析，令|λ₁(x,r)|≥|λ₂(x,r)|≥|λ₃(x,r)|，设计OOF血管增强算法f_OOF(x,r)＝‑(λ₁(x,r)+λ₂(x,r))结合OFA算法$s(x,r,\hat{\rho })=\frac{1}{4{\mathrm{\pi r}}^2}{\int }_{\partial B_r}(▿(G_{{\sigma }_0}*I\left(x+h\right))·\hat{\rho })dA$构建血管边界检测算法$f_{OFA}(x,r)=\sqrt{s^2(x,r,e_1(x,r\left)\right)+s^2(x,r,e_2(x,r\left)\right)}$为提高计算效率，所述多尺度OOF和OFA均采用快速傅里叶变换及其反变换进行优化，并对血管中心进行细化，计算过程包括：(1)设置最小半径r_min、最大半径r_max和半径增量r_step，初始化r＝r_min；(2)计算通量矩阵Q_r,x的元素i,j∈{1,2,3}，$q_{r,x}^{i,j}={\mathrm{FFT}}^{-1}\left\{FFT\left(I\right)4{\mathrm{\pi ru}}_i{u}_j{e}^{-2{\left(\pi r\right|\left|u\right|\left|\sigma \right)}^2}\frac{1}{\left|\right|u|{|}^2}\left(cos\right(2\pi r\left|\right|u\left|\right|)-\frac{sin\left(2\pi r\right|\left|u\right|\left|\right)}{2\pi r\left|\right|u\left|\right|})\right\}$其中，FFT(·)和FFT^‑1(·)分别是快速傅里叶变换和反快速傅里叶变换，σ为高斯函数的方差，u为频域的位置矢量；(3)计算通量矩阵Q_r,x的特征值λ₁(x,r)、λ₂(x,r)、λ₃(x,r)和特征向量e₁(x,r)、e₂(x,r)、e₃(x,r)；(4)计算半径为r的f_OOF(x,r)f_OOF(x,r)＝‑(λ₁(x,r)+λ₂(x,r))(5)计算半径为r的f_OFA(x,r)$f_{OFA}(x,r)=\sqrt{s^2(x,r,e_1(x,r\left)\right)+s^2(x,r,e_2(x,r\left)\right)}$其中，(6)结合OOF和OFA算法，细化血管中心M(x,r)：M(x,r)＝max(0,f_OOF(x,r)‑f_OFA(x,r))(7)如果r＜r_max，r＝r+r_step，重复(2)至(6)，否则，按照如下公式计算多尺度f_OOF(x)，f_OFA(x)和M(x)响应：$f_{OOF}\left(x\right)=\underset{r\in R}{max}\left(f_{OOF}\right(x,r\left)\right)$$f_{OFA}\left(x\right)=\underset{r\in R}{max}\left(f_{OFA}\right(x,r\left)\right)$$M\left(x\right)=\underset{r\in R}{max}\left(M\right(x,r\left)\right)$其中，R为半径r的集合；获取以体素x为圆心的最优血管截面半径此时，对应的通量矩阵的特征向量e₃的方向即为体素x的切线方向步骤3：根据CTA图像血管特征，自动选取血管中心线种子点并提取相应的中心线，结合血管的几何结构，构建肝脏血管树血管中心线种子点位于血管截面中心，考虑到细化图像中，血管中心的灰度值最大，首先选取M(x)大于阈值t_seed的体素作为候选种子点，然后，以每个候选种子点为中心，构建边长为r_seed的立方体区域，在该区域内选择M(x)最大的体素作为种子点x_i；血管中心线由一系列血管截面中心点连接而成，运用以下步骤提取，即：(1)沿种子点x_i的正切线方向搜索邻域体素；(2)如果该邻域体素的$M\left(x_i^n\right)$值大于阈值t_line，且与切线方向的夹角小于θ_line，将该邻域体素标记为中心线体素，并作为新的种子点；(3)重复(2)，当无法检测到符合条件的邻域体素，或该邻域体素已被标记为中心线体素时，循环结束；(4)沿种子点x_i的负切线方向搜索邻域体素；(5)如果该邻域体素的值大于阈值t_line，且与与负切线方向的夹角小于θ_line，将该邻域体素标记为中心线体素，并作为新的种子点；(6)重复(5)，当无法检测到符合条件的邻域体素，或该邻域体素已被标记为中心线体素时，循环结束；采用上述方法，遍历所有种子点，实现血管中心线的提取；基于肝脏血管的几何结构，构建静脉期肝脏血管树：首先，统计中心线的体素个数，去除体素个数小于n_voxels的孤立中心线，然后，以中心线端点x_e为中心，在边长为15～20的立方体区域内搜索其他中心线体素，计算x_e与其他中心线体素x_j的连接概率$conn(x_e,x_j)=\exp (-1.8{\theta }_{tree})\frac{{\hat{r}}_e}{{\hat{r}}_e+d}$其中，d为x_e到x_j的距离，为以x_e为圆心的最优血管截面半径，θ_tree为x_e的切线方向与的夹角；当同时满足以下条件时，将x_e与最大conn(x_e,x_j)的体素x_j连接起来：(a)conn(x_e,x_j)大于t_conn；(b)θ_tree小于(c)d小于20；最后，根据血管的连通性，选择体素数目最多的2条中心线作为肝静脉和肝门静脉血管树；步骤4：结合快速行进法初步分割肝脏血管并计算相应的血管和背景灰度直方图，构建图割能量函数，实现肝脏血管的准确分割首先，以肝脏血管树作为初始目标区域，采用快速行进法初步分割肝脏血管，分别统计目标和背景的灰度直方图；然后，利用血管的灰度分布和边界信息，构建图割能量函数，实现肝脏血管的准确分割：$E=\alpha \underset{p\in P}{\Sigma }{R}_p\left(L_p\right)+\beta \underset{\{p,q\}\in N}{\Sigma }{B}_{pq}·\delta (L_p,L_q)$其中，α和β分别为区域项和边界项的权重系数；L_p为体素p的标签，L_p的取值为1和0，R_p为区域项，h_obj和h_bkj分别为归一化的初步分割目标和背景灰度直方图，P表示体素p的集合；B_pq为边界项，dist(p,q)为体素p到邻域体素q的距离；N为体素p的邻域体素的集合，最后，采用优化算法最小化该能量函数，实现肝脏血管的准确分割。