一种基于图像特征提取的PET图像和谐化方法

一种基于图像特征提取的pet图像和谐化方法
技术领域
1.本发明属于图像处理技术领域，尤其涉及一种基于图像特征提取的pet图像和谐化方法。


背景技术：

2.正电子发射成像(positron emission tomography,pet)技术是一种能够通过示踪剂在人体内的富集程度判断组织或器官的代谢程度的功能性成像手段。目前已广泛应用于各类肿瘤患者的诊断、分级和预后工作以及人体脑部研究中。由于pet成像分辨率较低，pet成像通常与其他成像技术结合以弥补这一缺陷，例如pet/ct、pet/mr。理论上，只需要采用恰当的示踪剂，pet图像可以展现各种组织和器官的代谢特性。因此，研究人员通过示踪同位素与化学物质不同组合开发出针对不同疾病或者不同靶点的示踪剂以寻求病灶部位在图像中的特异性展示。与此同时，为了提高pet图像的重建效率和质量，更多的重建方法也被开发并应用到新一代扫描设备中，例如有序子集期望最大化(ordered subset expectation maximization,osem)、点扩散函数(point-spread function,psf)等。
3.标准摄取值(standard uptake value，suv)作为pet图像的计量标准，已被证实在不同扫描设备、示踪剂和重建方法的基础上获取的图像之间，suv值会存在一定的差异。除此之外，还有来自扫描设置等方面的影响。这种差异会使不同来源的数据之间的可比性较低，并严重限制了pet图像多中心联合研究的数据来源范围。尽管可以通过严格遵守的统一的操作步骤来提高数据之间的可比性，但这将导致现有pet图像数据的利用率大大降低。
4.关于pet图像的和谐化算法的作用在于对不同来源的pet图像进行处理以提高图像之间的可比性。和谐化算法能够提高现有pet图像数据的利用率，减轻多中心联合研究在数据获取上的负担，并尽可能避免因扫描过程中的操作差异所导致的suv值误差。目前国际上已有并广泛使用的和谐化算法包括德国西门子公司(siemens)提出的eqpet和日本佳能(canon)公司提出的gi-pet等。
5.现有的和谐化算法主要存在以下缺点：
6.现有算法的使用限制
7.eqpet和gi-pet都是商用的和谐化算法，在算法的适用范围上会有一定的限制。eqpet是西门子的pet扫描仪中的内置算法，他只能通过扫描开始时的设置对其进行修改。除此之外，它的算法处理过程是包含在图像重建过程中的，并没有相应的展示。因此，eqpet目前只适用于西门子所设计的扫描仪所获取的数据。gi-pet能够适用于大多数pet扫描仪，但其被限制于日本本土范围内使用，并且目前已停止销售。
8.现有算法的数据适用限制
9.eqpet和gi-pet都是基于三维高斯滤波的和谐化算法，且和谐化分别是以欧洲核医学协会和日本核医学协会所提出的suv参考范围作为最终标准。由于滤波算法只能将单向地将图像suv值降低，因此上述现有算法并不适用于低于参考范围的数据集，从而导致该部分数据无法得到利用。


技术实现要素：

