一种针对脑动静脉畸形的磁共振图像语义分割方法

本发明属于图像语义分割，具体涉及一种针对脑动静脉畸形的磁共振图像语义分割方法。

背景技术：

1、脑动静脉畸形(avm)是脑血管畸形疾病中的最常见类型，人群患病率约为16～23/10万，大部分发病人群为10～40岁，发病人群呈现年轻化，患者会伴有动脉瘤风险、癫痫和神经损伤症状等，严重情况下会导致血管破裂，出现颅内出血，带来较高的致残率和死亡率。动静脉畸形的病理特征是，动静脉之间没有毛细血管进行过渡连接，从而形成一团血管直径和厚度不均匀的畸形血管，经畸形血管的血液绕过周围组织，直接快速地从动脉进入到静脉，使得动静脉内壁承受极大压力，导致血管壁变薄变脆，最终可能导致脑血管破裂。由于avm结构的复杂性，难以区分滋养动脉和引流静脉的位置，不便于后续病情诊断、手术方案定制和术中指导，因此需要生成血管的三维可视化模型，从而解放医生基于经验在脑海中构造avm几何模型的工作，将其更多的精力放在观测血管真实三维结构、分析病情严重程度、详细制定诊治方案和术后疗效的预测。

2、血管分割任务是获得准确血管三维模型重建的重要环节，但是由于血管影像像素强度的不均匀性和血管结构本身的复杂性，脑部血管影像的分割很难达到令人满意的分割效果，同时存在过分割和欠分割的问题，该任务仍然具有极大的挑战。目前国内外已有大量学者对脑血管分割开展了相关研究，这些研究主要分为两大类：分别是基于模型和基于数据的分割方法，基于模型的分割方法可以总结为统计学方法、形变模型方法和追踪方法，基于模型的方法能够使用图像的统计特征或概率分布直接进行图像的自动分割，无需额外的模型训练，但其对图像的前提假设较多，对于细节丰富的分割结果准确性低，且计算复杂度较高，导致执行时间较长；基于数据驱动的分割方法主要是基于机器学习的分割方法，通过大量的图像数据和标签进行神经网络模型的训练，利用训练好的网络能够实现图像快速的语义分割，但该方法需要收集且标记大量的医学图像数据，由于avm病历的罕见性以及avm本身磁共振图像的复杂性，人工标注大量图像以训练模型是难以实现的，因此基于数据的分割方法亟待解决样本数量少及标签难获取问题。

技术实现思路

1、为解决上述技术问题，本发明提供了一种针对脑动静脉畸形的磁共振图像语义分割方法，融合统计学方法和深度学习方法，得到一种脑部血管语义分割算法，针对脑动静脉畸形血管，解决现有分割技术中计算资源消耗大、计算时间长、训练数据少和分割效果差的问题。

2、本发明采用的技术方案为：一种针对脑动静脉畸形的磁共振图像语义分割方法，具体步骤如下：

3、s1、对脑动静脉畸形的磁共振血管成像mra图像进行预处理；

4、s2、基于步骤s1预处理后的磁共振血管成像mra图像，构建低阶统计学模型，得到像素的类别似然概率，基于最大似然概率计算初始分割结果；

5、s3、基于步骤s2的初始分割结果，使用马尔科夫随机场和反向多尺度滤波算子构建高阶统计学模型，得到精确分割标签集和模糊分割标签集；

6、s4、基于步骤s3得到的精确分割标签集和模糊分割标签集，融合领域自适应方法、对比学习方法和弱监督学习方法，对双分支网络模型进行训练，得到针对脑动静脉畸形血管分割的双分支网络模型，完成磁共振图像语义分割。

7、进一步地，所述步骤s1具体如下：

8、首先对初始磁共振血管成像mra图像使用伽马校正调整图像亮度和对比度，凸显血管目标区域；然后通过中值滤波对图像进行降噪，同时保留图像细节；再使用fsl颅骨去除算法提取出脑部图像；最后对图像进行min-max归一化，将图像强度放缩至0至255的范围中，消除图像之间的强度差异。

9、进一步地，所述步骤s2具体如下：

10、基于步骤s1得到的预处理后的磁共振血管成像mra图像，构建低阶统计学模型，进行初步的统计学分割，获得初始分割结果。

11、其中，所述低阶统计学模型即基于图像像素强度直方图，组合瑞利分布和高斯分布构成的有限混合模型。

12、首先运用磁共振血管成像mra图像的先验知识，将磁共振血管成像mra图像的像素强度直方图按照像素强度分为三个部分：低像素强度区域对应于脑脊液、骨骼和背景空气；中间强度区域对应于脑组织，包括灰质和白质，以及部分眼睛；高强度区域对应于血管和皮下组织。

13、然后基于所述有限混合模型，使用em算法对像素强度分布进行高精度拟合。

14、其中，所述有限混合模型的高强度区域分布函数用于计算血管类别似然函数，所述有限混合模型的其他分布函数则被组合计算出背景类别似然函数。

