CN112927132B - PET image reconstruction method for improving spatial resolution uniformity of PET system - Google Patents

Abstract

The invention discloses a PET image reconstruction algorithm for improving the spatial resolution uniformity of a PET system based on deep learning, which fully utilizes the characteristic of nonuniform radial resolution in the FOV space of the PET system to solve the mutual depth effect, takes a high-resolution concentration distribution map reconstructed by projection data when phantom is positioned in the center of the FOV as a label, and improves the resolution of the reconstructed concentration distribution map of the same phantom positioned at the edge position of the FOV by means of neural network training without any novel detector or obtaining any additional information, such as DOI information or PSF information. The invention uses software means to replace the complex hardware method at the present stage, and solves the problem of nonuniform spatial resolution of the PET system.

Description

A PET Image Reconstruction Method for Improving Spatial Resolution Uniformity of PET System

technical field

The invention belongs to the technical field of biomedical image analysis, and in particular relates to a PET image reconstruction algorithm based on deep learning to improve the spatial resolution uniformity of a PET system.

Background technique

Positron emission tomography (PET) is a non-invasive functional imaging technique and one of the important imaging methods in nuclear medicine and molecular imaging. In the early stage of the disease, biochemical changes often precede anatomical changes. Because PET can use ¹⁵ O, ¹⁸ F and other radionuclide-labeled substances such as glucose and protein as tracers to participate in the normal physiological metabolism of organisms, it can be used at the molecular level. It dynamically and quantitatively reflects the pathophysiological changes and metabolic processes in animals, so it plays an irreplaceable role in the early diagnosis and treatment of heart diseases, brain diseases and malignant tumors. When these radionuclides decay in the body, the positrons generated by them will collide with the surrounding tissue and lose their kinetic energy almost immediately and combine with the negative electrons to undergo an annihilation reaction, and then radiate a pair of opposite directions with the same energy and 511kev gamma photons. According to the consistent principle, the in vitro PET detector can detect the decay site of γ photons in the body, obtain the original projection data (sinogram), and use different algorithms to process the sinogram to reconstruct the concentration of radionuclides in the organism. Distribution.

Undoubtedly, the resolution of the PET reconstructed concentration distribution map has a great influence on the accuracy of the doctor's diagnosis of the disease, so the uniform high-resolution image reconstruction occupies a particularly important position in the PET system. However, since the depth at which the annihilation photon interacts with the detector crystal is arbitrary, and even penetration occurs, there is a severe radial depth of interaction (DOI) effect in the defined lateral field of view (FOV) of the PET system Or radial parallax (parallax error), that is, the resolution of the reconstructed image at different radial positions is not uniform, especially when the object is close to the edge of the field of view, the resolution will be significantly reduced and the severe "smearing" effect will appear .

Different PET detector structures, such as ring PET and flat PET, all have radial DOI effects. At present, the hardware method is considered to be the main method to solve the radial DOI problem. It mainly designs new detectors based on different principles, so as to obtain accurate DOI information in the crystal. To affect the extent of the DOI problem, some use a phosphor sandwich or a discrete multilayer detector composed of a multi-layer scintillation crystal array with staggered peaks, and some couple two photoelectric conversion devices to the two ends of the crystal array. Take the signal output to judge the position where the Gamma photon acts on the crystal, and obtain continuous DOI information. Obviously, the above method requires a more complex detector structure and signal processing process, and the use of software means to replace the hardware DOI detector can achieve the effect of a certain degree of precision DOI detector, thereby solving the problem of uneven spatial resolution of the PET system. It is a worthy research direction to reduce the requirements for hardware equipment.

With the emergence and popularization of machine learning, researchers began to study machine learning to improve the uniformity of PET spatial resolution. For example, in the article Convolutional Neural Network for Crystal Identification and Gamma Ray Localization in PET, the hardware foundation of DOI encoder was proposed The deep learning method is introduced to more accurately locate the depth of interaction between Gamma photons and crystals. This technology reduces the requirements for the structural parameters of PET to a certain extent, but still requires a complex PET structure.

SUMMARY OF THE INVENTION

In view of the above, the present invention provides a PET image reconstruction algorithm based on deep learning to improve the spatial resolution uniformity of a PET system. Projection and Neural Network Solving.

A PET image reconstruction algorithm based on deep learning to improve the spatial resolution uniformity of a PET system, comprising the following steps:

(1) Inject PET radioactive tracer into the body membrane, place the body membrane at different positions radially from the center of the field of view (FOV) of the PET system for scanning, detect and count the coincident photons, and obtain the phase at different radial positions. corresponding projection data;

(2) According to the PET measurement equation, the reconstruction problem of projection data corresponding to different radial positions is divided into two sub-problems: sub-problem 1 is the reconstruction of the projection data corresponding to the position i=0, and sub-problem 2 is the position where i>0 Corresponding projection data reconstruction, i represents the radial distance from the center of the field of view;

(3) For sub-problem 1, the filtered back-projection algorithm (FBP) is used for reconstruction, and for sub-problem 2, the deep learning method is used for reconstruction;

(4) Build the ISTA-Net network model, use the projection data scanned at different radial positions with i>0 as the input of the model, use the PET image reconstructed from sub-problem 1 as the ground truth label of the corresponding output of the model, and obtain through iterative training The PET image reconstruction model is used to realize the reconstruction process of sub-problem 2.

Further, the scanning manner in the step (1) may be static scanning or dynamic scanning.

Further, the filtering back-projection algorithm in the step (3) includes two steps of frequency-domain filtering and back-projection.

算子经过软阈值算法后与一个

算子连接组成，所述

算子从输入至输出依次由卷积层A1、卷积层A2、Relu函数、卷积层A3依次连接组成；所述

算子与

算子的结构呈镜像对称，其从输入至输出依次由卷积层A3、Relu 函数、卷积层A2、卷积层A1依次连接组成，卷积层A1与A2之间插入有批归一化(BN)层，

算子的输出以及

算子的输入均经批归一化处理，

算子的输入与

算子的输出线性叠加后以残差的形式作为phase的最终输出结果。Further, the ISTA-Net network model is composed of multiple phases (phases) connected in sequence, and each phase consists of a

After the operator goes through the soft threshold algorithm, it is combined with a

The operator concatenates the composition, the

The operator is composed of convolutional layer A1, convolutional layer A2, Relu function, and convolutional layer A3 connected in sequence from input to output; the

operator and

The structure of the operator is mirror-symmetrical. It consists of convolutional layer A3, Relu function, convolutional layer A2, and convolutional layer A1 connected in sequence from input to output. Batch normalization is inserted between convolutional layers A1 and A2. (BN) layer,

the output of the operator and

The input of the operator is processed by batch normalization,

The input to the operator is the same as

The output of the operator is linearly superimposed in the form of residual as the final output result of the phase.

Further, the outputs of the convolutional layers A1 to A3 are batch standardized, and each convolutional layer is provided with 32 filters, the size of the convolution kernel used is 3×3, and the step size is 1.

Further, the specific process of training the ISTA-Net network model in the step (4) is: firstly initialize the model parameters, and divide the projection data obtained by scanning all i>0 different radial positions as samples into a training set and a test. The training set samples are input into the network model one by one, and the output results of the network model are obtained through forward propagation calculation; then the loss function L between each output result of the model and the corresponding true value label is calculated. According to the partial derivative of the loss function L The parameters in the network model are iteratively optimized by the Adam algorithm until the loss function L converges, and finally the PET image reconstruction model is obtained after the training is completed.

Further, the expression of the loss function L is as follows:

为第i个样本输入网络模型后其中第n个phase的输出结果，A为x_i的像素尺寸大小，n为ISTA-Net网络模型中的 phase数量，

为网络模型中第k个phase中的

算子，

为网络模型中第k个phase中的

算子，‖‖₂表示L2范数，

为第i个样本输入网络模型后其中第k个phase的输出结果，B为训练集中的样本数量。Where: x _i is the true value label corresponding to the ith sample,

is the output result of the n-th phase after the i-th sample is input to the network model, A is the pixel size of x _i , n is the number of phases in the ISTA-Net network model,

is the kth phase in the network model

operator,

is the kth phase in the network model

operator, ‖‖ ₂ represents the L2 norm,

is the output result of the kth phase after the ith sample is input to the network model, and B is the number of samples in the training set.

Further, the input of the ISTA-Net network model is 2D PET scan data. If the collected projection data is in a 3D form, then methods such as SSRB or FORB need to be used to convert it into 2D PET scan data. Different methods Different effects on the radial resolution of the same cross section.

Further, before training the ISTA-Net network model, the input projection data and its corresponding ground truth label need to be normalized for a single frame. The specific calculation formula is as follows:

Among them: X is any data value before normalization in the projection data or its corresponding truth label, X _norm is the data value after X corresponding normalization, and X _max is the projection data or its corresponding truth label. The maximum value of , X _min is the minimum value in the projected data or its corresponding ground-truth label.

The invention makes full use of the feature of the non-uniform radial resolution in the FOV space of the PET system to solve the mutual depth effect, and proposes a software method based on deep learning to improve the uniformity of the resolution of the PET system. The phantom is located in the center of the FOV. The high-resolution concentration distribution map reconstructed from the projection data at the time is used as the label, and the neural network training method is used to improve the resolution of the reconstructed concentration distribution map of the same phantom located at the edge of the FOV, without any new detectors or additional acquisition. Any information, such as DOI information or PSF information. The present invention replaces the complicated hardware method at the present stage by means of software, and solves the problem of uneven spatial resolution of the PET system.

Description of drawings

FIG. 1 is a schematic diagram of the implementation flow of the PET image reconstruction algorithm of the present invention.

FIG. 2 is a schematic structural diagram of the ISTA-Net network model of the present invention.

Figure 3 is a comparison chart of the PET reconstruction results of derenzo phantom at different radial positions (0.70mm, 0.95mm, 1.40mm, 1.95mm, 2.40mm, 2.80mm). The reconstruction map under the 4-layer DOI detector, the reconstruction map of the phantom located in the center of the field of view without DOI information, the reconstruction map of the phantom with no DOI information at a distance of 75mm from the center of the field of view, and the reconstruction map obtained by using the algorithm of the present invention.

Detailed ways

In order to describe the present invention more specifically, the technical solutions of the present invention will be described in detail below with reference to the accompanying drawings and specific embodiments.

As shown in FIG. 1 , the image reconstruction algorithm of the present invention for improving the spatial resolution uniformity of the PET system by deep learning specifically includes the following steps:

(1) Collect data. The PET radiotracer was injected into the donor film (phantom), and the phantom was placed at different positions radially from the center of the field of view (FOV) of the PET equipment for scanning, and the coincident photons were detected and counted, and the results obtained from different radial positions were obtained. The corresponding original projection data matrix Yi at _i .

(2) According to the principle of PET imaging, establish the measurement equation model:

Y=GX+R+S

Among them: G is the system matrix, X is the real tracer concentration distribution map, R is the number of random photons in the measurement process, and S is the number of scattered photons in the measurement process.

Due to the depth of interaction in the PET system, its spatial resolution is not uniform, that is, images located in the center of the same cross-sectional field of view have high resolution. The reconstruction problem of projection data Y _i at different radial positions is divided into two sub-problems, sub-problem 1 is to deal with the original projection data Y ₀ at the center of the field of view, and sub-problem 2 is to deal with different radial positions i (i≠0) The corresponding raw projection data matrix Y _i .

Sub-problem 1 directly uses filtered back projection (FBP) for image reconstruction, and sub-problem 2 is based on the reconstructed image of sub-problem 1 and is solved by means of deep learning (ISTA-Net).

经过软阈值算法后和其对称的算子

组成，其中算子

包括2个卷积层、一个Relu函数和另一个卷积层依次连接组成，网络中的每一个卷积神经网络层的输出均依次经过批量标准化处理；

输出的结果经过软阈值算法，而

则是与

完全对称的结构，先对数据进行批归一化处理再经过一个卷积层、一个Relu函数和另外两个卷积层，且两个卷积层之间插有一个批归一化(BN)层，并将

最后输出数据同输入的投影数据进行线性叠加以残差的形式输出。本实施方式的网络中一共设置了9个phase，每个phase里的所有的卷积层均设置有32个滤波器，每个卷积核的尺寸为3*3，步长为1。The structure of the network model adopted by the present invention is shown in Figure 2, which consists of multiple identical phases connected in sequence, and each phase consists of an operator from input to output

After the soft threshold algorithm and its symmetric operator

composition, where the operator

It consists of 2 convolutional layers, a Relu function and another convolutional layer connected in sequence, and the output of each convolutional neural network layer in the network is batch standardized in turn;

The output result goes through a soft threshold algorithm, while

is with

Completely symmetrical structure, batch normalization is performed on the data first, and then a convolutional layer, a Relu function and two other convolutional layers are passed, and a batch normalization (BN) is inserted between the two convolutional layers. layer, and will

Finally, the output data is linearly superimposed with the input projection data and output in the form of residual. A total of 9 phases are set in the network of this embodiment, all convolution layers in each phase are set with 32 filters, the size of each convolution kernel is 3*3, and the step size is 1.

(3) Training stage.

First, the input data (sinogram) and the corresponding label need to be normalized for a single frame:

Where: X _min and X _max are the minimum and maximum values of a single frame of data, respectively.

Then initialize the parameters of ISTA-Net, divide the original projection data matrix Y _i of different radial positions into training set and test set, and input the sinogram of the training set into ISTA-Net, and calculate the output of each layer through the forward propagation formula, Then obtain the final output of ISTA-Net, and calculate the loss function between the output of ISTA-Net and the label:

where: x _i is the label of the ith sample, n is the total number of stages of ISTA-Net, B is the number of total training sets processed in ISTA-Net, and A is the pixel dimension of x _i .

Finally, the partial derivative of the loss function is obtained, the learnable parameters in ISTA-Net are updated through the Adam algorithm, and forward propagation and reverse derivation are repeated until the value of the loss function is small enough.

(4) Estimation stage.

In the estimation stage, the projection data closest to the FOV edge in the radial direction is input into the trained ISTA-Net, and a relatively high-quality reconstruction map is directly obtained, thereby greatly weakening the influence of the depth effect and improving the problem of uneven spatial resolution.

Below we conduct experiments based on Monte Carlo simulation data to verify the effectiveness of this implementation. The tracer for Monte Carlo simulation is ¹⁸ F-FDG, the phantom is derenzo phantom (0.70mm, 0.95mm, 1.40mm, 1.95mm, 2.40mm, 2.80mm), and the simulated scanner is a 4-layer DOI detector, The DOI accuracy is 5mm and the single-layer detector without DOI information, in which the data obtained from the single-layer detector without DOI information at different radial positions from the center of the field of view are randomly divided into training sets (1080 sinograms) and Test set (360 sinograms). The training set is used to learn the parameters of ISTA-Net, and the test set is used to check the performance of the trained ISTA-Net.

Figure 3 compares the reconstruction results of the present invention and the 4-layer DOI detector when the derenzo phantom is at the edge of the FOV. The columns from left to right are: the reconstructed image of the phantom at a distance of 75mm from the center of the field of view under the 4-layer DOI detector, The reconstruction map of the phantom without DOI information in the center of the visual field, the reconstruction map of the phantom with no DOI information at a distance of 75 mm from the center of the visual field, and the reconstruction map obtained by using the algorithm of the present invention. As can be seen from Figure 3, even if the phantom at the center has no DOI information, the resolution of its concentration distribution map is also very high, while the phantom at the edge position uses FBP (when there is no DOI information) to reconstruct the resolution of the concentration distribution map. Very low, even a point source with a diameter of 2.8mm is indistinguishable, but the image resolution after obtaining DOI information can be slightly higher than the reconstructed image at the center without DOI information. We use the present invention to reconstruct the sinogram at the edge position and find that its resolution has been significantly improved, and the effect can be close to the acquisition of 4-layer DOI information, so that the resolution at the edge and center of the FOV is closer to the same, and the PET system space is improved. Uniformity of resolution.

The above description of the embodiments is for the convenience of those of ordinary skill in the art to understand and apply the present invention. It will be apparent to those skilled in the art that various modifications to the above-described embodiments can be readily made, and the general principles described herein can be applied to other embodiments without inventive effort. Therefore, the present invention is not limited to the above-mentioned embodiments, and improvements and modifications made to the present invention by those skilled in the art according to the disclosure of the present invention should all fall within the protection scope of the present invention.

Claims (7)

1. A PET image reconstruction method for improving the spatial resolution uniformity of a PET system based on deep learning comprises the following steps:

(1) injecting a PET radioactive tracer into the body membrane, respectively placing the body membrane at different positions of the PET system in the radial direction from the center of the visual field for scanning, detecting coincident photons and counting to obtain corresponding projection data at different radial positions;

(2) the reconstruction problem of the projection data corresponding to different radial positions is split into two sub-problems according to a PET measurement equation: sub-problem 1 is projection data reconstruction corresponding to a position where i is 0, sub-problem 2 is projection data reconstruction corresponding to a position where i > 0, and i represents a radial distance from the center of the field of view;

(3) reconstructing the subproblem 1 by adopting a filtering back projection algorithm, and reconstructing the subproblem 2 by adopting a deep learning method;

(4) an ISTA-Net network model is built, projection data obtained by scanning different radial positions with i being larger than 0 are used as input of the model, a PET image obtained by reconstruction of a subproblem 1 is used as a truth label of corresponding output of the model, a PET image reconstruction model is obtained through iterative training to realize the reconstruction process of a subproblem 2, and the specific training process is as follows: firstly, initializing model parameters, taking projection data obtained by scanning different radial positions with i larger than 0 as samples to be divided into a training set and a testing set, inputting the training set samples into a network model one by one, and obtaining an output result of the network model through forward propagation calculation; calculating a loss function L between each output result of the model and the corresponding truth label, continuously performing iterative optimization on parameters in the network model through an Adam algorithm according to the partial derivative of the loss function L until the loss function L is converged, and finally obtaining a PET image reconstruction model after training is completed;

wherein: x is the number of_iThe true label for the ith sample,

the output result of the nth phase is input into the network model for the ith sample, and A is x_iN is the number of phases in the ISTA-Net network model,

for in the kth phase in the network model

The operator(s) is (are) selected,

for in the kth phase in the network model

Operator, | | | purple wind₂The norm of L2 is shown,

and B is the number of samples in the training set, wherein the output result of the kth phase is the ith sample after the ith sample is input into the network model.

2. The PET image reconstruction method according to claim 1, characterized in that: the scanning mode in the step (1) may be static scanning or dynamic scanning.

3. The PET image reconstruction method according to claim 1, characterized in that: the filtering back projection algorithm in the step (3) comprises two steps of frequency domain filtering and back projection.

4. The PET image reconstruction method according to claim 1, characterized in that: the ISTA-Net network model is formed by connecting a plurality of phases in sequence, and each phase is formed by one phase

The operator is subjected to a soft threshold algorithm and then is compared with one

Operator concatenation composition, said

The operator is formed by sequentially connecting a convolutional layer A1, a convolutional layer A2, a Relu function and a convolutional layer A3 from input to output; the described

Operator and

the structure of the operator is mirror symmetry, and the operator is formed by sequentially connecting a convolution layer A3, a Relu function, a convolution layer A2 and a convolution layer A1 from input to output, wherein a batch normalization layer is inserted between the convolution layers A1 and A2,

output of operator and

the inputs of the operators are all processed by batch normalization,

input of operator and

and the output of the operators is linearly superposed and then is used as the final output result of the phase in the form of residual error.

5. The PET image reconstruction method according to claim 4, characterized in that: the outputs of the convolutional layers A1-A3 are subjected to batch standardization, each convolutional layer is provided with 32 filters, the size of an adopted convolutional kernel is 3 multiplied by 3, and the step length is 1.

6. The PET image reconstruction method according to claim 1, characterized in that: the input of the ISTA-Net network model is 2D PET scanning data, if the acquired projection data is in a 3D form, the projection data needs to be converted into the 2D PET scanning data by adopting an SSRB or FORB method, and different methods have different influences on the radial resolution of the same cross section.

7. The PET image reconstruction method according to claim 1, characterized in that: before training the ISTA-Net network model, the input projection data and the corresponding truth label need to be normalized by a single frame, and the specific calculation formula is as follows:

wherein: x is any data value before normalization processing in projection data or corresponding truth value label thereof, and X is_normFor X to correspond to the normalized data value, X_maxFor the maximum, X, in the projection data or its corresponding truth label_minIs the minimum in the projection data or its corresponding truth label.

Images