用于医学图像重建参数寻优的客观图像质量评价方法
技术领域
本发明属于医学图像处理技术领域,尤其涉及一种用于医学图像重建参数寻优的客观图像质量评价方法。
背景技术
医学图像质量评价,长期以来一直被医生用于评价X-ray、计算机断层成像、核磁共振成像、超声成像等众多模态成像设备的性能优劣。这种评价通常是主观的,称之为主观图像质量评价。另外,单种成像模态还具有不同的重建算法,重建质量的好坏也需要医生进行主观评价,这在进行复杂的重建参数寻优过程中显得费时费力。在这些应用场景中,基于人工智能和机器学习的客观图像质量评价,因计算快速、自动化程度高而凸显了优势。
目前,针对自然图像的客观图像质量评价已经取得了丰硕的研究成果。这是因为自然图像的特点是内容丰富,即语义多种多样,可以有人物、风景、建筑物等;而且失真机理易于辨别,比如JPEG压缩失真具有明显的块效应,高斯模糊失真会产生视觉模糊现象等,据此可以抽象出特定的数学模型用来量化失真程度,从而反映客观质量。这些数学模型包括自然场景统计模型、人眼视觉系统模型等。客观图像质量评价按照对参考图像的依赖程度,可分为全参考、部分参考和无参考三类。部分参考实质是由全参考向无参考过渡的衍生品,这种过渡是一种趋势。一方面,因为精密光学仪器价格昂贵和传输带宽有限,限制了参考图像的获取和传输;另一方面,自然图像本身的复杂语义结构具有规律性,依靠自然场景统计和人眼视觉系统可以发现这种规律,进而形成先验知识,这为摆脱参考图像的束缚实现无参考评价提供了可能。
虽然面向自然图像的客观图像质量评价日趋成熟,但是这些理论和方法却不完全适用于医学图像的客观质量评价,而且面向医学图像的客观图像质量评价大都是无参考的。主要原因包括:一、医学图像的成像方法不同于自然图像,与通过自然光的反射或透射成像有本质区别,因此,自然场景统计模型失去了基础;二、医学图像的语义没有自然图像复杂,尤其是根据不同参数重建后的系列图像,它们的内容几乎一致,因此语义是相同的,失去了统计意义;三、医学图像的参考图像难以获取或者根本不存在,因为医学图像在前期大都以rawdata的形式存在,需要通过处理或重建得到可供显示的医学图像。综合来看,重建的医学图像不具有自然场景统计特性,语义单一,无参考图像。通常,重建算法需要根据重建参数生成一系列重建图像,通过辨别这一系列图像的质量等级便可得到最优参数。通过客观图像质量评价,解决具有相同语义的系列重建图像的质量等级标定问题,将是一种可行且有效的方法。
复杂的语义结构往往是从大量具有不同语义的图像集中训练得到的,但在评价具有相同语义的系列重建图像时却失去了有效性,不能形成用于质量客观评价所需要的先验知识。然而,值得注意的是,主观评价却不受此限制,人眼依然能分辨出系列图像的质量等级。由此可见,人类视觉系统不仅能形成来自不同语义的先验知识,还能形成来自相同语义的先验知识。后者的一种数学解释便是相同语义图像的自相似性,这种自相似性会因为图像受到不同的失真而产生不同程度的丢失。研究和量化自相似性的丢失程度,将有助于评价具有相同语义的系列重建图像的客观质量。
发明内容
本发明的目的是解决具有相同语义的系列重建图像的质量等级标定问题。提供一种用于医学图像重建参数寻优的客观图像质量评价方法。发明立足于人眼视觉系统,充分利用相同语义的图像具有自相似性的特点,提取自相似信息,具体操作过程中为了解决系列重建图像的质量等级标定问题,采用了循环置顶策略并结合质量因子排序方法,得出系列图像的质量等级,取得了良好的图像质量评价效果。
本发明采取的技术方案是:
首先,通过循环置顶的方法,将系列重建图像输入一个基于自相似信息对比的多输入全参考图像质量评价系统中;其次,多输入全参考图像质量评价系统将置顶图像视为失真图像,非置顶图像视为参考图像,得到置顶图像相对众多非置顶图像的一系列客观质量;接着,将获取的系列客观质量取平均值得到置顶图像的质量因子,以此类推得到其它置顶图像的质量因子;最后,通过冒泡排序算法对质量因子进行排序,得到系列重建图像的质量等级。
本发明解决其技术问题所采用的技术方案如下:
步骤(1).输入一组相同语义的重建图像,不失一般性假定有n幅,命名为Ii,i=1,2,…,n。
步骤(2).将重建图像按照命名下标排成双向循环队列,I1位于队头,In位于队尾,完成I1的置顶操作;
步骤(3).将置顶图像视为失真图像D,非置顶图像视为参考图像R,不失一般性地假定两幅图像的大小为M×N;
步骤(4).将步骤(3)的参考图像R和失真图像D输入评价系统,分别对其实施s个尺度的Daubechies小波变换,该小波通过调整函式和小波函式分别对输入信号进行低通滤波和高通滤波,并且对所得到的低通子带信号和高通子带信号继续进行滤波处理,实现多尺度分析。
调整函式φ(x)在多尺度分析的每一层都满足
其中,(a0,…,aN′-1)为有限长度为N′的实数数列,称作调整系数。
而用于高通滤波的小波函式ψ(x)用调整函式的线性组合表示:
其中,(b0,…,bM′-1)是有限长度为M′的实数数列,称作小波系数。
由于Daubechies小波具有正交性,调整系数满足即调整系数必须和位移偶数间隔后的调整系数互相垂直(内积为0),其中,m为位移次数,k为调整系数在数组中的位置,代表整数集,δm,0代表冲激响应函数,该函数仅在m=0时具有非零值。而小波系数满足bk′=(-1)k′aN′-1-k′,k′为小波系数在数组中的位置。
实际计算时根据Daubechies小波的消失动量A与调整系数的查找表设计低通滤波器,再根据正交性得到小波系数,设计高通滤波器。在第s次尺度分解中,原始信号LLs-1根据低通(L)和高通(H)的不同组合会得到四个子带,分别为LLs,LHs,HLs,HHs,LLs继续作为下一尺度分解的原始信号,其余三个子带具有方向性,依次为水平方向,垂直方向和对角方向,记为θ=0,π/2,3π/4。为表达方便起见,用CRs,θ和CDs,θ分别表示参考图像R和失真图像D的子带系数,由于θ=0,π/2,3π/4,共得到3s组与参考图像R和失真图像D对应的子带系数。
步骤(5).将步骤(4)得到的对应的子带CRs,θ和CDs,θ做协方差变换,得到各自的协方差矩阵和
其中,s为分解尺度,θ为子带方向。
由于不同尺度的子带大小为M/2s-1×N/2s-1,根据矩阵乘法规则得到不同子带的协方差矩阵大小为M/2s-1×M/2s-1;
步骤(6).分别对步骤(5)得到的故特征值分解,得到各自的特征向量和特征值:
其中,是对应参考图像R各子带的特征值,且 是与对应的特征向量;是对应失真图像D各子带的特征值,且 是与对应的特征向量,diag()是对角矩阵,T为矩阵转置操作。
步骤(7).利用特征向量空间的对称性,通过L2范数求得特征向量的自相似性
其中,||||是L2范数算符,对于p维向量x,
步骤(8).对步骤(7)得到的取均值,得到失真图像D各个子带对参考图像R的自相似性sss,θ,
步骤(9).对于同一方向不同尺度的sss,θ做聚合,得到同一方向上的自相似性ssθ,
其中,α1,α2和α3为聚合系数。
步骤(10).对于三个不同方向的ssθ做聚合,得到整体自相似性,
其中,β1,β2和β3为聚合系数。
步骤(11).并行执行n-1个步骤(4)-(10)所述基于自相似性对比(self-similarity comparison,SSC)的全参考图像质量评价算法,得到关于失真图像D的n-1个自相似性ssi(i=1,2,…,n-1)。
步骤(12).对通过步骤(11)得到的n-1个关于第一幅重建图像D1的整体自相似性ssi取均值,得到D1的质量因子Q1,
同理,通过双向队列循环移位,依次完成第2至第n幅重建图像D2,D3,…Dn的置顶操作,分别得到各自的质量因子Q2,Q3,…Qn。
步骤(13).根据步骤(12)得到n幅失真图像的质量因子,通过冒泡排序算法,将质量因子按照从大到小排列,并且记录相应的重建图像的编号,质量因子最大的重建图像具有最好的重建质量,因而其对应的重建参数便是最优参数。本发明的有益效果:
本发明充分利用相同语义的图像之间存在的自相似性度量图像质量,不需要额外提供参考图像,而是通过循环置顶的方法构建虚拟参考图像,从而可以利用全参考图像质量评价算法。而且,采用特征向量自相似算法可以很好地度量图像的自相似性。虽然算法中的多尺度小波变换和矩阵特征分解的计算复杂度相对较高,但本发明在单次循环置顶时可以并行计算自相似性,从而大大降低了计算复杂度。实验结果表明,基于本发明所提出方法的图像质量客观评价与主观评价具有很好的一致性。
附图说明
图1为本发明基于自相似性对比的全参考图像质量评价原理图。
图2为本发明利用循环置顶方法计算质量因子的结构框图。
具体实施方式
下面结合附图对本发明方法作进一步说明。
如图1和图2所示,用于医学图像重建参数寻优的客观图像质量评价方法,其具体实施步骤如下:
步骤(1).本实施例中,在Matlab环境下进行编程,输入通过磁敏感加权成像模态得到的n幅(本实施例中,n=5)重建图像。这些重建图像是根据回波时间不同重建得到的,图像大小为512×512,即M=N=512。
步骤(2).将重建图像按照命名下标排成双向循环队列,I1位于队头,In位于队尾,完成I1的置顶操作;
步骤(3).将置顶图像视为失真图像D,非置顶图像视为参考图像R,两幅图像的大小为M×N;
步骤(4).将步骤(3)的参考图像R和失真图像D输入评价系统,本实施例中,分别对其实施3个尺度的Daubechies小波变换。该小波通过调整函式和小波函式分别对输入信号进行低通滤波和高通滤波,并且对所得到的低通子带信号和高通子带信号继续进行滤波处理,实现多尺度分析。
调整函式φ(x)在多尺度分析的每一层都满足
其中,(a0,…,aN′-1)为有限长度为N′的实数数列,称作调整系数。
而用于高通滤波的小波函式ψ(x)用调整函式的线性组合表示:
其中,(b0,…,bM′-1)是有限长度为M′的实数数列,称作小波系数。
由于Daubechies小波具有正交性,调整系数满足即调整系数必须和位移偶数间隔后的调整系数互相垂直(内积为0),其中,m为位移次数,k为调整系数在数组中的位置,代表整数集,δm,0代表冲激响应函数,该函数仅在m=0时具有非零值。而小波系数满足bk′=(-1)k′aN′-1-k′,k′为小波系数在数组中的位置。
实际计算时根据Daubechies小波的消失动量A与调整系数的查找表设计低通滤波器,再根据正交性得到小波系数,设计高通滤波器。在第s次尺度分解中,原始信号LLs-1根据低通(L)和高通(H)的不同组合会得到四个子带,分别为LLs,LHs,HLs,HHs,LLs继续作为下一尺度分解的原始信号,其余三个子带具有方向性,依次为水平方向,垂直方向和对角方向,记为θ=0,π/2,3π/4。为表达方便起见,用CRs,θ和CDs,θ分别表示参考图像R和失真图像D的子带系数。本实施例中,由于s=1,2,3;θ=0,π/2,3π/4,因此共得到9组与参考图像R和失真图像D对应的子带系数。
步骤(5).将步骤(4)得到的对应的子带CRs,θ和CDs,θ做协方差变换,得到各自的协方差矩阵和
其中,s为分解尺度,θ为子带方向,且在本实施例中,θ=0,π/2,3π/4。
本实施例中,由于不同尺度的子带大小为512/2s-1×512/2s-1,根据矩阵乘法规则得到不同子带的协方差矩阵大小为512/2s-1×512/2s-1;
步骤(6).分别对步骤(5)得到的做特征值分解,得到各自的特征向量和特征值:
其中,是对应参考图像R各子带的特征值,且 是与对应的特征向里;是对应失真图像D各子带的特征值,且 是与对应的特征向量;diag()是对角矩阵,T为矩阵转置操作。
步骤(7).利用特征向量空间的对称性,通过L2范数求得特征向量的自相似性
其中,||||是L2范数算符,对于p维向量x,
步骤(8).对步骤(7)得到的取均值,得到失真图像D的各个子带对参考图像R的自相似性sss,θ,
步骤(9).对于同一方向不同尺度的sss,θ做聚合,得到同一方向上的自相似性ssθ,
其中,α1,α2和α3为聚合系数,本实施例中,α1=0.5,α2=0.3,α3=0.2。
步骤(10).对于三个不同方向的ssθ做聚合,得到整体自相似性,
其中,β1,β2和β3为聚合系数,本实施例中,β1=0.4,β2=0.4,β3=0.2。
步骤(11).本实施例中并行执行4个步骤(11)-(10)所述基于自相似性对比的全参考图像质量评价算法,得到关于失真图像D的4个自相似性ssi(i=1,2,…,4)。
步骤(12).对通过步骤(11)得到的4个关于第一幅重建图像D1的整体自相似性ssi取均值,得到D1的质量因子Q1,
同理,本实施例中通过双向队列循环移位,依次完成第2至第5幅重建图像D2,D3,…D5的置顶操作,分别得到各自的质量因子Q2,Q3,…Q5。
步骤(13).根据步骤(12)得到5幅失真图像的质量因子,通过冒泡排序算法,将质量因子按照从大到小排列,并且记录相应的重建图像的编号,质量因子最大的重建图像具有最好的重建质量,因而其对应的重建参数便是最优参数。