15、所述有限混合模型由一个瑞利分布fr(ys)和三个高斯分布fgi(ys)(i＝1,2,3)构成，其函数表达式如下：

16、

17、其中，wr和wgi(i＝1,2,3)分别表示控制各个分布函数对拟合效果的权重，且满足权重归一化约束ys表示像素点s的强度值，其瑞利分布和高斯分布的函数表达式如下：

18、

19、其中，权重参数和σi、μi和β分布参数通过直方图分布特征及k-means聚类算法进行初始化，并且使用最大期望em算法对像素分布进行拟合，经过设定轮次em算法的迭代，每一次迭代都在上一次迭代的基础上进行参数更新，最终获得对像素分布的准确近似。

20、再基于完成拟合的有限混合模型的分布函数计算得到各个像素点对血管类别xv和背景类别xb的似然概率p(ys|xs)，表达式如下：

21、

22、其中，xs∈{xv,xb}表示像素点s的可能分类类别，图像中所有像素的似然概率构成似然概率图，似然概率图中每个像素的值都是该像素的似然概率。

23、最后通过最大似然概率估计得到像素点s的类别作为低阶统计学模型的分类结果，即初始分割结果表达式如下：

24、

25、进一步地，所述步骤s3具体如下：

26、s31、使用马尔科夫随机场和反向多尺度滤波算子构建高阶统计学模型，将步骤s2得到的初始分割结果和似然概率图输入高阶统计学模型，得到精确分割标签集；

27、首先基于反向多尺度滤波算子，计算得到血管边缘强度响应，先将图像进行像素强度翻转，将背景设置为高强度区域，而血管设置为低强度区域：

28、i'(s)＝255-i(s)

29、其中，i(s)表示原始图像在s＝(x1,x2)t位置的像素强度值，t表示转置操作，i'(s)表示翻转后对应的像素强度值，将翻转处理后的图像进行多尺度滤波，得到其目标边缘的增强，其滤波响应的hessian矩阵h(s,σ)表达式如下：

30、

31、

32、其中，*表示卷积运算，σ表示高斯算子的尺度；g(s,σ)＝(2πσ2)-1exp(-sts/2σ2)表示尺度为σ的高斯核。

33、对于每个像素，hessian矩阵h(s,σ)能够被分解为两个特征值：|λ1,σ|≤|λ2,σ|，特征值关系用来区分不同的结构特征；则滤波算子响应值vσ表达式如下：

34、

35、其中，λρ,σ表示用于降低噪声敏感度的正则项，表达式如下：

36、

37、其中，κ表示预先设定的参数，用于控制正则强度，其值介于0和1之间。

38、当施加多个不同尺度的滤波算子时，选择最大响应作为像素点s的最终滤波响应值表示正实数集。

39、然后将反向多尺度滤波响应值与马尔科夫随机场融合，得到像素级势能函数，表达式如下：

40、e(xs,xr)＝wae1(xs,xr)+l*wbv(r)

41、其中，xs∈{xv,xb}表示中心像素点s可能的分类类别；xr表示基于低阶统计学模型获得的初始分割结果中，s像素点的空间26邻域η(s)像素点r处的分类结果，r∈η(s)；e1(xs,xr)表示中心像素点s的分类类别xs与其邻域像素分类类别xr之间的关系，若两者的标签类别一致则赋值-1，反之为1；l∈{-1,1}表示控制多尺度滤波响应v(r)对后续先验概率的抑制和增强能力，若xs＝xv则l＝1实现对后续先验概率的抑制，反之l＝-1实现对后续先验概率的增强；wa和wb分别表示权重，用于权衡马尔科夫随机场算子和多尺度滤波算子的影响。

42、将所有邻域相关的势能函数相加，得到中心像素点s的能量函数u(xs)，表达式如下：

43、

44、再基于hammersley-clifford定理，得到像素点s关于各个类别的先验概率p(xs)，表达式如下：

45、

46、其中，z表示归一化常数，通过加和所有xs类别的分子得到，用于将先验概率p(xs)归一化到(0，1)区间。

47、再基于上述低阶统计学模型和高阶统计学模型中得到的似然概率和先验概率，通过贝叶斯定理求得中心像素点s后验概率值，再通过最大后验概率算法，选择后验概率最大类别作为像素点s的最终分割结果表达式如下：

48、

49、其中，表示正比例关系；通过对比像素强度ys条件下血管类别xv和背景类别xb的后验概率大小，选择最大后验概率对应的类别标签作为最终分割结果对图像所有像素点进行最大后验概率估计得到精确分割标签集b。

50、s32、基于步骤s31中获得的精确分割标签集和步骤s2中获得的似然概率图，通过加权最大似然概率估计和布尔差运算，得到模糊分割标签集；

51、对图像中所有像素进行加权最大似然概率估计得到扩张的分割结果，即扩大的前景标签集a，表达式如下：

52、p(ys|xv)＞λ*p(ys|xb)

