CN109242850A - 采用鲁棒局部各向异性聚类模型的磁共振图像修正方法 - Google Patents
采用鲁棒局部各向异性聚类模型的磁共振图像修正方法 Download PDFInfo
- Publication number
- CN109242850A CN109242850A CN201811184812.0A CN201811184812A CN109242850A CN 109242850 A CN109242850 A CN 109242850A CN 201811184812 A CN201811184812 A CN 201811184812A CN 109242850 A CN109242850 A CN 109242850A
- Authority
- CN
- China
- Prior art keywords
- image
- model
- pixel
- value
- expressed
- 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.)
- Pending
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
- G06T7/0012—Biomedical image inspection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
-
- 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
- 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]
Landscapes
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Data Mining & Analysis (AREA)
- General Physics & Mathematics (AREA)
- Physics & Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- Health & Medical Sciences (AREA)
- Evolutionary Computation (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Artificial Intelligence (AREA)
- Life Sciences & Earth Sciences (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Radiology & Medical Imaging (AREA)
- Quality & Reliability (AREA)
- Magnetic Resonance Imaging Apparatus (AREA)
Abstract
本发明提供了一种采用鲁棒局部各向异性聚类模型的磁共振图像修正方法,包括:获取磁共振图像;对所述图像各像素点构造各向异性权重来获得各像素点的更新灰度值;构建带偏移场的图像修正模型;构建能量模型,通过对所述能量模型最小化,获得修正后的真实图像。就鲁棒性、分割准确性、Jaccard相似性系数和Dice系数而言,与现有流行的模型进行比较和定量评估,证明了RLAC在实际和合成脑MR图像上有着更好的性能。
Description
技术领域
本发明涉及图像处理领域,具体而言,涉及一种采用鲁棒局部各向异性聚类模型的磁共振图像修正方法。
背景技术
在临床诊断中,医学图像已成为必不可少的辅助手段。而在成像过程中不可避免的存在噪声和偏移场,需要进行图像影像的修正,目前常见的修正方法都需要进行图像分割处理。
医学图像分割的主要目的之一是将图像分成不同的部分,每个部分内的像素具有相同的物理特征,如形状,纹理等。磁共振(MR)图像分割是MR图像处理和分析中第一步也是最重要的一步,其结果将直接影响诊断和治疗。
MR图像分割的准确性对临床诊断具有重要的指导意义。然而,噪声和偏移场是脑MR图像中错误分割的主要原因,目前几乎没有一个模型能够对所有的脑MR图像具有同样高的分割精度。由于成像设备的影响,MR图像中通常存在噪声,这种情况在成像过程中是无法避免。偏移场是由射频场的强度不均匀性引起的,对于相同组织的像素,它看起来是缓慢而平滑的变化。因此,MR图像分割过程中的这些缺陷将直接导致错误分类。
在脑MR图像分割中,聚类模型主要包括两类分别是硬聚类和软聚类。在硬聚类模型中,K均值聚类是它的典型方法。软聚类模型可以再次被分为混合模型和模糊c均值聚类(FCM)。FCM模型的主要思想是假设图像中的每个像素属于多个类别,但值得注意的是FCM并没有考虑空间信息。
为了增强FCM的鲁棒性,许多模型将空间信息融入到聚类模型中。例如偏差修正FCM(BCFCM)模型,该模型引入了邻域像素与其聚类中心之间的欧氏距离作为模糊聚类目标函数的约束项。然而,BCFCM不是各向异性的,即邻域中的每个像素具有与中心像素相同的权重。此外必须在每次迭代期间重新计算邻域和聚类中心之间的距离。
为了减少总耗时和复杂性,首先在迭代操作之前对图像进行滤波,从而提出了两种使用均值滤波和中值滤波的改进模型,即FCM_S1和FCM_S2。这两种模型都可以减少能量最小化的时间,但是均值滤波或中值滤波会不可避免的模糊图像的边缘信息。同时,随着迭代步数的增加,分割结果的模糊程度将逐渐增大。在EnFCM模型中,聚类方法是基于灰度信息而不是MR图像中的像素来实现的。而快速广义的FCM(FGFCM)模型用来增强EnFCM并定义了一个控制参数Sij来控制邻域的大小。因此,FGFCM的邻域受到灰度信息和空间信息的共同约束,但是它引入了太多的控制参数,使得该模型在实际应用中效率低下。
以上所有模型都没有考虑如何修正由强度不均匀性产生的偏移场。在修正后的图像中仍然存在噪声,对医疗诊断存在较为严重的影响,因此,迫切需要寻找一种较为完善的模型,对MR图像进行处理以获得理想的影像。
发明内容
本发明的目的在于提供一种采用鲁棒局部各向异性聚类模型的磁共振图像修正方法,能够解决上述提到的至少一个技术问题。具体方案如下:
一种采用鲁棒局部各向异性聚类模型的磁共振图像修正方法,包括如下步骤:
S100、获取磁共振图像;
S200、对所述图像各像素点构造各向异性权重来获得各像素点的更新灰度值;
S300、构建带偏移场的图像修正模型;
S400、构建能量模型,通过对所述能量模型最小化,获得修正后的真实图像。
可选的,所述步骤S200包括:
S202、选取图像中的任一像素点,划分出该像素点8邻域区域,计算各所述邻域区域内的点与其他邻域区域内的点的灰度均方差;
S204、将所述均方差通过高斯核函数投影到所述模型核空间,并计算所述均方差在所述核空间中的投影值;
S206、重复上述步骤S202、S204,获得所述图像中每个像素点在所述核空间的投影值,并得到相应像素点的权重值并将所述权重值归一化处理;
S208、将所述邻域区域中相似的点划分成一类,并根据该邻域内占有较大权重的类来更新当前像素点的灰度值。
可选的,所述步骤S300包括:
S302、将磁共振图像中的所有像素视为由偏移场和真实图像组成的乘积,则观察到的磁共振图像为
I(x)=b(x)J(x)+n(x) (1)
其中,I(x)是当前像素点x的观测值,b(x)是引起图像强度不均匀的偏移场,J(x)是当前像素点的真实值,n(x)是零均值的加性高斯噪声;
S304、根据偏移场平滑且缓慢变化的特性,偏移场b(x)表示为
其中,g1,...,gM为基函数,w1,…,wM是基函数的权重系数,M为自然数;
对偏移场b(x)转换表示为
b(x)=wTG(x) (3)
其中,w=(w1,…,wM)T,w为w1,…,wM的一个列向量,G(x)=(g1(x),…,gM(x))T;
S306、将真实的图像分成N个区域,真实图像J(x)表示为
其中,ci是常数,N为自然数,ui是隶属度矩阵。
可选的,所述步骤S202中灰度均方差表示为
其中,Ir和Ir'分别表示点r的灰度值和邻域内其他点的灰度值,nx表示邻域Nx中的像素点个数,
所述步骤S204中投影值表示为
其中,γxr为相应均方差在核空间中的投影值,
所述步骤S206得到所述权重值归一化处理结果表示为
所述步骤S208中当前像素点的灰度值表示为
可选的,所述步骤S400包括:
能量模型表示为
F(b,J)=∫Ω|AW(x)-b(x)J(x)|2dx (9)
结合公式(1)-(9),获得能量模型表示为
当上述公式(10)表示的能量模型被最小化之后,可以得到最优化的w、c和u,根据公式(4)得到分割之后的真实图像。
可选的,所述步骤S400中通过对所述能量模型最小化包括:对所述能量模型收敛迭代,设置迭代停止条件为|c(n)-c(n-1)|<ε,其中ε=0.0001。
本发明实施例的上述方案与现有技术相比,至少具有以下有益效果:本发明提出了一种新的鲁棒局部的各向异性(RLAC)聚类能量最小化模型,用于分割由噪声和强度不均匀性占据的脑部MR图像。首先构造各向异性加权方案,并利用图像块替换原始图像中的每个像素,可以将局部空间信息结合到能量方程中以提高鲁棒性。然后,使用由真实图像和偏移场乘积组成的乘法框架去修正偏移场和分割脑MR图像。为了确保偏移场平滑且缓慢变化的性质,使用平滑正交基函数的线性组合来表示偏移场。因此,组织分割和偏移场修正通过在能量最小化迭代过程中被同时执行。就鲁棒性、分割准确性、Jaccard相似性系数和Dice系数而言,与现有流行的模型进行比较和定量评估,证明了RLAC在实际和合成脑MR图像上有着更好的性能。
附图说明
此处的附图被并入说明书中并构成本说明书的一部分,示出了符合本发明的实施例,并与说明书一起用于解释本发明的原理。显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。在附图中:
图1示出了根据本发明实施例的图像修正方法流程图;
图2示出了根据本发明实施例的被提出的脑MR图像分割模型示意图;
图3示出了根据本发明实施例的带有噪声的不同组织的灰度值的示意图;
图4示出了根据本发明实施例的分割结果比较示意图;
图5示出了根据本发明的实施例的在真实图像中RLAC和BCFCM,FLICM,MICO模型比较的分割结果示意图;
图6示出了根据本发明的实施例的BCFCM,FLICM,MICO和RLAC分割结果的示意图;
图7示出了根据本发明的实施例的使用RLAC对具有严重强度不均匀性的合成脑MR图像进行分割结果的示意图;
图8示出了根据本发明的实施例的原始图像和偏移场修正图像的比较示意图;
图9示出了根据本发明的实施例的三种模型获得具有不同强度不均匀性的脑MR图像分割结果示意图;
图10示出了根据本发明的实施例的比较LIC,MICO和RLAC模型的结果示意图;
图11示出了根据本发明的实施例的LIC,MICO和RLAC模型分割结果的Dice系数示意图。
具体实施方式
为了使本发明的目的、技术方案和优点更加清楚,下面将结合附图对本发明作进一步地详细描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其它实施例,都属于本发明保护的范围。
在本发明实施例中使用的术语是仅仅出于描述特定实施例的目的,而非旨在限制本发明。在本发明实施例和所附权利要求书中所使用的单数形式的“一种”、“所述”和“该”也旨在包括多数形式,除非上下文清楚地表示其他含义,“多种”一般包含至少两种。
应当理解,本发明中使用的术语“和/或”仅仅是一种描述关联对象的关联关系,表示可以存在三种关系,例如,A和/或B,可以表示:单独存在A,同时存在A和B,单独存在B这三种情况。另外,本发明中字符“/”,一般表示前后关联对象是一种“或”的关系。
应当理解,尽管在本发明实施例中可能采用术语第一、第二、第三等来描述……,但这些……不应限于这些术语。这些术语仅用来将……区分开。例如,在不脱离本发明实施例范围的情况下,第一……也可以被称为第二……,类似地,第二……也可以被称为第一……。
取决于语境,如在此所使用的词语“如果”、“若”可以被解释成为“在……时”或“当……时”或“响应于确定”或“响应于检测”。类似地,取决于语境,短语“如果确定”或“如果检测(陈述的条件或事件)”可以被解释成为“当确定时”或“响应于确定”或“当检测(陈述的条件或事件)时”或“响应于检测(陈述的条件或事件)”。
还需要说明的是,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的商品或者装置不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种商品或者装置所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的商品或者装置中还存在另外的相同要素。
在传统的分割模型中,可以在不受噪声和偏移场影响的情况下获得准确的分割结果。然而,这些图像中的缺陷或多或少都存在于真实的脑MR图像中。因此,应用于真实脑部MR图像的传统模型几乎不能获得准确的分割结果。针对上述缺点,提出相关模型。
模糊的局部信息C均值聚类(FLICM)
FLICM不需要手动输入大量的参数,直接对图像引入模糊参数,模糊参数Gki的表达式为
其中,dkj表示第k个像素和第j个像素坐标的欧氏距离,uij表示隶属度矩阵,xj表示邻域内的像素点,vi表示聚类中心。在引入了模糊参数之后,FLICM的目标函数定义为
其中,uki的表达式如下
其中,vk的表达式如下
FLICM有令人满意的鲁棒性,但是在实际应用中FLICM并不总是有效的,在不同组织的边界处可能会产生误分割。这种情况发生的原因是相邻的组织的灰度值相差很大,在这两种组织的Gki会得到一个错误值。而此时如果第三种组织的灰度值接近这两种组织的中间值,那么Gki会使得相邻两种组织边界处的像素点划分到第三种组织。也就是说相邻两个组织边界处的像素点会被误分到第三种组织当中。
偏移场修正嵌入模糊C均值聚类模型(BCEFCM)
FLICM模型只能去除噪声,但无法校正偏置场。为了快速准确地修正偏移场,并对由噪声和强度不均匀性破坏的图像进行分割,提出了BCEFCM模型,在这个模型中引入了计算统一设备架构(CUDA)加速非局部均值去噪(NLMD)模型用以去除噪声像素。NLMD方法可以有效地保留原始信息,不会模糊图像的边缘。它们不使用给定像素的邻域信息,而是对整个图像域中的相似像素求平均值,并且估计值由下式给出
其中,I(y)是噪声图像,根据y到x的相似性,他们使用定义为权重的W(x,y)来替换带有加权平均值I(y)的原始像素,它可以用下面的等式写出
其中,h代表过滤系数。
由于NLMD的计算时间过长,它们将其纳入CUDA。最后引入了一系列正加权系数λ1,λ2,...,λN去将图像分成N类。能量最小化函数被定义为
BCEFCM模型首先使用NLMD模型计算加权系数,该模型通过CUDA加速。该模型可以准确地去除噪声MR图像中的噪声,但BCEFCM估计的偏移场不满足平滑和缓慢变化的特性。因此,迫切需要更有效的模型来分割具有噪声和偏置场的MR图像。
鲁棒局部各向异性聚类模型(RLAC)。
在许多模型中都没有考虑到局部信息,因此分割结果往往对高噪声敏感。鲁棒局部各向异性聚类模型(RLAC)能够分割具有噪声和强度不均匀性的脑MR图像。如图2所示,被提出的脑MR图像分割模型的组成,RLAC模型由两个主要部分组成:各向异性加权方案和偏移场修正模型。从图2第一行的各向异性加权方案可以看出,在由白质,灰质和脑脊髓液组成的前景中原始图像具有严重的噪声和强度不均匀性(左),因此原始图像的相应直方图中不具有明显的峰值。相比之下,使用图像补丁(右)的噪声去除图像比原始图像更好,但是由于强度不均匀性,它在相应的直方图中仍然没有良好分离的峰值。在图2最下一行中的偏移场修正模型可以有效地对无噪声图像进行分类并且修正由正交勒让德多项式的线性组合表示的偏移场。
各向异性权重的构造
考虑上文所述的缺点我们在本发明提出新的分割模型。在8邻域中,每个像素点对中心点的影响不同,和中心点有相似灰度值的像素点会占有较大的权重,而和中心点灰度值相差很大的像素点则会占有较小的权重。因此我们利用邻域内每个像素点的像素值和权重值来更新当前像素点的值。
对于图像中的每一个像素点x,我们考虑这个像素点的8邻域信息。假设其邻域内的点为r,则r∈Nx。首先计算出邻域内各点r与其他点的灰度均方差,通过均方差我们能够得到邻域内各点灰度值之间差异,计算均方差公式如下
其中,Ir和Ir'分别表示点r的灰度值和邻域内其他点的灰度值。nx表示邻域Nx中的像素点个数,我们将上式通过高斯核函数投影到核空间
其中,γxr为相应均方差在核空间中的投影值,通过计算了每个像素点在核空间的投影值能够得到权重值,最后各向异性权重可以通过将权重归一化得到
通过以上的公式计算,我们可以将邻域中较为相似的点划分成一类,并根据该邻域内占有较大权重的类来更新当前像素点的灰度值。最后,当前像素点的灰度值由以下公式给出
灰度值AW(x)被命名为图像补丁,它被用来代替原始图像中的每个像素。图3带有2%噪声的不同组织的灰度值,(a)-(c)脑MR图像中的部分白质和灰质组织;(d)相应的像素值;(e)相应的权重值。从图3中可以看出,第一行显示了由2%噪声占据的白质组织的灰度值和相应的权重,我们可以观察到灰度值为192和188的白质组织中噪声像素的权重倾向于0在使用各向异性加权方案后。图3中的下面一行显示了边界区域中的灰质组织,我们可以观察到具有灰度值57,61,68,124和126的像素被视为噪声像素,并且其权重也趋向于0,因此,在使用各向异性权重后,我们可以清楚地看到各向异性加权方案可以有效地消除噪声的影响。
带偏移场的图像模型
为了更准确地反映来自观察到的MR图像中的偏移场和真实图像,我们假设观察到的图像I是连续域中的乘法模型。我们可以将MR图像中的所有像素视为由偏移场和真实MR图像组成的乘积。观察到的MR图像被定义为
I(x)=b(x)J(x)+n(x) (12)
其中,I(x)是当前像素点x的观测值,b(x)是引起图像强度不均匀的偏移场,它的变化应该被假设为非常缓慢。J(x)是当前像素点的真实值,它在成像过程中反映了不同组织的内在物理性质.因此在理想情况下,对于同一组织内的像素点x应该取相同的值。n(x)是零均值的加性高斯噪声。
在本发明中,我们将图像I假设为一个函数:在一个连续域Ω中根据等式(12)将图像分解成两个乘法因子b、J和一个加性噪声n(x)。我们通过寻求b的最优解同时完成偏移场估计和图像分割,从而实现能量最小化。根据偏移场平滑且缓慢变化的特性,我们将其表示为以下形式
其中,g1,...,gM为基函数,在本发明,多阶正交勒让德多项式函数和最佳系数线性组合的乘积被用来表示偏移场。通常,我们使用四阶多项式来修正偏移场,即M=15。四阶基函数定义为g1(x)=1,g2(x)=x1, g6(x)=x2,g7(x)=x1x2, 其中x1和x2是图像I(x)的方向分量。因此,我们可以使用列向量来描述基函数,它可以被表示为
G(x)=(g1(x),…,gM(x))T (14)
其中,(·)T是转置运算符,w1,…,wM是基函数的权重系数。我们将w1,…,wM描述为一个列向量w=(w1,…,wM)T。从能量最小化得到的w的最优解可以被用来修正偏移场。因此偏移场b(x)可以用以下向量形式表示
b(x)=wTG(x) (15)
根据上文所述,我们假设一个真实的图像Ω可以被分成N种不同的组织。Ωi表示第i个组织所在的区域,相应的区域Ω1,…,ΩN形成Ω的一个分区,即且在第i个组织中,真实图像J(x)可以被近似看成常数ci。而每个组织Ωi可以用其隶属度矩阵ui来表示。在理想情况下,每个像素只被一种组织所包含,所以隶属度矩阵ui是一个二元函数,并且下列等式是一个分段常数。根据以上的描述,x∈Ωi,ui(x)=1或者是ui(x)=0。真实图像J可以被看成是
RLAC的能量表达式
基于上面,我们提出了偏移场估计和图像分割的能量表达式。考虑到同时估计偏移场和分割图像的问题,我们应用如下能量模型,当这个能量模型被最小化的时候就可以得到偏移场和图像分割的结果
F(b,J)=∫Ω|AW(x)-b(x)J(x)|2dx (17)
其中,AW(x)为经过各向异性权重处理后的像素值。为了使这个能量方程最小化,我们将上文所述b和J的表达式带入到这个能量方程当中。所以能量方程可以被重新写为关于三个自变量的方程,分别是w=(w1,…,wN)T,c=(c1,…,cN)T,u=(u1,…,uN)T,这个能量等式被表示为
因此,可以通过最小化能量方程得到w,c,u的结果来得到偏移场b和真实图像J。由于在每一个Ωi中ui是一个二元隶属度函数并且交换积分和求和的次序不改变能量方程,因此,我们可以将能量方程改写成如下形式
当上述能量方程被最小化之后,我们可以得到最优化的w,c和u,根据公式b(x)=wTG(x)可以得到偏移场,根据可以得到分割之后的真实图像。
能量最小化
能量最小化是通过交替迭代实现的,F(w,c,u)的能量最小化中关于三个变量的描述如下。
关于w的能量最小化
固定c和u,通过最小化能量方程得到w=(w1,…,wN)T的最优解,可以通过偏导数实现
令w的最优解由以下给出
通过上述等式我们可以得到w的最优解,由此我们可以得到估计偏移场b的表达式
关于c的能量最小化
固定w和u,通过最小化能量方程得到c=(c1,…,cN)T的最优解可以通过如下公式来实现
关于u的能量最小化
在本发明中我们所使用的隶属度矩阵是二元函数。所以固定w和c,通过最小化能量方程得到u=(u1,…,uN)T的最优解可以通过如下公式来实现
其中,αi(x)=|AW(x)-wTG(x)ci|2。
通过对上文各向异性权重结构和分段常数c、隶属度矩阵u和偏移场矢量w的描述,我们可以将本发明的流程描述如图1所示。
在本发明中,设置迭代停止条件为|c(n)-c(n-1)|<ε,其中ε=0.0001~0.01。
具体的,根据本发明的具体实施方式,采用鲁棒局部各向异性聚类模型的磁共振图像修正方法,包括如下步骤:
S100、获取磁共振图像;
S200、对所述图像各像素点构造各向异性权重来获得各像素点的更新灰度值;
可选的,所述步骤S200包括:
S202、选取图像中的任一像素点,划分出该像素点8邻域区域,计算各所述邻域区域内的点与其他邻域区域内的点的灰度均方差;
S204、将所述均方差通过高斯核函数投影到所述模型核空间,并计算所述均方差在所述核空间中的投影值;
所述步骤S204中投影值表示为
其中,γxr为相应均方差在核空间中的投影值,
S206、重复上述步骤S202、S204,获得所述图像中每个像素点在所述核空间的投影值,并得到相应像素点的权重值并将所述权重值归一化处理;
所述步骤S206得到所述权重值归一化处理结果表示为
S208、将所述邻域区域中相似的点划分成一类,并根据该邻域内占有较大权重的类来更新当前像素点的灰度值。
所述步骤S208中当前像素点的灰度值表示为
S300、构建带偏移场的图像修正模型;
可选的,所述步骤S300包括:
S302、将磁共振图像中的所有像素视为由偏移场和真实图像组成的乘积,则观察到的磁共振图像为
I(x)=b(x)J(x)+n(x)
其中,I(x)是当前像素点x的观测值,b(x)是引起图像强度不均匀的偏移场,J(x)是当前像素点的真实值,n(x)是零均值的加性高斯噪声;
S304、根据偏移场平滑且缓慢变化的特性,偏移场b(x)表示为
其中,g1,...,gM为基函数,w1,…,wM是基函数的权重系数,M为自然数;
对偏移场b(x)转换表示为
b(x)=wTG(x)
其中,w=(w1,…,wM)T,w为w1,…,wM的一个列向量,G(x)=(g1(x),…,gM(x))T;
S306、将真实的图像分成N个区域,真实图像J(x)表示为
其中,ci是常数,N为自然数,ui是隶属度矩阵。
可选的,所述步骤S202中灰度均方差表示为
其中,Ir和Ir'分别表示点r的灰度值和邻域内其他点的灰度值,nx表示邻域Nx中的像素点个数,
S400、构建能量模型,通过对所述能量模型最小化,获得修正后的真实图像。
可选的,所述步骤S400包括:
能量模型表示为
F(b,J)=∫Ω|AW(x)-b(x)J(x)|2dx
结合上述公式,获得能量模型表示为
当上述公式表示的能量模型被最小化之后,可以得到最优化的w、c和u,得到分割之后的真实图像。
可选的,所述步骤S400中通过对所述能量模型最小化包括:对所述能量模型收敛迭代,设置迭代停止条件为|c(n)-c(n-1)|<ε,其中ε=0.0001。
本发明提出了一种新的鲁棒局部的各向异性(RLAC)聚类能量最小化模型,用于分割由噪声和强度不均匀性占据的脑部MR图像。首先构造各向异性加权方案,并利用图像块替换原始图像中的每个像素,可以将局部空间信息结合到能量方程中以提高鲁棒性。然后,使用由真实图像和偏移场乘积组成的乘法框架去修正偏移场和分割脑MR图像。为了确保偏移场平滑且缓慢变化的性质,使用平滑正交基函数的线性组合来表示偏移场。因此,组织分割和偏移场修正通过在能量最小化迭代过程中被同时执行。就鲁棒性、分割准确性、Jaccard相似性系数和Dice系数而言,与现有流行的模型进行比较和定量评估,证明了RLAC在实际和合成脑MR图像上有着更好的性能。
为了证明RLAC的鲁棒性和偏移场修正的能力,首先使用仅由高噪声占据的脑MR图像来证明RLAC的鲁棒性。然后我们使用RLAC来修正带有严重强度不均匀脑MR图像上的偏移场。最后,我们在噪声和偏移场同时占据的脑MR图像将RLAC与其他相关的模型进行比较。所有的实验结果都和目前流行的分割模型定性比较并且在Core(TM)i5-73002.50GHz CPU,4GB RAM,and Windows 7的操作系统上进行实验,除非另有说明,否则在本发明使用四阶正交勒让德多项式函数来表示偏移场。
RLAC的鲁棒性
为了充分证明RLAC的鲁棒性,我们首先应用从BrainWeb获得的T1加权的合成脑MR图像,无偏场,1mm的切片厚度和9%噪声作为实验图像并与BCFCM,FLICM和MICO比较分割结果。从图4中可以看出,上面一行显示切片67图像,切片96图像被显示在下面一行。图4(b)显示了在BCFCM的分割结果中明显存在噪声像素,因为它使用的局部空间信息不是各向异性的,因此8邻域中的噪声像素将导致明显的错误分割。FLICM(图4(c))模型利用邻域中的空间欧几里德距离,这将容易在相邻的不同区域和角落中产生过分割,因为距离太大而不能在每次迭代过程中计算。MICO(图4(d))没有考虑局部空间信息,因此对噪声敏感,所以会对高噪声的脑MR组织产生错误的分割。因此,我们可以清楚地看到分割结果中存在许多异常值。本发明提出的RLAC利用8邻域中每个像素的加权来更新当前中心像素并考虑局部空间信息,因此它(图4(e))可以识别邻域中的噪声像素,可以减少噪声的影响。从上述实验结果来看,RLAC的鲁棒性略优于FLICM,并且明显优于BCFCM和MICO。
如上所述,RLAC在合成脑MR图像中比其他三个模型鲁棒性更强。我们将RLAC应用于不同角度的轴向切片,矢状切片和冠状切片带颅骨的真实图像,并将其与BCFCM,FLICM和MICO进行比较。如图5所示,上面一行是轴向切片图像,中间一行是矢状切片图像,下面一行是冠状切片图像。我们可以观察到在真实脑MR图像中白质和灰质相互交织,三个真实的脑MR图像被低噪声所占据,BCFCM模型(图5(b))没有清楚地区分相邻区域的白质和灰质。由于过分割,FLICM模型(图5(c))不能很好的区分背景和脑脊液,因此当背景与脑脊液相邻时,它将会产生一个不正确的分割结果。MICO的分割结果(图5(d))仍然存在明显的噪声。在真实的脑MR图像中,RLAC(图5(e))还可以清晰地将图像分为4类,并有效地去除图像中存在的噪声而不受头骨的影响。
为了更客观、更准确的比较四种模型的分割准确度,我们可以采用定量评价的方法去评价四种模型的分割结果。BrainWeb同时也提供标准分割结果,它可以被用来定量评估分割精度。Jaccard相似系数(JS)可以被用来衡量两个集合的相似度,它被定义为
其中,S1和S2是需要比较的分割结果,从每种模型获得的结果需要和GT进行比较。|·|表示相应区域中的像素点个数。JS∈[0,1],JS越大则表示分割结果更加的准确。我们应用不同等级的噪声和不均匀程度不同强度的20幅图像作为实验数据。图6显示了从BCFCM、FLICM、MICO和RLAC的分割结果得到的白质和灰质的JS值的箱形图。我们可以清楚地观察到RLAC与其他三个模型相比具有最高的JS值,并且RLAC的JS值中不存在异常值,这表明所提出的模型具有强稳定性。FLICM的JS值略低于RLAC的JS值,但是它仅考虑噪声的影响而不能消除图像中存在的偏移场。MICO的JS值相对较长且低于其他模型,因为它没有考虑图像中存在的噪声。表1显示了用噪声等级越来越高的脑MR图像作为实验数据来比较BCFCM,FLICM,MICO和RLAC模型分割结果的分割精度,很明显,四种模型的分割精度随着噪声的增加而降低,但RLAC的分割精度比其他三个模型更好。结果表明,本发明提出的RLAC不仅具有更好的抗噪声鲁棒性,而且具有比其他三种模型更高的分割精度。总之,视觉和定量评估都表明RLAC模型对不同噪声水平更具鲁棒性,并且无论是合成图像还是真实脑MR图像,RLAC都可以获得比其他被比较的模型更准确的结果。
表1
带有不同噪声的BCFCM,FLICM,MICO和RLAC模型的分割精度
我们应用轴向切片,矢状切片和冠状切片并且具有严重强度不均匀和0%噪声的合成脑MR图像用来证明RLAC修正偏移场的能力。从图7中,我们可以观察到最上面一行是轴向切片图像,中间一行是矢状切片图像,最下面一排是冠状切片图像。从第二列可以看出,估计的偏移场满足平滑并且缓慢变化的性质,因为我们使用具有清晰的偏移场表达式的乘法模型,并且还用四阶勒让德多项式的线性组合来表示偏移场。从第三列的偏移场修正图像中,我们可以观察到使用RLAC模型后每个组织具有相同的强度。如最后一列所示,可以观察到白质,灰质,脑脊液和背景被很好的区分为四类,因此对于高强度不均匀性所占据的脑MR图像,可以特别得观察到RLAC的优点。我们还比较了原始图像和偏移场修正图像的直方图,以证明图像质量的改善。如图8所示的直方图,我们可以观察到的原始像素中白质,灰质和脑脊液的强度彼此重叠,这意味着我们不太可能准确地将原始图像分类为四类。由于严重的强度不均匀性,在相应的直方图(左)中没有明确定义的代表不同组织的峰值,但在偏移场修正图像的直方图中存在明显的峰值(右)。从上述比较中我们可以观察到RLAC可以在强度不均匀的脑MR图像上获得更好的性能。然后我们使用只有强度不均匀的脑MR图像来比较LIC,MICO和RLAC模型的平均JS值,并在图9中显示定量分割结果。我们可以注意强度不均匀性的范围从20%到100%变化。在JS图中,我们可以观察到LIC,MICO和RLAC模型的平均JS值随着强度不均匀性水平的增加而降低,但RLAC的平均JS值略高于LIC和MICO的平均JS值。使用具有噪声和强度不均匀性的脑MR图像的LIC,MICO和RLAC的视觉对比结果将在下一部分中显示。
在实验中,我们对从BrainWeb得到的T1加权,厚度为1mm,噪声等级为7%,含80%偏移场的切片图像进行分割。并且和两种知名的偏移场估计模型即LIC和MICO进行比较。图10展示了三种模型的分割结果,第一、二、三行分别是RLAC、LIC和MICO。第一列为原始图像,第二、三、四列分别为估计的偏移场、偏移场修正图像和分割结果。从图10我们可以看出LIC模型估计的偏移场不符合偏移场是缓慢并且光滑变化的性质,MICO模型不能有效的去除噪声对分割结果的影响,所以不能得到令人满意的结果。本发明提出的RLAC使用四阶平滑正交基函数作为偏移场的乘法因子,因此偏移场满足缓慢和平滑变化的特性,并且还可以获得准确的分割结果。从实验结果可以直观的看到应用本发明模型分割之后的图像能够清晰地分辨出各个组织之间的边界,并且很好的消除了偏移场对分割结果造成的影响,具有令人满意的鲁棒性。
在视觉上观察到的分割结果不能完全证明RLAC的性能,我们用定量评价标准去衡量分割结果的准确度,Dice系数描述了实验结果和标准分割图像之间的差异,其定义公式为
其中,S1和S2是需要比较的两幅被分割图像,我们使用每种模型得到的结果和GT进行比较。|·|表示相应区域中的像素点个数。DC∈[0,1],DC越大则表示分割结果更加的准确。我们使用从BrainWeb获得的具有相同缺陷的20个图像并计算它们的DC以定量评估它们的性能。图11显示了从RLAC、LIC和MICO得到的白质,灰质和脑脊液DC值的箱形图。我们可以观察到RLAC与LIC和MICO相比具有更高的Dice系数。可以看出,本发明提出的RLAC模型在存在噪声和偏移场的情况下仍具有高精度。我们最终分别计算了从BrainWeb数据集上的LIC,MICO和RLAC模型获得的白质,灰质和脑脊液的分割精度和JS,并在表2中显示了分割结果。在表2中噪声范围为3%至7%,强度不均匀性范围为40%至80%。值得注意的是,分割结果与GT图像之间的更好匹配将会有更高的值不仅对于分割精度而言还对于JS同样适用。我们可以清楚地观察到,所提出的RLAC模型比LIC和MICO产生更好的分割精度和JS。
表2在不同水平噪
声和强度不均匀性的情况下,LIC,MICO和RLAC模型的分割精度和JS
本发明提出的一种新的鲁棒局部各向异性聚类模型用于同时分割组织和修正偏移场。我们构造了一种各向异性加权方案,将局部空间信息结合到能量函数中以提高鲁棒性。在该方案中,原始像素被图像块替换,并且还可以保留噪声图像中的细节。然后,我们将观察到的图像表示为真实图像和偏移场的乘积,它可以准确的分割带有噪声的脑MR图像并同时修正偏移场。我们还使用平滑正交基函数的线性组合来表示偏移场,以确保其平滑和缓慢变化的性质。然后基于图像块,偏移场和真实图像的表达来定义能量函数。RLAC的优点是它可以有效地消除MR图像中存在的噪声,并且它具有比FLICM,LIC和MICO更高的性能。与现有分割模型的比较实验已经证明了RLAC在实际和合成脑MR图像上的鲁棒性,偏移场修正,分割精度,JS和Dice系数等方面的优势。我们还成功地将RLAC应用于具有噪声和偏移场的轴向切片,矢状切片和冠状切片图像,并获得了高质量的分割结果。此外,从实验结果可以看出,RLAC可以同时轻松处理脑MR图像中不同程度的偏移场和噪声。在不久的将来,我们将进一步提高RLAC的稳健性和修正偏移场的能力。我们还尝试应用更高阶基函数作为偏移场的乘法因子来处理包含多重干扰的脑MR图像。
以上实施例仅用以说明本发明的技术方案,而非对其限制;尽管参照前述实施例对本发明进行了详细的说明,本领域的普通技术人员应当理解:其依然可以对前述各实施例所记载的技术方案进行修改,或者对其中部分技术特征进行等同替换;而这些修改或者替换,并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。
Claims (6)
1.一种采用鲁棒局部各向异性聚类模型的磁共振图像修正方法,其特征在于,包括:
S100、获取磁共振图像;
S200、对所述图像各像素点构造各向异性权重来获得各像素点的更新灰度值;
S300、构建带偏移场的图像修正模型;
S400、构建能量模型,通过对所述能量模型最小化,获得修正后的真实图像。
2.根据权利要求1所述的方法,其特征在于,所述步骤S200包括:
S202、选取图像中的任一像素点,划分出该像素点8邻域区域,计算各所述邻域区域内的点与其他邻域区域内的点的灰度均方差;
S204、将所述均方差通过高斯核函数投影到所述模型核空间,并计算所述均方差在所述核空间中的投影值;
S206、重复上述步骤S202、S204,获得所述图像中每个像素点在所述核空间的投影值,并得到相应像素点的权重值并将所述权重值归一化处理;
S208、将所述邻域区域中相似的点划分成一类,并根据该邻域内占有较大权重的类来更新当前像素点的灰度值。
3.根据权利要求2所述的方法,其特征在于,所述步骤S300包括:
S302、将磁共振图像中的所有像素视为由偏移场和真实图像组成的乘积,则观察到的磁共振图像为
I(x)=b(x)J(x)+n(x) (1)
其中,I(x)是当前像素点x的观测值,b(x)是引起图像强度不均匀的偏移场,J(x)是当前像素点的真实值,n(x)是零均值的加性高斯噪声;
S304、根据偏移场平滑且缓慢变化的特性,偏移场b(x)表示为
其中,g1,...,gM为基函数,w1,…,wM是基函数的权重系数,M为自然数;
对偏移场b(x)转换表示为
b(x)=wTG(x) (3)
其中,w=(w1,…,wM)T,w为w1,…,wM的一个列向量,G(x)=(g1(x),…,gM(x))T;
S306、将真实的图像分成N个区域,真实图像J(x)表示为
其中,ci是常数,N为自然数,ui是隶属度矩阵。
4.根据权利要求3所述的方法,其特征在于,
所述步骤S202中灰度均方差表示为
其中,Ir和Ir'分别表示点r的灰度值和邻域内其他点的灰度值,nx表示邻域Nx中的像素点个数,
所述步骤S204中投影值表示为
其中,γxr为相应均方差在核空间中的投影值,
所述步骤S206得到所述权重值归一化处理结果表示为
所述步骤S208中当前像素点的灰度值表示为
5.根据权利要求4所述的方法,其特征在于,所述步骤S400包括:能量模型表示为
F(b,J)=∫Ω|AW(x)-b(x)J(x)|2dx (9)
结合公式(1)-(9),获得能量模型表示为
当上述公式(10)表示的能量模型被最小化之后,可以得到最优化的w、c和u,根据公式(4)得到分割之后的真实图像。
6.根据权利要求5所述的方法,其特征在于,所述步骤S400中通过对所述能量模型最小化包括:对所述能量模型收敛迭代,设置迭代停止条件为|c(n)-c(n-1)|<ε,其中,ε=0.0001。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811184812.0A CN109242850A (zh) | 2018-10-11 | 2018-10-11 | 采用鲁棒局部各向异性聚类模型的磁共振图像修正方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811184812.0A CN109242850A (zh) | 2018-10-11 | 2018-10-11 | 采用鲁棒局部各向异性聚类模型的磁共振图像修正方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN109242850A true CN109242850A (zh) | 2019-01-18 |
Family
ID=65053229
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811184812.0A Pending CN109242850A (zh) | 2018-10-11 | 2018-10-11 | 采用鲁棒局部各向异性聚类模型的磁共振图像修正方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN109242850A (zh) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109886977A (zh) * | 2019-02-19 | 2019-06-14 | 闽南师范大学 | 一种带有邻域约束的图像分割方法、终端设备及存储介质 |
CN110232663A (zh) * | 2019-05-16 | 2019-09-13 | 福建自贸试验区厦门片区Manteia数据科技有限公司 | 一种危机器官自动勾画模型的修正方法及装置 |
CN113160213A (zh) * | 2021-05-12 | 2021-07-23 | 西南科技大学 | 一种融合图像局部和非局部信息的像素隶属度估计新方法 |
CN116908760A (zh) * | 2023-07-11 | 2023-10-20 | 哈尔滨工业大学 | 一种基于深度解耦网络的mri偏置场修正方法及系统 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106204592A (zh) * | 2016-07-12 | 2016-12-07 | 东北大学 | 一种基于局部灰度聚类特征的图像水平集分割方法 |
CN106204562A (zh) * | 2016-07-04 | 2016-12-07 | 西安交通大学 | 一种基于fMRI与DTI融合的穹隆白质分割的方法 |
CN107392866A (zh) * | 2017-07-07 | 2017-11-24 | 武汉科技大学 | 一种光照鲁棒的人脸图像局部纹理增强方法 |
CN107507209A (zh) * | 2017-08-07 | 2017-12-22 | 西安理工大学 | 极化sar图像的素描图提取方法 |
CN108272533A (zh) * | 2017-12-26 | 2018-07-13 | 中国科学院苏州生物医学工程技术研究所 | 皮肤创面区域的皮肤建模方法 |
-
2018
- 2018-10-11 CN CN201811184812.0A patent/CN109242850A/zh active Pending
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106204562A (zh) * | 2016-07-04 | 2016-12-07 | 西安交通大学 | 一种基于fMRI与DTI融合的穹隆白质分割的方法 |
CN106204592A (zh) * | 2016-07-12 | 2016-12-07 | 东北大学 | 一种基于局部灰度聚类特征的图像水平集分割方法 |
CN107392866A (zh) * | 2017-07-07 | 2017-11-24 | 武汉科技大学 | 一种光照鲁棒的人脸图像局部纹理增强方法 |
CN107507209A (zh) * | 2017-08-07 | 2017-12-22 | 西安理工大学 | 极化sar图像的素描图提取方法 |
CN108272533A (zh) * | 2017-12-26 | 2018-07-13 | 中国科学院苏州生物医学工程技术研究所 | 皮肤创面区域的皮肤建模方法 |
Non-Patent Citations (5)
Title |
---|
CHUNMING LI等: "Multiplicative intrinsic component optimization (MICO) for MRI bias field estimation and tissue segmentation", 《MAGNETIC RESONANCE IMAGING》 * |
杨红: "基于模糊聚类理论的脑部图像分割", 《中国优秀硕士学位论文全文数据库 信息科技辑》 * |
纪则轩: "基于模糊聚类的脑磁共振图像分割技术研究", 《中国博士学位论文全文数据库 信息科技辑》 * |
纪则轩等: "各向异性权重的模糊C均值聚类图像分割", 《计算机辅助设计与图形学学报》 * |
赵悦: "基于能量极小化方法的脑影像分割算法研究", 《中国博士学位论文全文数据库 医药卫生科技辑》 * |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109886977A (zh) * | 2019-02-19 | 2019-06-14 | 闽南师范大学 | 一种带有邻域约束的图像分割方法、终端设备及存储介质 |
CN110232663A (zh) * | 2019-05-16 | 2019-09-13 | 福建自贸试验区厦门片区Manteia数据科技有限公司 | 一种危机器官自动勾画模型的修正方法及装置 |
CN110232663B (zh) * | 2019-05-16 | 2021-04-30 | 福建自贸试验区厦门片区Manteia数据科技有限公司 | 一种危及器官自动勾画模型的修正方法及装置 |
CN113160213A (zh) * | 2021-05-12 | 2021-07-23 | 西南科技大学 | 一种融合图像局部和非局部信息的像素隶属度估计新方法 |
CN116908760A (zh) * | 2023-07-11 | 2023-10-20 | 哈尔滨工业大学 | 一种基于深度解耦网络的mri偏置场修正方法及系统 |
CN116908760B (zh) * | 2023-07-11 | 2024-05-14 | 哈尔滨工业大学 | 一种基于深度解耦网络的mri偏置场修正方法及系统 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109242850A (zh) | 采用鲁棒局部各向异性聚类模型的磁共振图像修正方法 | |
Ji et al. | A modified possibilistic fuzzy c-means clustering algorithm for bias field estimation and segmentation of brain MR image | |
US8073233B2 (en) | Image processor, microscope system, and area specifying program | |
Sijbers et al. | Automatic estimation of the noise variance from the histogram of a magnetic resonance image | |
CN107507173A (zh) | 一种全切片图像的无参考清晰度评估方法及系统 | |
JP2017521779A (ja) | 画像解析を用いた核のエッジの検出 | |
CN108364297B (zh) | 血管图像分割方法、终端、存储介质 | |
CN103218788B (zh) | 一种肝脏磁共振r2*参数的测量方法 | |
JP2013111420A (ja) | 画像処理装置、画像処理方法、及び画像処理プログラム | |
TWI624807B (zh) | 醫學影像之迭代式分析法 | |
CN111179295B (zh) | 改进的二维Otsu阈值图像分割方法和系统 | |
CN108550145A (zh) | 一种sar图像质量评估方法和装置 | |
CN103544695B (zh) | 一种高效的基于博弈框架的医学图像分割方法 | |
Galimzianova et al. | Stratified mixture modeling for segmentation of white-matter lesions in brain MR images | |
CN116416252B (zh) | 一种勃姆石生产工艺废水沉降图像检测方法 | |
WO2006087526A1 (en) | Apparatus and method for processing of specimen images for use in computer analysis thereof | |
Cheng et al. | A spatially constrained skew Student’st mixture model for brain MR image segmentation and bias field correction | |
Yusoff et al. | Automated two-dimensional K-means clustering algorithm for unsupervised image segmentation | |
Tung et al. | Identifying pathological slices of gastric cancer via deep learning | |
CN112396618A (zh) | 一种基于图像处理的晶界提取及晶粒度测量方法 | |
CN105956587B (zh) | 一种基于形状约束的膝关节磁共振图像序列半月板自动提取方法 | |
CN114862799B (zh) | 一种用于flair-mri序列的全自动脑体积分割方法 | |
CN108447066B (zh) | 胆道图像分割方法、终端、存储介质 | |
La Boissonière et al. | Atom based grain extraction and measurement of geometric properties | |
CN113538295A (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 | ||
WD01 | Invention patent application deemed withdrawn after publication | ||
WD01 | Invention patent application deemed withdrawn after publication |
Application publication date: 20190118 |