10.为解决上述技术问题，本发明提出了一种基于图像特征提取的pet图像和谐化方法，通过伽马分布参数估计和空间自回归模型估计提取多组图像的数值分布特征以及空间邻域关系，再根据以上结果模拟得到处于统一数值水平的模拟图像数据。本方法不仅能够帮助pet图像多中心研究提高现有pet图像数据的利用率，还能帮助解决临床应用中设备更新过渡期所遇到的问题。
11.为实现上述目的，本发明提供了一种基于图像特征提取的pet图像和谐化方法，包括如下步骤：
12.获得多组不同来源的pet图像，并对其进行预处理；
13.对预处理后的各组所述pet图像进行伽马分布参数估计，获得各组估计的伽马分布特征；
14.对各组提取伽马分布特征后的伽马残差进行空间自回归模型系数估计，获得空间邻域关系；
15.基于各组所述估计的伽马分布特征以及所述空间邻域关系进行图像模拟，获得处于统一数值水平的多组模拟图像数据。
16.可选地，对所述pet图像进行预处理的过程包括：对所述pet图像进行suv校正、勾画感兴趣区域、数据归一化处理和滤波处理。
17.可选地，对所述pet图像进行suv校正的方法为：
18.对所述pet图像进行数值上校正和衰减校正；
19.以体重为基础，对数值上和衰减校正后的所述pet图像进行标准化。
20.可选地，对所述pet图像进行suv校正的计算公式为：
21.suvbw＝(图像像素值
×
校正斜率+校正截距)
×
suv校正因子，
22.其中，
23.可选地，对所述pet图像勾画感兴趣区域的方法为：
24.将标准化后的所述pet图像转换成椭圆柱形的数据。
25.可选地，对所述pet图像进行数据归一化处理的方法为：
26.将所述椭圆柱形的数据转换为总体均值为一的形式。
27.可选地，对所述pet图像进行滤波处理的方法为：
28.利用二维伽马滤波算法，对数据归一化处理后的所述椭圆柱形的数据进行降噪处理。
29.可选地，对预处理后的各组所述pet图像进行伽马分布参数估计的最大似然目标方程为：
[0030][0031]
其中，x
ijk
为每个体素所服从的长度为n的随机生成的伽马分布数据,μ
ijk
和φ
ijk
是
根据每个体素所服从的分布计算的。
[0032]
可选地，对各组提取伽马分布特征后的伽马残差进行空间自回归模型系数估计，获得空间邻域关系的方法包括：
[0033]
基于空间自回归模型，获得邻域之间的线性关系；
[0034]
基于所述邻域之间的线性关系建立似然目标方程，估计对应邻域的系数；
[0035]
利用非线性加权最小二乘估计法，将目标方程最小化，获得包含对应邻域系数的邻域关系模型；
[0036]
将所述邻域关系模型应用在伽马估计残差中，提取对应的邻域信息并计算残差；
[0037]
若所述残差显示仍有残留信息，则调整所述邻域关系模型重新进行估计；若满足预设要求则保存为最终的邻域关系模型；
[0038]
基于满足预设要求的所述邻域关系模型，获得空间邻域关系。
[0039]
与现有技术相比，本发明具有如下优点和技术效果：
[0040]
本发明是基于图像特征提取的和谐化算法，因此不需要外部的参考范围作为和谐化标准，也没有类似eqpet、gi-pet的单向和谐化的局限。在pet多中心研究中，应用本算法可以将多个来源的pet数据进行和谐化，将其转换到一个统一且可比的水平来对疾病进行研究，而不再局限于统一来源并严格遵从统一操作流程的数据，提高现有pet数据在多中心研究中的利用率。
[0041]
目前在临床应用中，因院内设备更新导致已有患者在治疗后在新pet扫描设备上进行的扫描得到的病灶suv值高于治疗前在旧pet扫描设备上获取的病灶suv值。因此，医生难以依据两次扫描的pet数据判断患者的治疗效果，对患者的预后工作带来困难。本发明通过和谐化将两次扫描的pet图像转换至统一的数值水平，降低扫描设备设备、示踪剂、重建方法等因素的不同所导致的suv值差异，令两次扫描的结果变得可以互相比较，解决临床应用中因设备更新过渡期所遇到的问题。
附图说明
[0042]
构成本技术的一部分的附图用来提供对本技术的进一步理解，本技术的示意性实施例及其说明用于解释本技术，并不构成对本技术的不当限定。在附图中：
[0043]
图1为本发明实施例一的一种基于图像特征提取的pet图像和谐化方法流程示意图。
具体实施方式
[0044]
需要说明的是，在不冲突的情况下，本技术中的实施例及实施例中的特征可以相互组合。下面将参考附图并结合实施例来详细说明本技术。
[0045]
需要说明的是，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。
[0046]
实施例一
[0047]
如图1所示，本发明提供一种基于图像特征提取的pet图像和谐化方法，包括如下步骤：
[0048]
获得多组不同来源的pet图像，并对其进行预处理；
[0049]
对预处理后的各组所述pet图像进行伽马分布参数估计，获得各组估计的伽马分布特征；
[0050]
对各组提取伽马分布特征后的伽马残差进行空间自回归模型系数估计，获得空间邻域关系；
[0051]
基于各组所述估计的伽马分布特征以及所述空间邻域关系进行图像模拟，获得处于统一数值水平的多组模拟图像数据。
[0052]
本实施例以两台图像数值水平不同的pet扫描仪作为例子并分为主要的四个步骤按处理顺序分开讲述。需要注意的是，本发明能够在多种不同的扫描仪中获得的数据中对数据进行和谐化，这里仅以两台扫描设备为例。由于前三个步骤(预处理、伽马参数估计、空间自回归模型系数估计)是每组pet数据分开进行的，在步骤中的说明将会只有一组数据被提及。
[0053]
预处理主要分成四个部分：suv校正、勾画感兴趣区域(roi)、数据归一化处理、滤波处理。
[0054]
本实施例中所采用的suv值标准为标准化摄取值体重(standard uptake value body weight,suvbw)。suv校正的目的是对图像进行数值上的校正以及衰减校正，并以体重为基础对数值进行标准化，以保证数据之间的可比性。suv值的校正方法按照以下等式计算：
[0055]
suvbw＝(图像像素值
×
校正斜率+校正截距)
×
suv校正因子其中，
[0056]
suvbw＝(图像像素值
×
校正斜率+校正截距)
×
[0057]
suv校正因子勾画roi是指用同一个椭圆圈出病灶周围的人体区域，将pet图像中的人体形的数据转换为椭圆柱形的数据并输入到后续的计算步骤中。这一步骤的目的，一方面是为了排除图像中人体之外的区域对后续估计的影响，另一方面是通过减少需要进行估计的体素数量降低算法的计算量，从而提高算法的效率。
[0058]
数据的归一化处理是将椭圆柱形roi内的数据转换为总体均值为一的形式以提高算法的收敛速度。
[0059]
滤波处理则是采用二维伽马滤波算法对数据进行降噪处理以降低噪声对算法造成的影响。
[0060]
伽马分布参数估计
[0061]
在本实施例中，考虑使用参数化的伽马分布来代表pet图像的重建活动值并假设每个体素服从各自相互独立的伽马分布。设{z
ijk
,i＝1,
…
,n；j＝1,
…
,n；k＝1,
…
,k}为每个体素所对应的suv值。其中n为每张图像的像素长度，k为该组图像的张数。每个体素服从相互独立的伽马分布，即：
[0062]
e(z
ijk
)＝μ
ijk
,var(z
ijk
)＝μ
ijk
φ
ijk
[0063]
其中μ
ijk
和φ
ijk
的根据每个体素所服从的分布计算的。对于各点的伽马分布的参
数估计，选择使用极大伽马似然估计方法和牛顿更新过程来优化得到最佳参数。从伽马分布的概率密度函数中获取目标优化方程：
[0064][0065]
其中，x
ijk
为每个体素所服从的长度为n的随机生成的伽马分布数据。
[0066]
具体的优化过程如下所述：
[0067]
随机生成每个体素所服从的伽马分布
[0068]
假设{block}代表目标体素及其周围相邻一圈的体素，计算{block}的均值{mean(block)}和方差{variance(block)}。则目标体素所服从的伽马分布的形状参数和尺度参数就可以分别确定为{[mean(block)]2/variance(block)}和{variance(block)/mean(block)}。需要注意的是，这里用{block}计算得到的伽马分布理论上是{block}的均值所服从的分布。在后续步骤中，通过所估计的信息模拟得到的的suv值相较于图像校正后的suv值会有所偏差。因此，在生成随机分布后将会计算并保存{block}的均值与原值的比例用于后续的数据模拟。
[0069]
初始化目标参数μ
ijk
和φ
ijk
[0070]
目标参数的初始化值是根据体素各自所服从的随机生成的伽马分布{x
ijk
}计算得到的。μ
ijk
的初始化值为mean(x
ijk
)，φ
ijk
的初始化值为variance(x
ijk
)/mean(x
ijk
)。鉴于我们对数据进行了归一化处理，μ
ijk
的初始化值应该都处于1附近。
[0071]
最小化目标方程
[0072]
最小化目标方程的实现是基于牛顿更新步骤的优化算法，通过按照先保持φ不变更新μ，再保持μ不变更新φ的顺序不断迭代直至参数变化低于容差时确定最佳的参数选择。目标方程对μ和φ的一阶{g(
·
)}和二阶{h(
·
)}偏导如下所示：
[0073][0074][0075][0076]
[0077]
以μ的优化为例，牛顿更新步骤是通过来计算更新后的参数，其中μ1为更新后的参数，μ0为更新前的参数，λ则是引入的位于(0,1]之间的步长因子。λ的选择的标准为能够满足μ1大于零的最小值。容差的计算则是遵从当且仅当μ和φ的容差均小于10-5
时，迭代才会停止。
[0078]
此外，本实施例在计算伽马分布估计的残差时额外增加一个概率转换的步骤，以避免可能出现的残差偏斜的情况。残差的计算按照以下等式进行：
[0079][0080]
其中，φ-1
(
·
)表示标准高斯分布的逆累积分布函数，f(
·
)表示体素各自所服从的伽马分布的累积分布函数。
[0081]
空间自回归(spatial auto regressive,sar)模型参数估计
[0082]
本实施例通过空间自回归模型分析伽马残差的三维协方差，以获取体素及其邻域之间残存的空间邻域关系。sar模型指定了一组邻域之间的线性关系：
[0083][0084]
其中u(n)表示目标体素的伽马估计残差值，u(n-k)表示目标体素的相邻体素，且k＝(k1,k2,k3)。∈(n)表示一组方差为σ2的高斯白噪声，而θk则表示对应邻域的系数。对于邻域系数θk的估计，我们考虑采用一种基于似然的方法。该方法对应的似然目标方程为：
[0085][0086]
其中，表示包括目标体素在内的整体所包含的体素，则是表示不包括目标体素在内的邻域所包含的体素。为一组线性差分算子，满足表示sar过程的三维谱密度(在下面的式子中以f
θ
(λ)的形式简化表示)，其中p
θ
(λ)＝∑kθ
keiλ
·k为sar模型系数的三维离散傅里叶变换。将目标方程对θ求导并令导数等于零，即得到：
[0087][0088]
其中，索引非零系数的对应体素。通过帕斯瓦尔定理(parseval’s relation)，方程可以被转换为以下形式：
[0089][0090]
其中，是对协方差的样本估计，(l'-l|θ)＝∫e-iλ
·
(l-l
′
)fθ
(λ)dλ为谱密度f
θ
(λ)的逆傅里叶变换，其给出了三维模型的自协方差。
[0091]
利用非线性加权最小二乘估计进行将目标方程最小化即可得到包含对应邻域系数的邻域关系模型。将获得的邻域模型应用在伽马估计残差中提取其邻域关系并计算残差。通过残差分析确定邻域模型的有效性，若残差显示仍有残留信息，则调整邻域模型并重新进行估计；若满足要求则保存为最终的邻域模型用于数据模拟。
[0092]
数据模拟及数值和谐化
[0093]
数据模拟的起点是高斯白噪声，因此本算法能够通过同时生成多组白噪声用于模拟来实现模拟图像数据的区间估计。首先，将邻域关系应用于白噪声中：
[0094][0095]
其中，和分别表示快速傅里叶变换(fft)及其逆变换(ifft)。将添加了邻域关系的模拟数据作为伽马估计残差，结合估计得到的伽马分布特征信息即可模拟得到归一化的suv值。
[0096]zijk
＝f-1
(φ(u(n))|μ
ijk
,φ
ijk
)
[0097]
其中，f-1
(
·
)表示体素各自所服从的伽马分布的逆累积分布函数。需要注意的是，在应用上述估计的伽马特征信息之前，两组不同来源的数据所估计得到的μ
ijk
和φ
ijk
需要用直方图匹配的方法将其调整到近似的分布。
[0098]
将模拟数据从归一化中恢复需要额外乘以两个系数。第一个系数是原值在预处理以及计算分布时损失比例的倒数，用于恢复在预处理过程中的产生的数值损失。第二个系数是两个数据在归一化时的除数的均值，用于将其恢复至统一的数值水平。
[0099]
最终得到的两组数据在数值上将会是处于统一水平的，此时病灶区域的数值差异被认为是病灶之间真实差异，显著降低扫描设备、示踪剂、重建方法等因素的不同所带来的差异。
[0100]
本发明提供了一种基于图像特征提取的pet图像和谐化方法，本发明的核心思路是不以外部suv参考范围为标准，通过基于图像特征提取的算法对从不同pet设备、示踪剂、重建方法等获取的pet图像进行和谐化，以达到提高不同来源数据之间的可比性的目的。
[0101]
该方法是基于图像特征信息提取的算法，通过获取的特征信息模拟得到处于统一数值水平的和谐化图像。与基于滤波的算法对比，本发明能够更好地保留图像中的特征信息，尽可能避免因和谐化算法对图像信息造成模糊或改变。另一方面，本发明是以图像的数值特征作为标准的，因此不需要外部的标准参考范围参与。
[0102]
以上，仅为本技术较佳的具体实施方式，但本技术的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本技术揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本技术的保护范围之内。因此，本技术的保护范围应该以权利要求的保护范围为准。