53、其中，λ表示加权权重，λ＜1则扩大血管类别区域、缩小背景类别区域。

54、然后将扩大的前景标签集a与步骤s31得到的精确分割标签集b进行布尔差运算，将布尔差运算结果(a-b)作为后续需要深度学习方法重新进行标签确定的模糊分割标签集。

55、进一步地，所述步骤s4具体如下：

56、s41、基于对比学习和领域自适应算法，使用精确分割标签集和公开数据集对双分支模型进行训练；

57、引入min-max归一化处理后的正常脑部磁共振血管成像mra图像作为领域自适应学习中的源域，该源域数据集来源于公开数据集，而avm数据集作为目标域，通过源域和目标域构建领域自适应框架。

58、再结合对比学习与领域自适应算法进行域不变特征的提取，通过源域和目标域数据集分别对两个双分支模型进行训练。

59、其中，所述双分支模型在u-net骨架网络的瓶颈层添加并行网络层，则瓶颈层由两个并行网络层构成，分别是独立层和共享层，共享层用于提取源域与目标域之间的域不变特征，两个双分支模型共同使用一层共享层进行域不变特征提取，独立层针对各自数据集的域特定特征进行提取，两个模型分别有各自的独立层。

60、在领域自适应学习过程中使用对比学习方法，设置独立层与共享层的输出互为负对，设置不同域下共享层的输出特征之间为正对，通过后续模型的训练最大化正对特征之间的相似度，最小化负对特征之间的相似度。

61、将对比损失函数、源域数据的损失函数和目标域数据的损失函数加权结合得到综合损失函数，以最小化综合损失函数为目标实现网络模型的参数更新，经过预先设定的训练次数之后输出训练过程中最优的双分支模型。

62、s42、基于步骤s41训练好的模型，使用迭代式伪标签更新方法，筛选高置信度伪标签来微调下一代的模型；

63、将目标域avm数据集的磁共振血管成像mra图像输入到步骤s41训练好的模型中进行预测分割，首先通过目标域的双分支网络模型，计算图像中各像素的类别可能性prob＝net(g)。

64、其中，net表示针对目标域数据集g训练的双分支网络模型。

65、在得到当前prob之后，根据上一轮迭代的不确定掩膜m和伪标签ps初始化当前迭代中的不确定掩膜mcur和伪标签pscur。

66、然后使用预设定的分割阈值进行高置信度像素标签的筛选，选择满足(prob＞阈值)的类别作为当前伪标签，并设置其不确定掩膜mcur中对应像素位置的值为0，表示该像素位置的标签确定，不满足条件的则设置不确定掩膜mcur对应像素值为1，表示该像素位置的标签不确定。

67、最后对比上一代和当前的伪标签，若伪标签类别发生变化，将不确定掩膜mcur对应像素值设置为1，表示该像素位置的标签不确定，得到最终的伪标签集pscur和不确定掩膜mcur，当前迭代的伪标签集和不确定掩膜作为下一轮迭代的输入条件ps和m。

68、s43、基于精确分割标签集和更新好的伪标签集，进行新一轮的模型微调，经过预设迭代次数的模型训练和伪标签更新后，得到脑动静脉畸形血管分割的深度学习模型；

69、其中，模型的微调采用弱监督学习方法，该弱监督学习方法即融合精确分割标签集和伪标签集进行模型的训练，且在损失函数的计算过程中将精确分割标签集和伪标签集作为金标准。

70、s44、将avm的磁共振血管成像mra图像输入训练好的针对目标域的双分支模型，完成图像的语义分割。

71、本发明的有益效果：本发明的方法首先对脑动静脉畸形的磁共振血管成像mra图像进行预处理，通过预处理后的图像构建低阶统计学模型，得到像素的类别似然概率，基于最大似然概率计算初始分割结果，再使用马尔科夫随机场和反向多尺度滤波算子构建高阶统计学模型，得到精细分割标签集和模糊分割标签集，最后融合领域自适应方法、对比学习方法和弱监督学习方法，对双分支网络模型进行训练，得到针对脑动静脉畸形血管分割的双分支网络模型，完成磁共振图像语义分割。本发明的方法融合统计学分割法和深度学习方法，实现病灶血管的高性能分割，解决了现有分割技术中计算资源消耗大、计算时间长、训练数据少和分割效果差的问题。与现有分割技术相比较，本发明的方法使用统计学分割方法初步获取精确分割标签集和模糊分割标签集，无需专业医生手动标记大量血管标签，降低了标注成本和要求，大大减轻了标注负担；同时结合对比学习方法和领域自适应方法解决了脑动静脉畸形病历数量少的问题，使用已有的正常血管数据集进行域不变特征的充分提取，提高模型的泛化能力和分割效果；另外基于伪标签更新的弱监督学习方法引入模糊区域信息至模型训练过程中，并且保证其伪标签的历史连续性和可靠性，提高了模型对难分割区域的特征理解。