一种基于波形信号互相关的PET飞行时间估计方法

一种基于波形信号互相关的pet飞行时间估计方法
技术领域
1.本发明属于pet成像技术领域，具体涉及一种基于波形信号互相关的pet飞行时间估计方法。


背景技术：

2.正电子发射断层成像(positron emission tomography，pet)是一种基于核物理的医学成像技术，它通过放射性核素标记的示踪剂来检测人体内的生理化学活动，获得人体内组织中代谢、血流等相关生理指标，已经为各类疾病的检测和预防提供了有效依据。pet的工作原理是：将放射性核素标记的药物通过注射或口服的方式输入进病人体内，这些放射性药物会随着人体代谢、血流等过程在人体内形成特定的浓度分布，放射性核素持续发生衰变并发出正电子，在正电子与人体内负电子相遇时会发生湮灭反应，产生一对方向相反能量为511kev的伽马光子，发射的伽马光子被pet探测器探测得到，将探测得到的结果进行校正、重建的过程，即可得到放射性示踪剂在人体内的分布图像，以此进行医学诊断。
3.然而目前pet图像的空间分辨率约为3～4mm，远高于ct、mri等医学检测手段的分辨率；为了获得更高的空间分辨率，除了需要在pet探测器的硬件层面做出优化，也需要在重建算法层面做出优化。其中tof(time-of-flight)-pet是记录了光子飞行时间的pet探测方式，相较于普通的pet探测数据，tof-pet记录了伽马光子到达两个探测器的时间，然而tof时间的分辨率受到探测器物理性能的限制，tof时间的探测精度越高，则可以给重建算法添加越强的约束，从而能获得越好的重建效果。因此，提高tof的分辨率是tof-pet亟待解决的首要问题。
4.目前，通过探测所得的波形信息提取tof时间信息的方式主要有两种：前沿鉴别与恒比定时鉴别。前沿鉴别方式将每个脉冲的到达时间确定为信号的前沿越过设定阈值的时间，此方法需人工设定一个信号强度阈值。恒比定时鉴别方式与前沿鉴别方式相似，区别在于其设定的不是信号强度阈值，而是信号强度的百分比阈值，脉冲的前沿到达此脉冲峰值信号的设定百分比的时间则被确定为脉冲的到达时间。相较于前沿鉴别方式，恒比定时鉴别器消除了每个波形因峰值信号大小不同而带来的影响。
5.上述两种方式均源于最基础的信号处理技术，因为这两种方式在历史上很容易通过电子电路的方式进行实现，此类方法暴力地将探测器探测所得到的波形信息压缩为两个时间点信息，并直接通过两个时间点相减而求得tof信息，忽略了波形信号中随机变化可能对tof时间精度造成的影响。然而，随着快速波形数字化仪的投入使用，现在可以利用更先进的信号处理技术来检测波形中可能包含的更多的定时信息。文献[berg e,cherry s r.using convolutional neural networks to estimate time-of-flight from pet detector waveforms[j].physics in medicine&biology,2018,63(2):02lt01]提出使用卷积神经网络来估计tof时间，通过将波形信号数字化，将信号序列与tof时间作为神经网络的输入与标签进行训练，通过机器学习的模式来学习波形信号中潜在的时间信息，以提高tof时间分辨率。
[0006]
对于上述方法，无论是传统的鉴别方法还是深度学习方法，均是对于单个脉冲信号的单独处理，而忽视了tof信号来源于两个相似脉冲信号的事实；前述方法均没有对两个脉冲信号之间的关系作处理，导致可能会丢失一些信号间的时间信息。


技术实现要素：

[0007]
鉴于上述，本发明提供了一种基于波形信号互相关的pet飞行时间估计方法，能够有效地提高tof时间探测的分辨率。
[0008]
一种基于波形信号互相关的pet飞行时间估计方法，包括如下步骤：
[0009]
(1)利用一对探测器对点源发生的湮灭事件进行探测，进而对探测得到的符合波形进行预处理并保存；
[0010]
(2)对于任一组符合波形，对其中的两个波形序列进行互相关估计，得到互相关波形序列；
[0011]
(3)对所述互相关波形序列进行无偏归一化处理；
[0012]
(4)通过步骤(1)～(3)的操作得到大量样本，每组样本包括对应点源位置探测到的符合波形以及无偏归一化后的互相关波形序列，进而将所有样本分为训练集和测试集；
[0013]
(5)构建双通道卷积神经网络模型，其由两个并行的卷积神经网络h1和h2以及后端的全连接神经网络组成；
[0014]
(6)以训练集样本中的符合波形和互相关波形序列分别作为卷积神经网络h1和h2的输入，对应点源位置与探测器对连线中点的偏移作为标签，对模型进行训练；
[0015]
(7)将测试集样本输入至训练好的模型中，即可输出预测出对应光源发生湮灭事件的位置与探测器对连线中心的距离，进而实现对pet飞行时间的估计。
[0016]
进一步地，所述步骤(1)中使用pet探测系统中的一对探测器相向放置，在两探测器连线上等间隔均匀放置多个点源，利用两探测器在每个位置上探测点源发生的湮灭事件，得到多组符合波形，每一组符合波形包含两个探测器对同一湮灭事件所发出的成对伽马光子响应得到的两条波形序列，两波形序列的时延即代表pet飞行时间。
[0017]
进一步地，所述步骤(1)中对符合波形进行预处理包括归一化与裁剪，归一化处理使波形最高点高度均为1，以去除掉不同符合事件波形序列高度大小不同所带来的影响；裁剪处理是取符合波形的两个波形序列中各自的最大值，并以在时间轴上比较靠后的最大值点为基准，向前采样3ns，向后采样2ns，组成共5ns长度的波形，以此保证裁剪后所得到的波形中大多为波形上升沿。
[0018]
进一步地，所述步骤(2)中采用以下表达式对两个波形序列进行互相关估计：
[0019][0020]
其中：符合波形中的两个波形序列分别表示为x和y，为互相关波形序列中第m个元素值，x
n+m
为波形序列x中第n+m个元素值，yn为波形序列y中第n个元素值，m为自然数且1-n≤m≤n-1，n为波形序列的长度。
[0021]
进一步地，所述步骤(3)中采用以下表达式对互相关波形序列进行无偏归一化处理；
[0022][0023]
其中：为互相关波形序列中第m个元素值，为无偏归一化后的元素值，m为自然数且1-n≤m≤n-1，n为波形序列的长度。
[0024]
进一步地，所述卷积神经网络h1由5个卷积层c1～c5依此连接组成，其输入为符合波形；所述卷积神经网络h2由6个卷积层c6～c11依此连接组成，其输入为互相关波形序列；卷积层c1为二维卷积层，其中卷积核大小为2
×
2，用于将符合波形中的两条波形序列合并为一维序列数据，其余卷积层均为一维卷积层，其中c2、c4、c5、c7、c8、c10、c11中的卷积核大小为2，c3、c6、c9中的卷积核大小为3，所有卷积层的步长均设置为2，且每一卷积层输出均需经过激活函数leakyrelu处理；h1与h2的输出拼接后作为全连接神经网络，全连接神经网络的输出即为模型的预测结果。
[0025]
进一步地，所述步骤(6)对模型进行训练的过程如下：
[0026]
6.1初始化模型参数，包括可学习参数、卷积层参数、学习率、优化算法以及最大迭代次数；
[0027]
6.2将训练集样本中的符合波形和互相关波形序列分别输入至卷积神经网络h1和h2，再将h1与h2的输出拼接后输入至全连接神经网络，正向传播后即得到模型的预测结果；
[0028]
6.3计算模型预测结果与对应标签之间的损失函数，并根据损失函数采用梯度下降法对模型参数进行反向传播更新；
[0029]
6.4根据步骤6.2和6.3利用优化算法对模型参数不断进行迭代更新，直至损失函数最小化收敛或达到最大迭代次数，即训练完成。
[0030]
进一步地，所述损失函数采用均方误差(mse)函数。
[0031]
本发明通过信号互相关的方式将两波形信号先进行融合，以提取两波形间的联系；波形信号互相关的方法常用于两个相似但有时延的信号，用于检测两个信号的时延信息，而pet探测所得到的两个波形信号也是相似且时延的波形，tof时间信息即是波形的时延时间，因此本发明波形互相关方法能够更加准确地求出tof时间信息。此外，本发明采用深度学习方法能够自动地提取得到输入信号中的潜在信息，在使用深度学习方法的过程中，除去输入两波形信号外，同时输入两波形信号的互相关信号，即给予深度神经网络以更多的输入信息，并将指导两波形联系的互相关信号也输入至神经网络，以提高神经网络的学习能力。
附图说明
[0032]
图1为本发明pet探测器对及点源的分布示意图。
[0033]
图2(a)为点源放置在80mm处时探测器所探测到的一组波形示意图。
[0034]
图2(b)为点源放置在-60mm处时探测器所探测到的一组波形示意图。
[0035]
图3(a)为点源放置在80mm处时探测到的波形经互相关得到的波形序列图。
[0036]
图3(b)为点源放置在-60mm处时探测到的波形经互相关得到的波形序列图。
[0037]
图4为本发明双通道卷积神经网络的结构示意图。
具体实施方式
[0038]
为了更为具体地描述本发明，下面结合附图及具体实施方式对本发明的技术方案进行详细说明。
[0039]
本实施例中我们使用模拟的方式来得到所需的波形数据，实验使用geant4蒙特卡罗应用软件包与gossip模拟软件包结合来对探测器所产生的波形进行模拟，geant4主要用于模拟放射源放射出粒子输运的物理过程，gossip主要用于模拟探测器后端的硅光电倍增管sipm对于可见光子的响应并产生出实验所需的波形。
[0040]
如图1所示，整个模拟过程中设置两个相同探测器相向放置，两探测器前表面的距离为400mm，即探测器系统半径为200mm，探测器轴向的长度与宽度均为50mm，每个探测器由50
×
50个lyso晶体组成的阵列与后端的硅光电倍增管sipm构成。
[0041]
首先将放射性点源放置在两探测器中心点的连线上，固定点源位置，在此位置上进行模拟测量得到多组符合波形信号，每个信号的起始位置均为此次模拟过程的开始时间点，这些信号被成对地存储在700
×
2的矩阵中，每个波形的长度为700，代表从模拟开始后35ns内探测器的响应，即波形矩阵中每两个元素的间隔为0.05ns，这对应数字示波器20gs/s的采样速率。随后移动点源的位置，移动方式为在连线上每隔1mm移动一次，移动范围为[-100mm,100mm]，即共有201个点源位置，在每一个点源位置上均采样约50组数据。
[0042]
得到采样的原始数据后，需要对数据进行预处理，由于在波形上与tof相关的信息在波形的上升沿过程中。而对于lyso晶体而言，其上升沿时间往往在2～3ns左右，故而探测的波形数据在宽度上远大于上升沿持续时间，这对于后续的时间检测不利。我们分别求取每一组符合波形中两个波形的最大值对应位置，并且以两个位置中的较大位置为基准位置，以此基准位置向波形左侧选取60个元素(对应3ns)，向波形右侧选取39个元素(对应2ns)，与基准位置共同组成连续具有100个元素的序列，对应时长为5ns，这样即可使我们所关注的波形更多为上升沿波形。
[0043]
此外，探测器在探测过程中会产生随机过程，导致实际上的探测器产生的各个波形的幅值有较大差异，这也会对后续对波形信息的提取造成影响，故需对所有波形进行归一化处理，即波形中每一个元素除以此波形序列中最大元素值，使得最大元素值归为1。
[0044]
如图2(a)和图2(b)所示，分别为点源放置在80mm与-60mm位置上时的经过预处理后的波形，点源位置不同导致两波形的起始位置不同，两波形的时延即代表tof时间。
[0045]
而后则需要对预处理后的波形进行互相关，得到互相关波形，计算表达式如下：
[0046][0047]
表达式中的x与y即为上述经过预处理后的两个波形序列，其长度n均为100，m的取值范围为[-99,99]，则互相关序列的长度为199。
[0048]
而后对互相关序列进行无偏的归一化估计，计算表达式如下：
[0049][0050]
图3(a)和图3(b)分别为两组符合波形所对应的经无偏归一化后的互相关序列，可见两序列具有比较平滑的特性，主要区别在于最高峰值点位置的不同，最高峰值点位置是对于互相关两波形时延信息的反映。
[0051]
而后我们需要构建训练用的卷积神经网络，如图4所示，该网络由两个通道的卷积神经网络与后端的全连接网络组成，其中两个卷积网络通道的结构相似，第一个卷积通道的输入为符合波形对，其中第一个卷积层c1为二维卷积层，其目的是合并两波形为一维序列数据，而后的4个卷积层c2～c5均为一维卷积层；第二个卷积通道的输入为波形互相关序列，此通道有c6～c11共6个卷积层。所有卷积层的步长stride设置为2，卷积核大小根据输入数据长度的奇偶性设置为2或3，在每一层卷积层后会连接一层leakyrelu层作为激活函数，其参数设置为0.1；最后将两个通道的输出连接到一起，共同输入至三个连续的全连接层神经网络中。
[0052]
神经网络的输出为一个元素值，代表神经网络预测的符合事件位置距离探测器中心的距离，神经网络的损失函数为均方误差损失(mse)函数。由于距离与tof时间为线性关系，线性系数为光速，故而可以轻松地将预测距离转化为对tof的预测。
[0053]
为了对比本发明方法的有效性，我们设计了其它三种实验来求取tof值：
[0054]
(1)使用恒比定时鉴别器对波形进行计算，得到两波形对应的百分比时间点后相减得到tof信息。恒定比例设置为75％。
[0055]
(2)提取互相关序列的峰值点与零位置的时间差作为tof值。
[0056]
(3)使用单通道神经网络，只对符合波形对进行学习，得到tof预测值。
[0057]
表1为实验结果，为进行对比，结果为其它方法相较于恒比定时鉴别器时间分辨率的提升比率。
[0058]
表1
[0059][0060][0061]
表1中的互相关峰值点方法即为仅仅使用波形互相关方式所得到的结果，比恒比定时的方式提升20.35％，可见波形互相关的方式能够提取更为准确的tof时间信息；单通道网络方法即为使用深度学习的方式对两波形信号进行学习，以得到更好的结果，相较于恒比定时方式，提升了68.03％；互相关双通道网络即为本发明所提出的方法，相较于恒比定时方式，提升了72.82％，是几种方法中提升最大的方案。
[0062]
由以上实验结果可以看出，本发明基于波形信号互相关的pet飞行时间估计方法有效地提升了tof时间估计的分辨率。
[0063]
上述对实施例的描述是为便于本技术领域的普通技术人员能理解和应用本发明，熟悉本领域技术的人员显然可以容易地对上述实施例做出各种修改，并把在此说明的一般
原理应用到其他实施例中而不必经过创造性的劳动。因此，本发明不限于上述实施例，本领域技术人员根据本发明的揭示，对于本发明做出的改进和修改都应该在本发明的保护范围之内。