CN108921860B - 一种前列腺磁共振图像全自动分割方法 - Google Patents

一种前列腺磁共振图像全自动分割方法 Download PDF

Info

Publication number
CN108921860B
CN108921860B CN201810749052.7A CN201810749052A CN108921860B CN 108921860 B CN108921860 B CN 108921860B CN 201810749052 A CN201810749052 A CN 201810749052A CN 108921860 B CN108921860 B CN 108921860B
Authority
CN
China
Prior art keywords
equation
prostate
image
parameter
asymmetric
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810749052.7A
Other languages
English (en)
Other versions
CN108921860A (zh
Inventor
魏榕
朱翊
张珏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Peking University
Original Assignee
Peking University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Peking University filed Critical Peking University
Priority to CN201810749052.7A priority Critical patent/CN108921860B/zh
Publication of CN108921860A publication Critical patent/CN108921860A/zh
Application granted granted Critical
Publication of CN108921860B publication Critical patent/CN108921860B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/11—Region-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G06F18/232—Non-hierarchical techniques
    • G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • G06F18/23213—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions with fixed number of clusters, e.g. K-means clustering
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/20—Image enhancement or restoration using local operators
    • G06T5/30—Erosion or dilatation, e.g. thinning
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/70—Denoising; Smoothing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G06T2207/10088—Magnetic resonance imaging [MRI]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30004—Biomedical image processing
    • G06T2207/30081—Prostate

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Biology (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Probability & Statistics with Applications (AREA)
  • Image Analysis (AREA)
  • Image Processing (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)

Abstract

本发明提供了一种前列腺磁共振图像全自动分割方法，涉及医学图像处理领域，利用了关于前列腺解剖形状的先验知识，通过聚类算法提取DWI图像上前列腺的位置，再使用非对称的超椭圆模型对前列腺所在区域进行拟合，实现了一种前列腺磁共振图像全自动分割方法，相对传统的多图谱分割算法，本方法先验知识丰富，选取的特征较为显著，对于前列腺的分割能够具有更好的稳定性和精确性。

Description

一种前列腺磁共振图像全自动分割方法
技术领域
本发明涉及医学图像处理领域，尤其是一种多序列磁共振图像的医学图像分割方法。
背景技术
T2加权像(T2WI)主要反映的是组织问质子密度弛豫时间差别。在T2WI图像中，图像信噪比相对较高，能够反映前列腺的轮廓以及纹理特征。而弥散加权成像(DWI)图像用于衡量水分子受限扩散的难易程度。在DWI图像中前列腺和前列腺癌天然的显示为高信号，因此更适用于区分有严重病灶的前列腺图像和其他的组织。两个成像序列可以用于检测前列腺组织的异常，进行手术规划，辅助临床诊断有着重要意义。自动分割前列腺能够帮助衡量前列腺的大小，从而能够辅助医生进一步确定前列腺的病变程度。
然而，随着临床数据的急剧增加，对于医学图像的分割来说，手动分割已经满足不了需求，并且由于医生的经验以及主观因素，分割结果都会因人而异，因此对于感兴趣区域(ROI)的自动分割技术成为近年来医学图像领域的研究热点。
发明内容
为了克服现有技术的不足，本发明提出了一种前列腺磁共振图像全自动分割方法。
该方法参考了DWI图像中前列腺天然显示为高信号的特点，并通过自适应聚类将其作为前列腺分割的形状和位置先验。并且，受前列腺的解剖形状呈栗子型的外形特点启发，本文结合了DWI图像提供的先验星系，加上T2WI提供的边缘以及纹理信息，提出一种一种前列腺磁共振图像全自动分割方法，并通过匹配和拟合的策略进行前列腺分割。
本发明解决其技术问题所采用的技术方案具体包括以下步骤：
步骤一：首先要收集前列腺磁共振T2WI和DWI扫描数据，为了降低图像噪声的影响，其中的DWI图像中b值大于600的J0图像进行各项异性滤波前处理，各向异性滤波选取的迭代次数为N次，通常设定为10～30，为了平衡滤波结果和计算性能，本发明选取N＝15，传导系数为P，通常设定为20～100，本发明选用优选的P＝1，稳定系数为S，通常设定为0.00～0.25，，本发明选用优选的0.05，滤波完成后，得到去噪后的图像J1；
步骤二：对去噪后的DWI图像J1进行自适应聚类，考虑到前列腺腺体在J1中能够天然的显示为高亮信号，因此，得到的聚类中心灰度值最高的那个聚类结果作为前列腺初始外轮廓J2，其中，自适应聚类方法为k-means算法，选取的聚类特征为图像灰度值，聚类个数为K类，通常设定为3～15，根据经验，本发明中为优选的k＝4，此时聚类中心灰度值最低的那个类别对应图像中的空气，第二低的对应骨头和肌肉等组织，第三低的为膀胱等富含水的组织，第四低(最高)的即为前列腺组织，然后使用简单的形态学操作：填洞和腐蚀即可得到前列腺的初始轮廓J2；
步骤三：为了避免DWI低信噪比的导致J2存在锯齿状的边缘，利用非对称超椭圆模型对J2进行最小二乘拟合，这里使用的非堆成超椭圆模型可以由下列6个方程组成，分别为：标准椭圆参数方程、平移方程、旋转方程、凸化方程、凹化方程和非对称化方程；非对称超椭圆模型共含有9个待定参数，分别为：椭圆的半长轴参数ax,椭圆的半短轴参数ay,椭圆的形状参数ε,在x轴上的平移参数lx,在y轴上的平移参数ly,旋转角度r,变凸参数t,变凹参数b,非对称参数ω；θ为标准椭圆方程的控制参数，范围为0～2π，
标准椭圆参数方程：
Figure BDA0001725087100000011
平移方程：
Figure BDA0001725087100000021
旋转方程：
Figure BDA0001725087100000022
凸化方程：
Figure BDA0001725087100000023
凹化方程：
Figure BDA0001725087100000024
其中
Figure BDA0001725087100000025
非对称化方程：
Figure BDA0001725087100000026
得到分割后前列腺的外轮廓J3；
步骤四：根据分割后前列腺的外轮廓J3，按照医疗数位影像传输协定(DICOM)规定的内容，将DWI图像和T2WI图像根据头文件提供的物理坐标匹配起来，然后，将前列腺T2WI图像中外轮廓J3以外的区域置0，得到F0；
步骤五：考虑到DWI的图像信噪比相对T2WI较低，并且在DWI图像上区分出前列腺的外周带和移形带，因此，引入更加清晰的T2WI图像，考虑T2WI的边缘和纹理特征进行外周带的分割，这里使用3×3大小的Prewitt算子，遍历整幅图像F0，提取T2WI图像中的纹理信息，从而得到边缘特征图F1；
步骤六：通过局部二值法得到F0的局部纹理特征图F2，具体采用如下的方法：对于图像F0，选择3×3大小的窗口，以窗口中心像素灰度值为阈值，将8相邻的每个像素的灰度值与阈值进行比较，若大于阈值，则该像素点的位置被标记为1，否则为0；从左上角像素开始经顺时针比较可产生8位二进制数，并用这个8位二进制数来表示该点的纹理信息；使用局部二值法遍历整幅图像F0，得到F2；
步骤七：为了同时考虑边缘信息和纹理信息，这里将F1和F2加权求和得到前列腺外周带初始轮廓F3，其中，加权求和方式为F1+δ×F2＝F3，权重δ的变化范围为0.0～1.0，本文根据实际经验选δ＝0.6；
步骤八：利用公式(1)(2)(3)(4)(5)(6)定义的非对称超椭圆模型对F3进行最小二乘拟合，虽然和步骤三一样使用非对称超椭圆进行拟合，然而，此时由于基于不同的先验知识，所以选取的参数的求解范围稍有不同，因此能够得到符合前列腺外周带外观特征的分割结果F4；
本发明的有益效果是由于利用了关于前列腺解剖形状的先验知识，并且同时参考了T2WI和DWI图像的不同信息，实现了一种前列腺磁共振图像全自动分割方法，相对传统的多图谱分割算法，本方法先验知识丰富，选取的特征较为显著，对于前列腺的分割能够具有更好的稳定性和精确性。
附图说明
图1是本发明一种前列腺磁共振图像全自动分割方法的基本框图。
具体实施方式
下面结合附图和实施例对本发明进一步说明。
整个分割过程全部由软件MATLAB编写完成。
步骤一：首先要收集前列腺磁共振T2WI和DWI扫描数据，将收集到的DWI图像中b值大于600的J0图像进行各项异性滤波前处理，各向异性滤波选取的迭代次数为N次，通常设定为10～30，为了平衡滤波结果和计算性能，本发明选取N＝15，传导系数为P，通常设定为20～100，本发明选用优选的P＝1，稳定系数为S，通常设定为0.00～0.25，，本发明选用优选的0.05，滤波完成后，得到去噪后的图像J1；
步骤二：对去噪后的DWI图像J1进行自适应聚类，讲得到的聚类中心灰度值最高的那个聚类结果作为前列腺初始外轮廓J2，其中，自适应聚类方法为k-means算法，选取的聚类特征为图像灰度值，聚类个数为K类，通常设定为3～15，根据经验，本发明中为优选的k＝4，然后使用简单的形态学操作：填洞和腐蚀即可得到前列腺的初始轮廓J2，这里选取的腐蚀操作使用MATLAB提供的3×3的圆形作为腐蚀的结构元素；
步骤三：利用非对称超椭圆模型对J2进行最小二乘拟合，这里使用的非堆成超椭圆模型可以由下列6个方程组成，分别为：标准椭圆参数方程、平移方程、旋转方程、凸化方程、凹化方程和非对称化方程；非对称超椭圆模型共含有9个待定参数，分别为：椭圆的半长轴参数ax,椭圆的半短轴参数ay,椭圆的形状参数ε,在x轴上的平移参数lx,在y轴上的平移参数ly,旋转角度r,变凸参数t,变凹参数b,非对称参数ω；θ为标准椭圆方程的控制参数，范围为0～2π，
标准椭圆参数方程：
Figure BDA0001725087100000031
平移方程：
Figure BDA0001725087100000032
旋转方程：
Figure BDA0001725087100000033
凸化方程：
Figure BDA0001725087100000034
凹化方程：
Figure BDA0001725087100000035
其中
Figure BDA0001725087100000036
非对称化方程：
Figure BDA0001725087100000037
这里，为了加速超椭圆的拟合，根据解剖学的先验知识以及图像的实际大小，设定每个参数的求解范围分别为0≤ax≤40,0≤ay≤40,1≤ε≤2,20≤lx≤200,20≤ly≤200,-45≤r≤45,0≤t≤1,0≤b≤1,-1≤ω≤1，得到分割后前列腺的外轮廓J3；
步骤四：根据分割后前列腺的外轮廓J3，按照医疗数位影像传输协定(DICOM)规定的内容，读取DWI图像和T2WI图像头文件提供的Image Position参数，将DWI和T2WI根据这个参数的数值匹配起来，然后，将前列腺T2WI图像中外轮廓J3以外的区域置0，得到F0；
步骤五：使用3×3大小的Prewitt算子，遍历整幅图像F0，提取T2WI图像中的梯度信息，从而得到边缘特征图F1；
步骤六：通过局部二值法得到F0的局部纹理特征图F2，具体采用如下的方法：对于图像F0，选择3×3大小的窗口，以窗口中心像素灰度值为阈值，将8相邻的每个像素的灰度值与阈值进行比较，若大于阈值，则该像素点的位置被标记为1，否则为0；从左上角像素开始经顺时针比较可产生8位二进制数，并用这个8位二进制数，可以将这个八进制数换算程十进制数，即为一个0～255之间的数字，来表示该点的纹理信息；使用局部二值法遍历整幅图像F0，即可得到关于F0的纹理映射图F2；
步骤七：将F1和F2加权求和得到前列腺外周带初始轮廓F3，其中，加权求和方式为F1+δ×F2＝F3，权重δ的变化范围为0.0～1.0，本文根据实际经验选δ＝0.6；
步骤八：利用公式(1)(2)(3)(4)(5)(6)定义的非对称超椭圆模型对F3进行最小二乘拟合，虽然和步骤三一样使用非对称超椭圆进行拟合，然而，此时基于前列腺外周带类似马蹄形的先验知识，选取的参数的求解范围稍有不同，选取的参数范围分别为：0≤ax≤40,0≤ay≤40,1≤ε≤2,20≤lx≤200,20≤ly≤200,-45≤r≤45,-1≤t≤0,0≤b≤1,-1≤ω≤1，得到符合前列腺外周带外观特征的分割结果F4；
整个算法中，利用自适应的聚类以及前列腺的物理特征保证了前列腺初始轮廓的稳定性以及精确性，随后，利用超椭圆进行拟合，得到符合前列腺解剖形状特点且相对平滑的前列腺外轮廓，最后，提取T2WI图像中的边缘以及纹理信息进行超椭圆的拟合，得到符合前列腺外周带解剖信息的外周带分割结果。

Claims (1)

1.一种前列腺磁共振图像全自动分割方法，其特征在于，由以下步骤组成：
S1：对前列腺弥散加权成像(diffusion welghted imaging,DWI)图像中b值大于600的图像J0进行各项异性滤波进行前处理，得到去噪后的图像J1；
S2：对去噪后的图像J1使用K-means算法进行聚类，得到前列腺的初始外轮廓J2；
S3：根据由6个方程组和9个控制参数组成的非对称超椭圆模型M对初始外轮廓J2使用最小二乘法拟合匹配的参数，得到精细外轮廓J3；所述的非对称超椭圆模型M包含标准椭圆参数方程、平移方程、旋转方程、凸化方程、凹化方程和非对称化方程；所述的非对称超椭圆模型M包含椭圆的半长轴参数ax,椭圆的半短轴参数ay,椭圆的形状参数ε,在x轴上的平移参数lx,在y轴上的平移参数ly,旋转角度r,变凸参数t,变凹参数b,非对称参数ω；标准椭圆方程的控制参数θ；具体方程和参数如下：
标准椭圆参数方程：
Figure FDA0003175212690000011
平移方程：
Figure FDA0003175212690000012
旋转方程：
Figure FDA0003175212690000017
凸化方程：
Figure FDA0003175212690000013
凹化方程：
Figure FDA0003175212690000014
其中
Figure FDA0003175212690000015
非对称化方程：
Figure FDA0003175212690000016
S4:根据分割后前列腺的精细外轮廓J3，将前列腺T2加权成像(T2 weighted imaging,T2WI)图像中精细外轮廓J3以外的区域置0，得到图像F0；
S5：通过使用3×3大小的Prewitt算子，遍历整幅图像F0，得到边缘特征图F1；
S6:对于图像F0，选择3×3大小的窗口，以窗口中心像素点Q的灰度值为阈值，将8个相邻像素点Pk(k＝1,2,3,…,8)的灰度值与阈值进行比较，若大于阈值，则Pk(k＝1,2,3,…,8)所在的位置被标记为1，否则Pk(k＝1,2,3,…,8)所在的位置被标记为0；从左上角像素点开始经顺时针比较可产生8位二进制数，并用这个8位二进制数来表示Q的纹理信息；使用局部二值法遍历整幅图像F0，得到F2；
S7:将F1和F2加权求和得到前列腺外周带初始轮廓F3，其中，加权求和方式为F1+δ×F2＝F3，权重δ的变化范围为0.0～1.0；
S8：根据非对称超椭圆模型M对F3使用最小二乘法拟合匹配的参数，得到前列腺的外周带的分割结果F4。
CN201810749052.7A 2018-07-10 2018-07-10 一种前列腺磁共振图像全自动分割方法 Active CN108921860B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810749052.7A CN108921860B (zh) 2018-07-10 2018-07-10 一种前列腺磁共振图像全自动分割方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810749052.7A CN108921860B (zh) 2018-07-10 2018-07-10 一种前列腺磁共振图像全自动分割方法

Publications (2)

Publication Number Publication Date
CN108921860A CN108921860A (zh) 2018-11-30
CN108921860B true CN108921860B (zh) 2021-09-10

Family

ID=64411426

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810749052.7A Active CN108921860B (zh) 2018-07-10 2018-07-10 一种前列腺磁共振图像全自动分割方法

Country Status (1)

Country Link
CN (1) CN108921860B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115861600B (zh) * 2022-12-20 2023-09-05 西北民族大学 一种spect图像的roi区域识别方法及系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102737382A (zh) * 2012-06-22 2012-10-17 刘怡光 一种前列腺超声图像自动精确分割方法
CN104504720A (zh) * 2015-01-07 2015-04-08 四川大学 一种新的前列腺超声图像分割技术
CN106846349A (zh) * 2017-02-26 2017-06-13 哈尔滨理工大学 一种基于水平集的前列腺磁共振图像分割方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9269156B2 (en) * 2012-07-24 2016-02-23 Siemens Aktiengesellschaft Method and system for automatic prostate segmentation in magnetic resonance images

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102737382A (zh) * 2012-06-22 2012-10-17 刘怡光 一种前列腺超声图像自动精确分割方法
CN104504720A (zh) * 2015-01-07 2015-04-08 四川大学 一种新的前列腺超声图像分割技术
CN106846349A (zh) * 2017-02-26 2017-06-13 哈尔滨理工大学 一种基于水平集的前列腺磁共振图像分割方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
先验形状约束的超声前列腺图像分割方法;刘维平等;《微型电脑应用》;20120620(第06期);第5-9页 *
椭球先验约束的前列腺磁共振图像分割;李雪丽等;《南方医科大学学报》;20170331(第03期);第71-77页 *

Also Published As

Publication number Publication date
CN108921860A (zh) 2018-11-30

Similar Documents

Publication Publication Date Title
Chang et al. Segmentation of breast tumor in three-dimensional ultrasound images using three-dimensional discrete active contour model
CN107292312B (zh) 肿瘤ct图像处理方法
CN110120048B (zh) 结合改进U-Net和CMF的三维脑肿瘤图像分割方法
CN109753997B (zh) 一种ct图像中的肝脏肿瘤自动精确鲁棒分割方法
CN111667467B (zh) 基于聚类算法的下肢血管钙化指数多参数累积计算方法
WO2019000455A1 (zh) 图像分割的方法及系统
CN100561518C (zh) 基于感兴趣区域的自适应医学序列图像插值方法
CN102332161B (zh) 基于图像的血管内中膜厚度自动提取方法及系统
CN107292835B (zh) 一种眼底图像视网膜血管自动矢量化的方法及装置
CN112308846B (zh) 血管分割方法、装置及电子设备
WO2016191870A1 (en) Surface modeling of a segmented echogenic structure for detection and measurement of anatomical anomalies
CN109902682A (zh) 一种基于残差卷积神经网络的乳腺x线图像检测方法
US20030099384A1 (en) Detection and analysis of lesions in contact with a structural boundary
CN116503607B (zh) 一种基于深度学习的ct图像分割方法和系统
CN109191468B (zh) 一种血管提取的方法、装置及存储介质
CN115830016B (zh) 医学图像配准模型训练方法及设备
CN110782434A (zh) 一种脑结核mri图像病灶智能标记定位装置
CN111583286B (zh) 基于Flow-XDoG算子的腹部MRI图像轮廓提取方法
CN117893539B (zh) 一种乳腺图像识别处理方法
CN110610491B (zh) 一种腹部ct图像的肝脏肿瘤区域分割方法
Dabass et al. Effectiveness of region growing based segmentation technique for various medical images-a study
CN108921860B (zh) 一种前列腺磁共振图像全自动分割方法
JP2003515368A (ja) 画像が向上した視覚表示を得るための類似データのたたみこみフィルタリング。
Kitasaka et al. Automated extraction of aorta and pulmonary artery in mediastinum from 3D chest X-ray CT images without contrast medium
CN108629780B (zh) 基于颜色分解和阈值技术的舌图像分割方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant