X-ray diagnostic system, image processing apparatus, and program
By training a deep learning network with physical model information and analytical reconstruction images, the method enhances X-ray CT image quality and reduces computational time, overcoming the limitations of existing iterative reconstruction methods.
Patent Information
- Application Number
- JP2020098817
- Authority / Receiving Office
- JP · JP
- Patent Type
- Patents
- Current Assignee / Owner
- Priority Date
- 2019-07-12
- Filing Date
- 2020-06-05
- Publication Date
- 2025-07-30
- Estimated Expiration
- 2040-06-05
AI Technical Summary
Existing X-ray CT reconstruction methods, such as iterative reconstruction (IR) and model-based iterative reconstruction (MBIR), require significant computational resources and time, limiting their widespread adoption despite their potential for improved image quality.
A deep learning (DL) network is trained using physical model information and analytical reconstruction images to filter initial CT images, combining fast analytical reconstruction with DL to achieve image quality comparable to MBIR, by incorporating physical model data as input to enhance the DL network's ability to produce high-quality images.
The proposed method significantly reduces computational time while achieving image quality comparable to MBIR, addressing the trade-off between speed and quality in CT image reconstruction.
Smart Images

Figure 0007715492000013 
Figure 0007715492000014 
Figure 0007715492000015
Abstract
Description
Technical Field
[0001] The embodiments disclosed in this specification and the drawings relate to an X-ray diagnostic system, an image processing apparatus, and a program.
Background Art
[0002] As image reconstruction methods in an X-ray CT apparatus, there are an iterative reconstruction (IR) method and a model-based iterative reconstruction (MBIR) method. These reconstruction methods can improve image quality as compared with analytical reconstruction methods such as a filtered back projection (FBP) method.
[0003] However, it is known that the IR method and the MBIR method require a large amount of calculation time.
Prior Art Documents
Patent Documents
[0004]
Patent Document 1
Patent Document 2
Patent Document 3
Patent Document 4
Non-Patent Documents
[0005]
Non-Patent Document 1
[0006] One of the problems to be solved by the embodiments disclosed in this specification and the drawings is to improve the image quality. However, the problems to be solved by the embodiments disclosed in this specification and the drawings are not limited to the above problems. The problems corresponding to the effects of each configuration shown in the embodiments described later can also be regarded as other problems. [Means for Solving the Problems]
[0007] The X-ray diagnostic system according to the embodiment includes an acquisition unit, a reconstruction unit, and a generation unit. The acquisition unit acquires X-ray projection data. The reconstruction unit reconstructs a first image from the X-ray projection data based on a first CT method. The generation unit inputs the first image and the data of the physical model to a neural network that has been trained using a training data set and physical model information and generates the data of the physical model based on at least one of the projection data or the first image, thereby generating a second image. The training data set includes input data that is an image reconstructed using the first CT method and target data that is an image reconstructed using a second CT method, and the physical model information is generated in the same manner as the generation of the data of the physical model.
Brief Description of the Drawings
[0008]
Figure 1
Figure 2A
Figure 2B
Figure 3
Figure 4
Figure 5
Figure 6
Modes for Carrying Out the Invention
[0009] Hereinafter, embodiments of an X-ray diagnostic system, an image processing apparatus, and a program will be described in detail with reference to the drawings.
[0010] First, the overall configuration of the X-ray diagnostic system according to the embodiment will be described with reference to FIG. 6.
[0011] FIG. 6 shows a non-limiting example of a CT scanner. As shown in FIG. 6, the X-ray imaging gantry 500 is shown in a side view, and further includes an X-ray tube 501 as an X-ray source, an annular frame 502, and a multi-row or two-dimensional array type X-ray detector 503. The X-ray tube 501 and the X-ray detector 503 are mounted diametrically across the imaging object S on the annular frame 502 rotatably supported around the rotation axis RA.
[0012] The multi-slice X-ray CT apparatus further includes a high voltage generator 509 that generates the tube voltage applied to the X-ray tube 501 via the slip ring 508, and the X-ray tube 501 generates X-rays. The X-rays are radiated toward the imaging object S whose cross-sectional area is represented by a circle. For example, the X-ray tube 501 has an average X-ray energy during the first scan that is smaller than the average X-ray energy during the second scan. Therefore, two or more scans corresponding to different X-ray energies can be acquired. The X-ray detector 503 is arranged on the opposite side of the X-ray tube 501 across the imaging object S in order to detect the transmitted X-rays of the imaging object S. The X-ray detector 503 further includes individual detector elements or units. That is, the X-ray detector 503 as a detector array is arranged diametrically opposite to the X-ray tube 501 as the X-ray source across the opening of the gantry, and detects the X-rays emitted from the X-ray tube 501 to generate projection data.
[0013] The CT apparatus further includes another apparatus for processing the detection signal from the X-ray detector 503. The data acquisition circuit or data acquisition system (Data Acquisition System: DAS) 504 converts the signal output from the X-ray detector 503 of each channel into a voltage signal, amplifies this signal, and further converts this signal into a digital signal.
[0014] The above data is sent via the contactless data transmitter 505 to the preprocessing circuit 506 housed in the external console of the X-ray gantry 500. The preprocessing circuit 506 performs several corrections such as sensitivity correction of the raw data. The storage device 512 stores the resulting data, which is also called projection data at the stage immediately before the reconstruction process. The storage device 512 is connected to the system controller 510 via the data / control bus 511 together with the acquisition circuit 518, the generation circuit 517, the reconstruction circuit 514, the input interface 515, and the display device 516. The system controller 510 controls the current regulator 513, which limits the current to a sufficient level to drive the CT system.
[0015] The storage device 512 can store measurement values representing the X-ray irradiance at the X-ray detector 503. Further, the storage device 512 can store a dedicated program for executing Method 10.
[0016] The reconstruction circuit 514 can execute various steps of Method 10. Further, the reconstruction circuit 514 can execute image processing before the reconstruction process, such as volume rendering processing and image difference processing, as necessary.
[0017] The acquisition circuit 518 acquires projection data from the preprocessing device 506.
[0018] The generation circuit 517 performs learning using, for example, a neural network, and executes the learned model to generate an image. The detailed processing of the generation circuit 517 will be described in detail.
[0019] The pre-reconstruction preprocessing of the projection data executed by the preprocessing circuit 506 can include, for example, corrections for detector calibration, detector non-linearity, and polar effect.
[0020] The post-reconstruction processing executed by the reconstruction circuit 514 can include, as necessary, filtering and smoothing of images, volume rendering processing, and image difference processing. The image reconstruction process can implement various steps of Method 10. The reconstruction circuit 514 can use memory to store, for example, projection data, reconstructed images, calibration data and parameters, and computer programs.
[0021] Various circuits (e.g., the reconstruction circuit 514, the preprocessing circuit 506, the acquisition circuit 518, and the generation circuit 517) can include a CPU (processing circuit) that can be implemented as individual logic gates, an application specific integrated circuit (ASIC), a field programmable gate array (FPGA), or other complex programmable logic device (CPLD). The implementation of the FPGA or CPLD can be coded in VHDL, Verilog, or other hardware description languages, and the code can be stored directly in the electronic memory within the FPGA or CPLD or stored as separate electronic memory. Further, the storage device 512 can be non-volatile memory such as ROM, EPROM, EEPROM, or flash memory. The storage device 512 can also be volatile, such as static or dynamic RAM, and a processor such as a microcontroller or microprocessor can be provided to manage the interaction between the electronic memory and the FPGA or CPLD and the memory.
[0022] Note that the reconstruction circuit 514, the acquisition circuit 518, and the generation circuit 517 are each an example of a reconstruction unit, an acquisition unit, and a generation unit, respectively.
[0023] In one embodiment, the reconstructed image can be displayed on a display device 516. The display device 516 can be an LCD display, a CRT display, a plasma display, an OLED, an LED, or any other display known in the art.
[0024] Subsequently, the background of the embodiment will be described.
[0025] Compared with the CT analysis reconstruction method, the iterative reconstruction (IR) method can improve the image quality by investigating the statistical characteristics of the measurement. The model-based iterative reconstruction (MBIR) method can generate better image quality than the IR method by incorporating the physical model of the CT system and the scanned object to be measured. However, since the MBIR method has a slow reconstruction speed, the computational amount increases with an accurate physical model, and the iterative image reconstruction becomes slow, the MBIR method is not widely adopted.
[0026] In the trade-off between computational complexity / speed and image quality, one approach to reach the optimal compromise point is to first reconstruct the initial CT image using the analytical reconstruction method, and then use a deep learning (DL) network (also called an Artificial Neural Network (ANN) or a Convolutional Neural Network (CNN)) to filter the initial CT image to generate a final CT image. Since the DL network can be trained using training data where the input data is the analytical reconstruction image and the target image is the IR reconstruction image, in response to the input of the analytical reconstruction image, the DL network learns to output an image approximating the image quality of the IR reconstruction image. Since the DL network is fast compared to the IR method and the MBIR method, the combination of the analytical reconstruction method and the DL network is fast and can provide image quality comparable to that of the MBIR method at least theoretically. However, in practice, it is difficult to achieve image quality comparable to that of the IR CT method when only the analytical reconstruction image is provided as the input. However, the situation changes when physical model information / data is also used as the input to the DL network.
[0027] Therefore, for the input to the DL network, the method described herein uses information from one or more physical models in addition to the analytical reconstruction image to improve the ability of the DL network to achieve high image quality.
[0028] Referring now to the drawings, in which like reference numerals refer to like or corresponding parts throughout several views, FIG. 1 shows a flowchart of a non-limiting example of a method 10 for performing image processing on a CT reconstructed image (referred to herein by abbreviation as a CT image) using a DL network 170 for training and use. The method 10 shown in FIG. 1 learns a method of optimally using physical model information 117 using the DL network 170 to filter an input CT image 112 reconstructed from raw data 105 using a fast analytical reconstruction algorithm. The method 10 includes two parts, namely, (i) an offline training process 150 and (ii) a medical imaging process 100. That is, the offline training process 150 trains the DL network 170, and the medical imaging process 100 filters the input CT image 112 using the trained DL network 170, thereby generating a final image 135 having improved image quality compared to the input CT image 112. In some embodiments, step 130 can be omitted.
[0029] The network 170 is trained using an offline DL training process 160. In the offline DL training process 160, the loss function is minimized by repeatedly adjusting the parameters of the DL network 170 (for example, the parameters of the DL network 170 can include weight coefficients connecting network layers and activation functions / potentials of nodes within the layer). Optimization of the network parameters continues until a stopping criterion (for example, the stopping criterion can be whether the value of the loss function has converged to a predefined threshold) is met, resulting in a trained DL network 170.
[0030] The loss function compares the target data 153 with the output resulting from applying the target data 153 as an input to the current version of the DL network 170 with the input data 157 and the physical model information 155. For example, the input data can be a CT image reconstructed using an analytical reconstruction algorithm (e.g., the same analytical reconstruction algorithm used in step 110 of the medical imaging process 100 is preferred).
[0031] The target data can be respective images generated using an IR method (or preferably an MBIR method) that operates on the same raw data used to generate the CT image of the input data.
[0032] The CT image of the input data can be referred to as a low-quality CT image, and the CT image of the target data can be referred to as a high-quality CT image. For example, the low-quality CT image may have low image quality because it has more noise and artifacts than the high-quality CT image, or because it has a lower resolution than the high-quality CT image. The relative image quality can be based on image quality criteria such as (i) the peak signal-to-noise ratio, (ii) the Structural Similarity (SSIM) index, the Visual Information Fidelity (VIF), the Blind / Referenceless Image Spatial Quality Evaluator (BRISQUE), the Natural Image Quality Evaluator (NIQE), and the Perception Based Image quality Evaluator (PIQE).
[0033] The physical model information 155 is of the same type as the physical model information 117 used in the medical imaging process 100. In the case of the X-ray CT image reconstruction method, the physical model can include X-ray scattering estimated from an initial patient model as a scattering object. Additional details of the physical model information 155 are described below.
[0034] In the training data, each low-quality CT image of the input data forms a pair with the corresponding high-quality CT image of the target data. Imaging scans for obtaining the low-quality CT images for the input data 157 and the high-quality CT images for the target data 153 can be performed, for example, on a phantom.
[0035] By applying the low-quality CT images from the input data to the current version of the DL network 170, an output from the network is generated that is likely to match as closely as possible the corresponding high-quality CT image from the target data. If the loss function indicates that the output from the DL network 170 does not match the target data 153, the weighting coefficients of the DL network 170 are repeatedly adjusted until the output from the DL network 170 matches the target data 153.
[0036] That is, the DL network 170 is trained by repeatedly adjusting the network coefficients within the DL network 170 to minimize the difference between the filtered CT image output from the DL network 170 and the high-quality CT image from the target data 153. The training of the DL network 170 is determined to be complete when the difference between the target data and the output of the DL network 170 is minimized. The question of whether this difference has been sufficiently minimized is resolved based on one or more predefined stopping criteria of the offline DL training process 160. Once the stopping criteria are met, the trained DL network 170 can be saved and then called for use in the medical imaging process 100.
[0037] In another embodiment, the DL network 170 is implemented as a Residual Network (ResNet). In this case, the method described herein can filter an image by treating the difference between a low-quality and a high-quality CT image as an additive residual that can be directly removed from the low-quality CT image. Thus, when a low-quality CT image is applied to the neural network, the network outputs an image corresponding to the difference image. Next, a corrected CT image can be generated by subtracting the output of the network (i.e., the residual) from the low-quality CT image, thereby generating a corrected / final CT image. That is, the DL network 170 is a residual network, and the generation circuit 517 generates a second image, which is the final image 135, by subtracting the output result of the DL network 170, which is a neural network, from the input image 112, which is the first image.
[0038] The medical imaging process 100 is performed by acquiring raw data 105, for example, by performing an X-ray CT scan to generate a sinogram (e.g., X-ray CT projections at a series of view angles). That is, the acquisition circuit 518 acquires X-ray projection data as the raw data 105.
[0039] In step 110 of the medical imaging process 100, a CT image is reconstructed from the noise-removed CT image. That is, the generation circuit 517 reconstructs a first image from the X-ray projection data based on a first CT method, such as an analytical reconstruction method.
[0040] Preferably, a fast analytical reconstruction method (e.g., Filtered Back Projection (FBP)) is used to reconstruct the input CT image 112. However, various methods can be used to reconstruct a CT image from projection data, including Filtered Back Projection (FBP) and Statistical Iterative Reconstruction (IR) algorithms.
[0041] Examples of analysis methods that can be used to reconstruct the input CT image 112 include: (i) the Feldkamp Davis Kress (FDK) method, (ii) the generalized FDK method, (iii) the rebinning FBP method, (iv) the n-Pi method, (v) the Pi-slant method, the exact method of Katsevich, (vi) the Adaptive Multiple Plane Reconstruction (AMPR) method, (vii) the Advanced Single-Slice Rebinning (ASSR), (viii) the weighted FBP method, and (ix) the Adaptive Iterative Dose Reduction 3D (AIDR 3D) method.
[0042] Compared with the FBP reconstruction method, the image quality is improved in the IR method. An example of the IR method is executed by performing an optimization search to find an argument f that minimizes an objective function C(f) (also called a cost function), as shown in Equations (1) and (2).
[0043]
Number
[0044]
Number
[0045] Here, f * is the optimal image to be reconstructed, l is the projection data representing the logarithm of the X-ray intensity of the projection images taken at a series of projection angles, and f is the reconstructed image of the X-ray attenuation for the voxels / volume pixels in the image space (or the 2D pixels of the 2D reconstructed image). In the system matrix A, each matrix value a ij (i is the row index, j is the column index) is the volume corresponding to the voxel f j and the projection value l irepresents the overlap with the X-ray trajectory corresponding to f. When the forward projection A of the reconstructed image f provides a good approximation of all the measured projection images l, the data fidelity term ||Af-l || 2 W is minimized. The data fidelity term therefore involves solving the system matrix equation A=l, which represents the Radon transform (i.e., projection) of various rays from the source through the imaged object S in the space represented by f to the x-ray detector, which yields a value of l (e.g., x-ray projection through the 3D imaged object S to a 2D projection image l).
[0046] Notation ||g|| 2 w is g T means a weighted dot product of the form Wg, where W is a weighting matrix (e.g., representing the reliability of the projection data based on the signal-to-noise ratio per pixel). In other embodiments, the weighting matrix W can be replaced by an identity matrix. When the weighting matrix W is used in terms of data fidelity, the above IR method is called the Penalized Weighted Least Squares (PLWS) approach.
[0047] The function U(p) is a regularization term, which aims to impose one or more constraints (e.g., a Total Variation (TV) minimization constraint), which often has the effect of smoothing or denoising the reconstructed image. The value β is a regularization parameter, which weights the relative contributions of the data fidelity and regularization terms.
[0048] Additional image domain denoising is performed in step 130 of the medical imaging process 100. This step is optional and may be omitted in some embodiments.
[0049] Examples of noise removal methods include linear smoothing filters, anisotropic diffusion, non-local means, or non-linear filters. A linear smoothing filter removes noise by convolving the original image with a convolution kernel representing a low-pass filter or a smoothing operation.
[0050] For example, a Gaussian convolution kernel is composed of elements determined by a Gaussian function. Through this convolution, the value of each pixel more closely matches the values of adjacent pixels. Anisotropic diffusion removes noise while maintaining sharp edges by evolving the image under a smoothing partial differential equation similar to the heat equation. The median filter is an example of a non-linear filter, and if properly designed, a non-linear filter can also preserve edges and prevent blurring. The median filter is an example of a Rank-Conditioned Rank-Selection (RCRS) filter, and when applied, it can remove salt-and-pepper noise without producing significant blurring artifacts in the image.
[0051] In addition, if the assumption of uniformity across large regions where the imaged area is separated by sharp boundaries between uniform regions is supported, a filter using a Total Variation (TV) minimization regularization term can be applied. The TV filter is another example of a non-linear filter. Furthermore, non-local means filtering is an exemplary method for determining the pixel with noise removed by using the weighted average of similar patches within the image.
[0052] Finally, a final image (reconstructed image) 135 with good image quality is output, and this final image (reconstructed image) 135 is either displayed to the user or saved for later use.
[0053] Figures 2A and 2B show how the DL network 170 is used in step 120 and the offline DL training process 160, respectively. In step 120, the CT image and the physical model information 117 are inputs to the DL network 170, and the filtered image 122 is the output from the DL network 170. If step 130 is omitted, the filtered image 122 is the final image (output image) 135. That is, the generation circuit 517 inputs the input image 112, which is the first image, and the physical model information 117, which is the data of the physical model generated in step 115, to the DL network 170, which is a neural network that has been trained, to generate the second image, which is the final image.
[0054] In the offline DL training process 160, the inputs to the DL network 170 are the input data 157 (i.e., the low-quality image) and the physical model information 155, and the output from the DL network 170 is the filtered image 162. The loss function combines the filtered image 162 and the high-quality image of the target data to generate a value representing how closely these two images match. To train the DL network 170, the parameters of the DL network 170 are adjusted to optimize (e.g., minimize) the value of the loss function.
[0055] Regarding the physical model information 155 and the physical model information 117, many different physical models can be used to provide the physical model information. Also, various methods can be used to extract information from the physical model. Generally, the physical model information is organized / formatted in the same way when used in the offline training process 150 to train the DL network 170 and when used in the medical imaging process 100. To show the generation of the physical model information 155 and the physical model information 117, non-limiting examples are shown. In this non-limiting example, the physical model information is provided as the gradient of the objective function of the MBIR method.
[0056] The MBIR method includes a physical model, and information of the physical model is provided by the gradient of the objective function. Next, the DL network 170 can be trained to improve the image quality based on the physical model using the physical model information presented / encoded in the form of the gradient of the objective function.
[0057] Thus, in the offline DL training process 160, the learning of the DL network 170, which is a neural network, is performed using the training dataset and the physical model information 155. The training dataset includes the input data 157, which is an image reconstructed using a first CT method such as the analytical reconstruction method, and the target data 153, which is an image reconstructed using a second CT method such as an IR method (Iterative Reconstruction) typified by the MBIR method. The physical model information 155 is generated in the same manner as the generation of the data of the physical model in step 115 described later.
[0058] As described above, IR reconstruction can be solved by minimizing the objective function C(f). In an example of the objective function C(f) shown above, the system matrix A is included in the data fidelity term. In some model-based methods, a physical model can be incorporated into the system matrix. Further, the projection data term l can be modified as Pl in order to incorporate the physical model. For example, the system matrix can be modified to include beam hardening correction, and the projection data l can be modified to correct for scatter. That is, the objective function C(f) can be modified to be Equation (3). s As described above, IR reconstruction can be solved by minimizing the objective function C(f). In an example of the objective function C(f) shown above, the system matrix A is included in the data fidelity term. In some model-based methods, a physical model can be incorporated into the system matrix. Further, the projection data term l can be modified as Pl in order to incorporate the physical model. For example, the system matrix can be modified to include beam hardening correction, and the projection data l can be modified to correct for scatter. That is, the objective function C(f) can be modified to be Equation (3).
[0059]
Equation
[0060] Here, A tot =A fp A bh is the total system matrix, A fp is the forward projection matrix, and A bhis beam hardening correction, P s is scatter correction. Therefore, the objective function C(f) can be modified to incorporate various physical models, resulting in more accurate image reconstruction at the cost of increased complexity in the calculation of each iteration of the IR method. When the IR method is modified to incorporate a physical model, the method is referred to as the MBIR method.
[0061] Given an objective function incorporating one or more physical models, information on the physical model can be generated by estimating the gradient of the objective function given by Equation (4).
[0062]
Equation
[0063] That is, the generation circuit 517 generates data of the physical model using the MBIR method, for example, using an estimated value of the gradient of the objective function of the MBIR method. In some embodiments, an adjoint function-based gradient estimation method can be used to estimate the gradient of the objective function (or the data fidelity term). In some embodiments, the gradient of the objective function can be estimated using one or more iterations of the MBIR method or using only the data fidelity term of the MBIR method. That is, the generation circuit 517 updates an image based on, for example, the data fidelity term of the objective function of the MBIR method and generates the updated image as data of the physical model.
[0064] For example, in many IR CT methods, f (n+1) = f (n) - λ∇C(f (n) ) is repeatedly calculated to obtain the final image f * . Here, λ and ∇C(f (n) ) are the step size and gradient of the objective function, respectively, and include the gradients of the data fidelity term and the regularization term. That is, the gradient of the objective function is obtained, for example, by performing iterations of the IR method and using the step direction f (n+1) - f (n)can be estimated using. Therefore, step 115 executes one (or more) iterations of the MBIR method, and the estimated gradient ∇C ∝ f (1) -f (0) by generating model information using. The initial estimated value of the reconstructed image f (0) can be (but does not have to be) the input CT image 112.
[0065] In this way, in step 115, the generation circuit 517 generates data of the physical model, for example, using the MBIR method, based on at least one of the projection data or the first image.
[0066] Note that the process of generating model information can be further generalized by noting that the objective function is not limited to the above functional form. More generally, model-based reconstruction can be performed by solving the optimization problem given by Equation (5).
[0067]
Equation
[0068] Here, J(·, ·) is the data conformity term (also called the data fidelity term), and G(·) and P(·) are general physical models. Here, J(·, ·) acts on both the pre- and post-log measurement data and covers different evaluation criteria for different estimations including least squares estimation, maximum likelihood estimation, and maximum a posteriori probability estimation. As described above, one way to generate physical model information for input is to utilize gradient information, i.e., ∂J / ∂f.
[0069] In a specific case where the objective function takes the functional form described above, for example, as given by Equation (6), a completely physical model reconstruction can be performed by minimizing the following regularized weighted least squares cost function.
[0070]
Equation
[0071] Here, the matrices A and P are related to the physical models applied to the image domain and the projection domain, respectively. As described above, the final image f * is such that f (n+1) = f (n) - λ∇C(f (n) ) can be obtained by iterative calculation. One way to generate the physical model information 117 is to utilize the updated image IM generated using only the data fidelity term given by Equation 7.
[0072]
Equation
[0073] Here, f (0) is the input CT image 112 and is calculated through analytical reconstruction.
[0074] In some embodiments, the image f (0) may be different from the input CT image 112 (e.g., f (0) may be reconstructed using a reconstruction kernel different from that used to reconstruct the input CT image 112). In some embodiments, the first iteration of the data fidelity term IM is used as the physical model information 117.
[0075] In other embodiments, the physical model information 117 can be an estimate of the gradient (e.g., IM - f (0) , f (1) - f (0) , f (n) - f (0) , f (n) - f (n-1) etc.). Generally, the same method used to generate the physical model information 117 is also used to generate the physical model information 155.
[0076] The physical model information 155 can be generated using only one (or several) iterations of the MBIR method, but the MBIR method can be run until convergence to generate the target data. Also, the physical model information 155 can be generated using only the data fidelity term, but the MBIR method can be run using both the data fidelity term and the regularization term to generate the target data.
[0077] As described above, various physical models can be incorporated into the objective function, and these physical models are the above-mentioned scatter correction P s and beam hardening correction A bh and are not limited thereto. Multiple factors, such as a CT system (gantry, X-ray source, detector, etc.) and physical mechanisms (scatter in the patient, polychromaticity of radiation, etc.), can affect the quality of CT images.
[0078] Therefore, model-based reconstruction can improve imaging performance by considering both advanced physical models and the integration of various physical mechanisms. These heterogeneous physical models can be integrated into the data fidelity term of the objective function in both the image domain (e.g., the physical model can be incorporated into the system matrix A) and the projection domain (e.g., the physical model can be incorporated into the matrix P). More generally, operations in the image domain can be represented by the function G(f), and operations in the projection domain can be represented by the function P(l).
[0079] In certain embodiments, the physical models can include a deterministic radiative transfer scatter model for total scatter correction P s as described in U.S. Patent Application Nos. 15 / 210,657, 15 / 436,043, 16 / 392,177, and 16 / 252,392, which are hereby incorporated by reference in their entirety. For example, using the initial reconstructed image f (0) the scatter contribution to each pixel value of the projection data l can be estimated, and since the scatter correction can remove the contribution due to scatter, P sl becomes the main signal, which is obtained by subtracting the scattered signal from the total signal. That is, as an example of the physical model, an X-ray scattering model can be cited.
[0080] In some embodiments, the physical model includes a multi-material beam hardening model for beam hardening correction (A bh ) as discussed in the document entitled "Three-dimensional two material based beam hardening correction for iterative reconstruction" (I. Hein, Z. Yu, S. Nakanishi, "Three-dimensional two material based beam hardening correction for iterative reconstruction", Proc. 4th Int. Conf. Image Formation X-Ray Comput. Tomogr., pp. 455-458, 2016), which is incorporated herein by reference in its entirety. For example, the initial reconstructed image f (0) can be segmented into material components (e.g., by using material discrimination or by mapping the Hounsfield Units in each voxel to respective material components, such as bone, water, combinations of bone and water), and then forward projected using material-dependent spectral attenuation, taking into account beam hardening correction. Additionally, when generating the projection data l using dual-energy or spectral CT, material discrimination can be used to separate the attenuation of each voxel into the material component(s), and then this is used for forward projection to take into account beam hardening correction. That is, as an example of the physical model, a beam hardening model can be cited.
[0081] In certain embodiments, the physical model includes blurring of the source and detector (A psfcan include simulated and / or experimental system response simulation algorithms for []. Further, the physical model can include simulated and / or experimental system responses regarding the non-linearity of the detector, as described in U.S. Patent Application No. 14 / 593,818 (Patent No. 9,977,140), which is incorporated herein by reference in its entirety. That is, detector correction A psf can account for the point spread function resulting from the spatial resolution limitations of the detector elements (e.g., diffractive optical photons, charge sharing in direct detectors), and in some embodiments, detector correction A psf can account for the non-linearity of the detector due to k escape, energy resolution limitations, pile-up, and Compton scattering in gamma-ray detectors. That is, as an example of the physical model, a spatial resolution detection model can be mentioned.
[0082] Also, as an example of the physical model, a forward projection model or a system geometry model can be mentioned. For example, in some embodiments, the physical model can include an accurate forward projection (A fp ) that takes into account the system geometry, such as the advanced footprint method described in the document entitled "Distance-driven projection and backprojection in three dimensions" (D. De Man and s. Basu, “Distance-driven projection and backprojection in three dimensions,” Phys. Med. Biol., vol. 49, issue 11, page 2463 (2004)), which is incorporated herein by reference in its entirety. For example, using an accurate forward projection model A fp can use distance-driven projection and backprojection that provide a highly sequential memory access pattern. The forward projection model A fp can be executed by mapping the pixel boundaries and detector boundaries to a common axis and applying a kernel operation to map the data from one set of boundaries to another.
[0083] Distance-driven projection can be better understood in the context of pixel-driven backprojection and ray-driven projection. Pixel-driven backprojection works by connecting the lines from the focal point through the centers of the pixels of interest to the detector. When the positions of the intersections of the detector are calculated, values are obtained from the detector (usually linearly) by interpolation, and the results are accumulated in the pixel. Ray-driven projection works by connecting the lines from the focal point through the centers of the detector cells of interest through the image. For all image rows (or columns), the positions of the intersections are calculated, values are usually obtained from the image rows by linear interpolation, and the results are accumulated in the detector cell.
[0084] Distance-driven projection / backprojection combines the advantages of the ray-driven and pixel-driven methods. First, all views (or source positions) define a bijection between positions on the detector and positions within an image row (or column) (i.e., all points within an image row are uniquely mapped to points on the detector and vice versa). This makes it possible to define the length of the overlap between each image pixel and each detector cell. To calculate this overlap length, it would be possible to map all pixel boundaries in the image row of interest to the detector, or all detector cell boundaries to the center line of the image row of interest.
[0085] In practice, both boundary sets are mapped onto a common line. This is achieved by connecting all pixel boundaries and all detector cell boundaries to the source and calculating the x-intercepts. Based on these boundaries, the overlap length is calculated between each image pixel and each detector cell, and this overlap length is used to normalize the weights used in the projection and backprojection. This corresponds to applying a distance-driven kernel operation to the mapped boundary positions, which is achieved by looping over all boundaries (e.g., starting at the boundary at the minimum x-intercept and stopping at the boundary at the maximum x-intercept). Normalization consists of dividing by the pixel width (in the case of FBP) or the detector width (in the case of simulation).
[0086] The total system matrix is obtained by combining each modification. For example, A Tot =A geo A psf A fp A bh is. Finally, the operators (A and P) of the data fidelity term of the objective function for iterative reconstruction, shown in Equation (8), are the simple system geometry model A geo When combined, A geo A psf A fp A bh and P s become.
[0087]
Number
[0088] Next, the training of the DL network (for example, the offline DL training process 160) will be described in more detail. Here, the target data 153 is a high-quality CT image reconstructed using the IR method as described above, and the input data 157 is a low-quality CT image reconstructed using the fast analysis method.
[0089] Figure 3 shows a flowchart of the offline DL training process 160 according to an embodiment. In the offline DL training process 160, the input data 157 and the target data 153 are used as training data for training the DL network 170. As a result, the DL network 170 is output from step 319 of the offline DL training process 160. The offline DL training process 160 uses a large number of reconstructed CT images of the input data 157 paired with the corresponding physical model information 155 and the reconstructed CT image of the target data 153 for training the DL network 170, and generates a filtered CT image similar to the target data 153, which is the target CT image, from the input CT image that is the input data 157.
[0090] In the offline DL training process 160, a set of training data is acquired, and the DL network 170 is repeatedly updated to reduce the error (e.g., the value generated by the loss function). The DL network infers the mapping implied by the training data, and the loss function generates an error value related to the discrepancy between the CT image of the target data 153 and the result generated by applying the current form of the DL network 170 to the CT image of the input data 157.
[0091] For example, in some embodiments, the loss function can use the mean squared error to minimize the mean squared error. In the case of a multilayer perceptrons (MLP) neural network, the backpropagation algorithm can be used to train the network by minimizing the mean squared error-based loss function using (stochastic) gradient descent.
[0092] Thus, the generation circuit 517 performs the learning of the DL network 170, which is a neural network, by repeatedly adjusting the weight coefficients of the DL network 170 to minimize the value of the loss function, which is a measure of the discrepancy between the input data 157 and the target data 153. Here, for example, the training data set includes the input data 157, which is an image reconstructed using a first CT method such as the analytical reconstruction method, and the target data 153, which is an image reconstructed using a second CT method such as the IR method.
[0093] In step 316 of the offline DL training process 160, an initial guess for the coefficients of the DL network 170 is generated. For example, the initial guess can be based on prior knowledge of the region to be imaged or one or more exemplary noise removal methods, edge detection methods, and / or blob detection methods. Additionally, the initial guess can be made based on any of LeCun initialization, Xavier initialization, and Kaiming initialization.
[0094] Steps 316 to 319 of the offline DL training process 160 provide non-limiting examples of optimization methods for training the DL network 170.
[0095] An error is calculated (e.g., using a loss function or a loss function) to represent a measure of the difference (e.g., a distance measure) between the CT image of the target data 153 (i.e., the ground truth) and the CT image of the input data 157 after applying the current version of the DL network 170. This error can be calculated using any known loss function including the above loss function, or a distance measure between the image data. Further, in certain embodiments, the error / loss function can be calculated using one or more of the hinge loss and the cross-entropy loss.
[0096] In some embodiments, the loss function can be the l p norm of the difference between the target data and the result of applying the input data to the DL network 170. Different values of "p" of the l p norm can be used to emphasize different aspects of the noise. Further, a weighting mask (e.g., based on the attenuation coefficient of the signal intensity) can be applied pixel-by-pixel to the difference between the target data and the result generated from the input data. In some embodiments, instead of minimizing the l p norm of the difference between the target data and the result from the input data, the loss function can represent similarity (e.g., using the Peak Signal-To-Noise Ratio (PSNR) or the Structural Similarity (SSIM) index).
[0097] That is, the loss function may include the peak signal-to-noise ratio, the structural similarity index, or the l p norm of the difference between the image obtained from the input data and the image of the target data.
[0098] In certain embodiments, the training is performed by minimizing the loss function represented by the following equation (9).
[0099]
Number
[0100] Here, θ is an adjustable weight coefficient of the DL network 170, h is a non-adjustable parameter (e.g., a parameter selected by the user such as the selection of a reconstruction kernel), y (n) represents the n-th input CT image, and y’ (n) represents the n-th target CT image. The number N is the total number of training projections. In some embodiments, a weighted mean absolute error loss function L represented by the following equation (10) is used.
[0101]
Number
[0102] Here, d j is a weight having the following form as shown in equation (11).
[0103]
Number
[0104] Note that p is a scalar quantity.
[0105] The selection of this weight is triggered by the statistical mean estimation method, and d here j is often selected to be the reciprocal of the data noise variance. To address the problem of overfitting, an additional regularization R of h is used, which is R(h)=Σ j h j given by. The regularization strength can be adjusted via the parameter β.
[0106] In some embodiments, the DL network 170 is trained using backpropagation. Backpropagation can be used for training neural networks and is used in combination with the gradient descent optimization method. During the forward pass, the algorithm calculates the network's predictions based on the current parameters θ. These predictions are then input into a loss function, which thereby compares them to the corresponding ground truth labels (i.e., high-quality target data 153). During the backward pass, the model calculates the gradient of the loss function with respect to the current parameters, and then the parameters are updated by taking a step size of a predetermined size in the direction of the minimized loss (in acceleration methods such as Nesterov's momentum method and various adaptive methods, the step size can be selected to converge more quickly and optimize the loss size).
[0107] In the optimization method that performs backprojection, one or more of gradient descent, batch gradient descent, stochastic gradient descent, and mini-batch stochastic gradient descent can be used. The forward pass and the backward pass can be performed step by step through each layer of the network. In the forward pass, the execution starts by supplying the input through the first layer, thereby activating the output for the subsequent layer, and so on. This process is repeated until the loss function of the last layer is reached. During the backward pass, the last layer calculates the gradient with respect to its own trainable parameters (if any) and its own input that acts as the upstream derivative of the previous layer. This process is repeated until the input layer is reached.
[0108] Returning to FIG. 3, step 317 of the offline DL training process 160 can calculate a function of the changes in the network (e.g., error gradient) to determine the change in the error, and use this change in the error to select a direction and step size for subsequent changes to the weights / coefficients of the DL network 170. Calculating the gradient of the error in this way is consistent with some embodiments of the gradient descent optimization method. As will be appreciated by those skilled in the art, in some other embodiments this step may be omitted and / or replaced with another step that follows another optimization algorithm (e.g., a non-gradient descent optimization algorithm such as the simulated annealing method or genetic algorithm).
[0109] In step 317 of the offline DL training process 160, a new set of coefficients is determined for the DL network 170. For example, the weights / coefficients can be updated using the changes calculated in step 317 in the same way as they are updated in the gradient descent optimization method or the over-relaxation acceleration method.
[0110] In step 318 of the offline DL training process 160, a new error value is calculated using the updated weights / coefficients of the DL network 170.
[0111] In step 319, a predetermined stopping criterion is used to determine whether the training of the network is complete. For example, the predetermined stopping criterion can evaluate whether the new error and / or the total number of executed iterations exceeds a predetermined value. For example, when the new error falls below a predetermined threshold or the maximum number of iterations is reached, the stopping criterion can be satisfied. If the stopping criterion is not satisfied, the training process executed in the offline DL training process 160 returns to step 317 using the new weights and coefficients, and repeats this step to return to the start of the iterative loop (the iterative loop includes steps 317, 318, and 319). When the stopping criterion is satisfied, the training process executed in the offline DL training process 160 is completed.
[0112] Figures 4 and 5 show two examples of the interconnections between layers in the DL network 170. The DL network 170 can include fully connected convolutional layers and pooling layers, all of which will be described below. In a particular preferred embodiment of the DL network 170, the convolutional layers are arranged near the input layer, while the fully connected layers that perform high-level inferences are arranged further downstream in the architecture towards the loss function. Inserting a pooling layer after the convolutional layer has proven to reduce the spatial extent of the filters and thus the amount of learnable parameters. Activation functions are also incorporated into the various layers to introduce non-linearity and also to enable the network to learn complex prediction relationships. The activation function is a saturation activation function (e.g., sigmoid or hyperbolic tangent activation function), or a normalization activation function (e.g., the Rectified Linear Unit (ReLU) applied in the first and second examples above). The layers of the DL network 170 can also incorporate batch normalization, as illustrated in the first and second examples described above.
[0113] Figure 4 shows an example of a general artificial neural network (ANN) with N inputs, K hidden layers, and three outputs. Each layer is composed of nodes (also called neurons), and each node performs a weighted sum of the inputs and compares the result of this weighted sum with a threshold to generate an output. The ANN constitutes a class of functions, for which members of the class are obtained by varying architectural details such as the threshold, the weights of the connections, or the number of nodes and / or their connections. The nodes of the ANN are called neurons (or neuron nodes), and the neurons can have interconnections between different layers of the ANN system. Synapses (i.e., connections between neurons) store values called "weights" (also called "coefficients" or "weight coefficients") that manipulate data in calculations. The output of the ANN depends on three types of parameters, namely, (i) the interconnection pattern between different layers of neurons, (ii) the learning process for updating the weights of the interconnections, and (iii) the activation function that converts the weighted input of the neurons into their output activation.
[0114] Mathematically, the neural network function m(x) of a neuron is defined as the composition of other function(s) n i (x), which can in turn be defined as the composition of other functions. This can be conveniently represented as a network structure using arrows to represent the dependency relationships between variables, as shown in Figure 4. For example, the ANN can use a non-linear weighted sum, such as that represented by Equation (12).
[0115]
Number
[0116] Here, K (generally referred to as the activation function) is a predefined function, such as the hyperbolic tangent.
[0117] In FIG. 4 (and similarly in FIG. 5), neurons (i.e., nodes) are represented as circles around a threshold function. In the non-limiting example shown in FIG. 4, inputs are depicted as circles surrounding a linear function, and the arrows indicate directional connections between neurons. In some embodiments, DL network 170 is a feed-forward network.
[0118] FIG. 5 shows a non-limiting example where DL network 170 is a convolutional neural network (CNN). A CNN is a type of ANN that has properties beneficial for image processing and is thus particularly relevant for image noise removal applications. A CNN uses a feed-forward ANN, where the connection pattern between neurons can represent the convolution of image processing. For example, a CNN can be used for image processing optimization by using multiple layers of small neuron collections called receptive fields that process a part of the input image. The outputs of these collections can then be arranged and displayed overlapping to obtain a better representation of the original image. This processing pattern can be repeated across multiple layers having alternating convolutional and pooling layers.
[0119] After the convolutional layer, a CNN can include local and / or global pooling layers, which combine the outputs of the neuron clusters in the convolutional layer. Additionally, in certain embodiments, various combinations of fully connected layers with convolutional layers can be included in the CNN, and point-wise non-linearity is applied at or after the end of each layer.
[0120] Note that the embodiments are not limited to the above examples, and an image processing apparatus including the acquisition circuit 518, reconstruction circuit 514, and generation circuit 517 in FIG. 6 may be configured as an independent image processing apparatus.
[0121] Such an image processing apparatus includes an acquisition circuit 518, a reconstruction circuit 514, and a generation circuit 517. The acquisition circuit 518 acquires X-ray projection data as raw data 105. The reconstruction circuit 514 reconstructs a first image, which is an input image 112, from the X-ray projection data based on a first CT method. Also, the generation circuit 517 generates data of a physical model in step 115 based on at least one of the projection data or the first image. The generation circuit 517 inputs the first image and the data of the physical model to a DL network 170, which is a neural network trained using a training data set and physical model information 155, to generate a second image, which is a final image 135. The training data set includes input data 157, which is an image reconstructed using the first CT method, and target data 153, which is an image reconstructed using a second CT method, and the physical model information 155 is generated in the same manner as the generation of the data of the physical model in step 115.
[0122] Also, a program executed on a computer may cause the computer to execute steps of acquiring X-ray projection data as raw data 105, reconstructing a first image, which is an input image 112, from the X-ray projection data based on a first CT method, generating data of a physical model in step 115 based on at least one of the projection data or the first image, and inputting the first image and the data of the physical model to a DL network 170, which is a neural network trained using a training data set and physical model information 155, to generate a second image, which is a final image 135. Such a computer may include, for example, a CPU (processing circuit) that can be implemented as individual logic gates, an application specific integrated circuit (ASIC), a field programmable gate array (FPGA), or other complex programmable logic device (CPLD).
[0123] According to at least one embodiment described above, the image quality can be improved.
[0124] Although several embodiments have been described, these embodiments are presented by way of example and are not intended to limit the scope of the invention. These embodiments can be implemented in various other forms, and various omissions, replacements, changes, and combinations of embodiments can be made without departing from the gist of the invention. These embodiments and their modifications are included in the scope and gist of the invention, as well as in the invention described in the claims and its equivalent scope.
Description of Reference Numerals
[0125] 514 Reconfiguration circuit 517 Generation circuit 518 Acquisition circuit
Claims
1. An acquisition unit that acquires X-ray projection data, A reconstruction unit that reconstructs a first image from the X-ray projection data based on the first CT method, Generates data of a physical model based on at least one of the X-ray projection data or the first image A generation unit that generates a second image by inputting the first image and the data of the physical model to a neural network that has been trained using a training data set and physical model information, The training data set includes input data that is an image reconstructed using the first CT method and target data that is an image reconstructed using the second CT method, The physical model information is generated in the same way as the generation of the data of the physical model, and the physical model includes at least one of an X-ray scattering model, a beam hardening model, a spatial resolution detector model, a forward projection model, or a system geometry model, The first CT method is an analytical reconstruction method, The second CT method is an MBIR (Model-Based Iterative Reconstruction) method, The second image has better image quality than the first image, The reconstruction unit reconstructs the first image using the analytical reconstruction method, The generation unit generates the data of the physical model using the MBIR (Model-Based Iterative Reconstruction) method, and performs the learning of the neural network by repeatedly adjusting the weight coefficients of the neural network to minimize the value of a loss function that is a measure of the discrepancy between the output of the neural network when the input data and the physical model information are input and the target data, The generation unit generates an estimated value of the gradient of the objective function of the MBIR method as the data of the physical model by executing one iteration of the MBIR method, an X-ray diagnostic system.
2. The X-ray diagnostic system according to claim 1, wherein the generation unit updates an image based on a data fidelity term of the objective function of the MBIR method.
3. The neural network is a residual network, The generation unit generates the second image by subtracting the output result of the neural network from the first image, The X-ray diagnostic system according to claim 1.
4. The loss function includes a peak signal-to-noise ratio, a structural similarity index, or an l p norm of a difference between an image obtained from the input data and an image of the target data, the X-ray diagnostic system according to claim 1.
5. An X-ray source that emits X-rays, A detector array that is arranged diametrically opposite to the X-ray source across the opening of the gantry and detects the X-rays emitted from the X-ray source to generate the X-ray projection data The X-ray diagnostic system according to claim 1, further comprising.
6. An acquisition unit that acquires X-ray projection data; Reconstructing a first image from the X-ray projection data based on a first CT method, Generating data of a physical model based on at least one of the X-ray projection data or the first image A generation unit that generates a second image by inputting the first image and the data of the physical model to a neural network that has been trained using a training data set and physical model information, The training data set includes input data that is an image reconstructed using a first CT method and target data that is an image reconstructed using a second CT method, The physical model information is generated in the same manner as the generation of the data of the physical model, and the physical model includes at least one of an X-ray scattering model, a beam hardening model, a spatial resolution detector model, a forward projection model, or a system geometry model, The first CT method is an analytical reconstruction method, The second CT method is an MBIR (Model-Based Iterative Reconstruction) method, The second image has better image quality than the first image, The generation unit reconstructs the first image using the analytical reconstruction method, generates the data of the physical model using the MBIR (Model-Based Iterative Reconstruction) method, and performs the learning of the neural network by repeatedly adjusting the weight coefficients of the neural network to minimize the value of a loss function that is a measure of the discrepancy between the output of the neural network when the input data and the physical model information are input and the target data, The generation unit generates an estimated value of the gradient of the objective function of the MBIR method as the data of the physical model by executing one iteration of the MBIR method, An image processing apparatus.
7. A step of acquiring X-ray projection data; A step of reconstructing a first image from the X-ray projection data based on a first CT method A step of generating data of a physical model based on at least one of the X-ray projection data or the first image A program for causing a computer to execute a step of generating a second image by inputting the first image and the data of the physical model to a neural network trained using a training data set and physical model information, wherein the training data set includes input data that is an image reconstructed using a first CT method and target data that is an image reconstructed using a second CT method, the physical model information is generated in the same manner as the generation of the data of the physical model, and the physical model includes at least one of an X-ray scattering model, a beam hardening model, a spatial resolution detector model, a forward projection model, or a system geometry model, the first CT method is an analytical reconstruction method, the second CT method is a MBIR (Model-Based Iterative Reconstruction) method, the second image has better image quality than the first image, the reconstructing step reconstructs the first image using the analytical reconstruction method, the generating step generates the data of the physical model using the MBIR (Model-Based Iterative Reconstruction) method, and performs the learning of the neural network by repeatedly adjusting the weight coefficients of the neural network to minimize the value of a loss function that is a measure of the discrepancy between the output of the neural network when the input data and the physical model information are input and the target data, the step of generating the second image generates an estimated value of the gradient of the objective function of the MBIR method as the data of the physical model by executing one iteration of the MBIR method, Program.
Citation Information
Patent Citations
Iterative reconstruction with system optics modeling using filters
US20160300369A1
Apparatus and method of iterative image reconstruction using regularization-parameter control
US20170294034A1
Deep learning based acceleration for iterative tomographic reconstruction
US20180197317A1
System and method for image conversion
US20190035118A1