Free breathing rapid three-dimensional heart magnetic resonance imaging method and system

Abstract

The invention discloses a free-breathing rapid three-dimensional heart magnetic resonance imaging method and system, and belongs to the field of three-dimensional imaging, and the method comprises the steps: collecting under-sampling K space data of a heart in different breathing states; constructing a motion solving network and an image reconstruction network based on implicit neural expression; according to sampling K space data in different breathing states, iterative updating and iterative reconstruction are carried out on a motion solving network and an image reconstruction network through a loss function which comprises regular term constraints in multiple aspects of a motion compensation data fidelity term, motion field smoothness, image smoothness and image edge sharpness; a breathing motion field estimation result and a heart nuclear magnetic resonance motion compensation reconstruction image are obtained; according to the method, the respiratory motion field and the three-dimensional heart magnetic resonance image are directly reconstructed from the under-sampled K space data by means of implicit neural expression, multiple regular term constraints are combined, the imaging time is shortened, the signal-to-noise ratio of the image is improved, the image reconstruction quality can be guaranteed, and extra training data are not needed.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of three-dimensional cardiac magnetic resonance imaging, and particularly relates to a free-breathing fast three-dimensional cardiac magnetic resonance imaging method and system, and particularly relates to a free-breathing fast three-dimensional cardiac magnetic resonance imaging method based on implicit neural representation. BACKGROUND

[0002] Three-dimensional cardiac magnetic resonance imaging can obtain high-resolution three-dimensional images of the heart, and further realize non-invasive assessment of the whole heart anatomical structure and function, thereby providing a basis for diagnosis and treatment of heart diseases. However, compared with conventional two-dimensional imaging, three-dimensional cardiac magnetic resonance imaging faces some key challenges that hinder its wide clinical application. First, limited by the slow acquisition speed of magnetic resonance and the large sampling matrix of three-dimensional imaging, three-dimensional high-resolution cardiac magnetic resonance scanning requires a long scanning time, which is challenging for patients who cannot maintain stillness for a long time. High-oversampling acquisition is required in clinical practice to shorten the scanning time. Second, respiratory motion during three-dimensional cardiac scanning often leads to significant respiratory artifacts and reduces image quality. To deal with respiratory motion in free-breathing magnetic resonance scanning, diaphragm navigation technology and self-navigation technology have been proposed. Among them, the self-navigation technology uses data collected at different respiratory states to achieve a scanning efficiency close to 100% and a predictable scanning time. Because the collected data is scattered in different motion states, image registration is generally needed in the image reconstruction process to estimate the respiratory motion field between different motion states, which is then used for subsequent motion compensation image reconstruction. However, due to high-oversampling acquisition, the images at different respiratory states often contain significant oversampling artifacts, which reduces the registration accuracy and affects the motion compensation reconstruction quality.

[0003] The methods for motion compensation image reconstruction of under-sampled three-dimensional magnetic resonance imaging data containing respiratory motion proposed in the field at the present stage mainly include two categories. The first category is the traditional iterative method, including NR-SENSE, NR-PROST and the like. This category of method first solves the three-dimensional non-rigid respiratory motion field between different respiratory states by using a free deformation non-rigid registration algorithm based on the smoothness constraint of the motion field; then the solved motion field is substituted into the loss function, and the motion compensation reconstruction image is obtained by minimizing the loss function through the iterative algorithm. This kind of method performs well on low-fold under-sampled data, and as the under-sampling multiple increases, the motion estimation accuracy decreases, and the reconstructed image appears obvious artifacts and noise amplification. The second category is the motion compensation reconstruction method based on deep learning, such as MoCo-MoDL. The network model constructed by this category of method includes a motion estimation network and an image reconstruction network, which first estimates the three-dimensional non-rigid respiratory motion field from the different respiratory state images containing artifacts by using the motion estimation network, and then substitutes the non-rigid respiratory motion field into the image reconstruction network to output the motion compensation reconstruction result. This kind of method can realize higher under-sampling multiple than the traditional iterative method by training on a large amount of magnetic resonance data. However, in most actual situations, it is difficult to collect a large-scale three-dimensional cardiac magnetic resonance training data set; and when the test set and the training set are inconsistent, the trained neural network may have generalization problems, resulting in a significant decrease in reconstruction quality. In order to avoid using large-scale training data, some self-supervised image reconstruction methods are proposed, such as DIP, DIP-CS and the like. However, this kind of method can only deal with the under-sampled data reconstruction problem without motion, that is, it can only reconstruct the under-sampled K-space data in the same respiratory motion state, and cannot utilize the data collected in all respiratory states for motion compensation image reconstruction, which has low data utilization and is difficult to realize high-fold acceleration. Therefore, a scheme for high-fold under-sampling motion compensation reconstruction without training data is needed to improve the clinical usability of three-dimensional cardiac magnetic resonance imaging. SUMMARY

[0004] To solve the above technical problems, the present application provides a free breathing fast three-dimensional cardiac magnetic resonance imaging method and system to solve the problems existing in the prior art.

[0005] To achieve the above purpose, the present application provides a free breathing fast three-dimensional cardiac magnetic resonance imaging method and system, comprising: Acquiring under-sampled K-space data of the heart in different respiratory states, wherein the under-sampled K-space data is under-sampled Fourier space cardiac magnetic resonance data; Constructing a motion solving network and an image reconstruction network based on implicit neural representation; Constructing a loss function for joint motion estimation and motion compensation image reconstruction; According to the undersampled K-space data in different respiratory states, the motion solving network and the image reconstruction network are iteratively updated and iteratively reconstructed through a loss function, to obtain a final respiratory motion field estimation result and a motion compensation reconstructed image. The loss function includes different loss terms of data fidelity, data smoothing and image detail constraint.

[0006] Optionally, the acquisition process of the undersampled K-space data includes: The magnetic resonance imaging data of the free-breathing heart is acquired through a self-navigation imaging sequence; and for the magnetic resonance imaging data, the data of the center point of the K-space in the magnetic resonance imaging data in each cardiac cycle acquisition window is obtained as navigation data. A respiratory curve is generated according to the navigation data, the respiratory curve is divided according to the respiratory amplitude of the respiratory curve to obtain a respiratory phase, wherein each respiratory phase corresponds to a respiratory state from the end of inspiration to the end of expiration. The undersampled K-space data in different respiratory states is obtained by integrating the respiratory phase with the magnetic resonance imaging data.

[0007] Optionally, the loss function includes a motion compensation data fidelity term, a motion field smoothing regular term, an image smoothing regular term and an image edge sharpness regular term; wherein the motion compensation data fidelity term is a loss term between the true undersampled K-space data and the K-space data synthesized by the motion compensation reconstructed image, the motion field smoothing regular term is a smoothing measurement loss term of the respiratory motion field, the image smoothing regular term is a smoothing measurement loss term of the motion compensation reconstructed image, and the image edge sharpness regular term is an edge sharpness loss term of the motion compensation reconstructed image.

[0008] Optionally, the loss function is a weighted sum result of the motion compensation data fidelity term, the motion field smoothing regular term, the image smoothing regular term and the image edge sharpness regular term.

[0009] Optionally, the motion solving network adopts a multi-layer perception machine, wherein the motion solving network includes 3 layers of hidden layers, the hidden neurons of each hidden layer are 32, 32 and 3, the first layer of hidden layer and the second layer of hidden layer adopt an Elu activation function, the length of the output layer is 3, and the input is the encoding information of the respiratory motion field, wherein the encoding information of the respiratory motion field is obtained by hashing encoding the spatial coordinates of the voxel points of the respiratory motion field and the respiratory state, and the output result is the respiratory motion field estimation result, wherein the respiratory motion field estimation result is the displacement amount of each voxel point in different directions.

[0010] Optionally, the image reconstruction network employs a multilayer perceptron, wherein the image reconstruction network includes 6 hidden layers, the input is the encoded information of the reconstructed image, which is obtained by hashing the spatial coordinates of the voxel points of the reconstructed image; the output is the real and imaginary parts of the voxel signals of the voxel points of the reconstructed image.

[0011] Optionally, the process of iteratively updating and reconstructing the motion solving network and the image reconstruction network includes: The loss function is optimized by a semi-quadratic splitting method, and by introducing auxiliary variables, the optimization steps of the loss function for the motion solving network and the image reconstruction network are split into a first sub-problem and a second sub-problem. For the first subproblem, the spatial coordinates of voxel points in the reconstructed image and motion field are obtained. The motion solving network and the image reconstruction network are used to reconstruct the motion field based on the spatial coordinates of voxel points in the reconstructed image and motion field, thus obtaining the breathing motion field estimation result and the motion-compensated reconstructed image for the first subproblem. The motion solving network and the image reconstruction network are then optimized based on the breathing motion field estimation result and the motion-compensated reconstructed image for the first subproblem. For the second subproblem, the motion-compensated reconstructed image of the first subproblem is denoised and sharpened to obtain the processed image of the second subproblem; the motion solving network and the image reconstruction network are optimized based on the processed image of the second subproblem. The optimization steps of the first and second subproblems are executed iteratively and alternately to obtain the final respiratory motion field estimation results and motion compensation reconstruction images.

[0012] Optionally, a loss function is used to optimize the motion solving network and the image reconstruction network for the first subproblem. for: Where z represents an auxiliary variable, This indicates the respiratory phase, or the respiratory state. Indicates the first K-space data collected from each respiratory state. Indicates the first K-space sampling mask corresponding to each breathing state Indicates Fourier transform, Indicates coil sensitivity. The first part represents the estimation of the motion solution network. The respiratory motion field corresponding to each breathing state This represents the motion-compensated reconstructed image estimated by the image reconstruction network. Represents the finite difference operator. , , denote weight coefficients corresponding to different terms.

[0013] Optionally, the loss function for optimizing the motion solving network and the image reconstruction network for the second sub-problem is: wherein, denotes an image edge sharpness loss term, denotes corresponding weight coefficients.

[0014] In another aspect, the application provides a free breathing fast three-dimensional cardiac magnetic resonance imaging system for performing the above method.

[0015] Compared with the prior art, the application has the following advantages and technical effects: 1. The application jointly solves the respiratory motion field and the motion compensation three-dimensional cardiac magnetic resonance image directly from the undersampled K-space data by means of implicit neural representation, avoids additional reconstruction and registration of high-undersampled images of each respiratory state, prevents the propagation of motion solving errors, and improves the quality of high-undersampling motion compensation reconstruction; 2. The application combines multiple regularization constraints, improves the signal-to-noise ratio of the reconstructed image, and effectively preserves the details and image edges of the reconstructed image; 3. No other training data is needed except for the undersampled K-space data itself, which can be generalized to different imaging modalities and undersampling rates. BRIEF DESCRIPTION OF DRAWINGS

[0016] The accompanying drawings, which form a part of this application, are included to provide a further understanding of the application and are incorporated in and constitute a part of this application. The embodiments of this application and their Figure 1 is an implementation flowchart of the free breathing fast three-dimensional cardiac magnetic resonance imaging joint motion estimation and motion compensation image reconstruction method based on implicit neural representation of the embodiment of the application. Figure 2 is a technical principle diagram of the free breathing fast three-dimensional cardiac magnetic resonance imaging joint motion estimation and motion compensation image reconstruction method based on implicit neural representation of the embodiment of the application. Figure 3 is a result schematic diagram of the free breathing fast three-dimensional cardiac magnetic resonance imaging joint motion estimation and motion compensation image reconstruction method based on implicit neural representation of an embodiment of the embodiment of the application. DETAILED DESCRIPTION

[0017] It should be noted that the embodiments in the present application and the features in the embodiments can be combined with each other without conflict. The present application will be described in detail below with reference to the accompanying drawings and in combination with the embodiments.

[0018] It is noted that the steps shown in the flowcharts of the drawings can be performed in a computer system such as a set of computer-executable instructions, and although a logical order is shown in the flowcharts, in some cases the steps shown or described can be performed in an order different from that shown here.

[0019] The application provides a free breathing fast three-dimensional cardiac magnetic resonance imaging method and system, comprising: acquiring undersampled K-space (Fourier space) data under each breathing state by using a self-navigation imaging sequence; constructing a loss function including motion compensation data fidelity term, motion field smoothness regularization term, image smoothness regularization term and image edge sharpness regularization term constraint as an optimization target; constructing an iterative reconstruction algorithm based on implicit neural representation to minimize the loss function, and jointly solving the breathing motion field and the motion compensation reconstructed image. The application directly solves the breathing motion field and the three-dimensional cardiac magnetic resonance image from the undersampled K-space data by means of implicit neural representation, avoids the additional reconstruction and registration of each breathing state image in the traditional method, and combines various regularization term constraints to reduce the imaging time and improve the image signal-to-noise ratio, so as to accelerate the three-dimensional cardiac magnetic resonance imaging data acquisition while ensuring the image reconstruction quality, and without additional training data.

[0020] In order to describe the above technical solutions in detail, the following embodiments are provided for detailed description: Embodiment 1: With reference to Figure 1 and Figure 2 , the application provides a free breathing fast three-dimensional cardiac magnetic resonance imaging method, which comprises the following steps: Step S1: acquiring undersampled K-space data under each breathing state by using a self-navigation imaging sequence; Step S1.1: acquiring free breathing three-dimensional cardiac magnetic resonance imaging data by using a self-navigation sequence; Step S1.2: generating a breathing curve according to the self-navigation signal of the acquired data; Step S1.3: dividing a plurality of breathing phase according to the breathing amplitude of the breathing curve, each breathing phase corresponding to a breathing state from the end of inspiration to the end of expiration; Step S1.4: extracting all the data acquired under each breathing phase respectively to obtain the undersampled K-space data under the breathing state.

[0021] Step S2: constructing a loss function including motion compensation data fidelity term, motion field smoothness regularization term, image smoothness regularization term and image edge sharpness regularization term constraint as an optimization target.

[0022] Step S3: constructing an iterative reconstruction algorithm based on implicit neural representation to minimize a loss function, jointly solving a respiratory motion field and a motion compensation reconstructed image; using a Half-Quadratic Splitting algorithm to optimize the loss function, introducing an auxiliary variable z to split the optimization problem into two sub-problems; using implicit neural representation to solve sub-problem 1 (first sub-problem), as follows: Step S3.1: constructing a motion solving network and an image reconstruction network, the motion solving network being used to output a motion field matrix corresponding to the spatial coordinates of the voxel points of the motion field, and the image reconstruction network being used to output an image matrix corresponding to the spatial coordinates of the voxel points of the three-dimensional cardiac magnetic resonance image; Step S3.2: defining the spatial coordinates of the voxel points of the motion field and the magnetic resonance image in advance; Step S3.3: taking the optimization objective of sub-problem 1 as a loss function for network training, and training the motion solving network and the image reconstruction network based on the defined spatial coordinates; Step S3.4: inputting the spatial coordinates of the voxel points of the motion field and the magnetic resonance image into the trained motion solving network and image reconstruction network respectively to obtain the respiratory motion field estimation result of sub-problem 1 and the motion compensation reconstructed image.

[0023] The solution of sub-problem 2 (second sub-problem) is obtained by sequentially performing image denoising and image sharpening on the motion compensation reconstructed image of sub-problem 1. By alternately solving sub-problem 1 and sub-problem 2, the respiratory motion solving result and the motion compensation reconstructed image of step S3 are obtained.

[0024] The application also provides a free-breathing fast three-dimensional cardiac magnetic resonance imaging system, which can be realized by performing the flow steps of the free-breathing fast three-dimensional cardiac magnetic resonance imaging method, that is, the free-breathing fast three-dimensional cardiac magnetic resonance imaging method can be understood by those skilled in the art as a preferred embodiment of the free-breathing fast three-dimensional cardiac magnetic resonance imaging system.

[0025] Example 2: The application also provides a free-breathing fast three-dimensional cardiac magnetic resonance imaging system, which comprises the following modules: Module M1: acquiring under-sampled K-space data in each respiratory state using a self-navigation imaging sequence; Module M1.1: acquiring free-breathing three-dimensional cardiac magnetic resonance imaging data using a self-navigation sequence; Module M1.2: generating a respiratory curve according to the self-navigation signal of the acquired data; Module M1.3: dividing a plurality of respiratory phase according to the respiratory amplitude of the respiratory curve, each respiratory phase corresponding to a respiratory state from the end of inspiration to the end of expiration; Module M1.4: extracting all data collected in each respiratory phase to obtain the undersampled K-space data in the respiratory state.

[0026] Module M2: constructing a loss function including a motion compensation data fidelity term, a motion field smoothness regularization term, an image smoothness regularization term and an image edge sharpness regularization term as an optimization objective.

[0027] Module M3: constructing an iterative reconstruction algorithm based on implicit neural representation to minimize the loss function, and jointly solving the respiratory motion field and the motion compensation reconstructed image; The above loss function is optimized by using a Half-Quadratic Splitting algorithm, and by introducing an auxiliary variable z, the optimization problem is split into two sub-problems; sub-problem 1 is solved by using implicit neural representation, including the following modules: Module M3.1: constructing a motion solving network and an image reconstruction network, the motion solving network is used to output a motion field matrix corresponding to the spatial coordinates of the voxel points of the motion field, and the image reconstruction network is used to output an image matrix corresponding to the spatial coordinates of the voxel points of the three-dimensional cardiac magnetic resonance image; Module M3.2: defining the spatial coordinates of the voxel points of the motion field and the magnetic resonance image in advance; Module M3.3: taking the optimization objective of sub-problem 1 as the loss function of network training, and training the motion solving network and the image reconstruction network based on the defined spatial coordinates; Module M3.4: inputting the spatial coordinates of the voxel points of the motion field and the magnetic resonance image into the trained motion solving network and the image reconstruction network respectively to obtain the respiratory motion field estimation result of sub-problem 1 and the motion compensation reconstructed image.

[0028] The solution of sub-problem 2 is obtained by sequentially performing image denoising and image sharpening on the motion compensation reconstructed image of sub-problem 1. By alternately solving sub-problem 1 and sub-problem 2, the respiratory motion solving result and the motion compensation reconstructed image of the module M3 are obtained.

[0029] Embodiment 3: An embodiment of the application applied to three-dimensional cardiac magnetic resonance coronary imaging will be described in detail below with reference to the accompanying drawings. Figure 1 The flowchart of the embodiment of the application is as follows, Figure 2 The algorithm principle diagram of the application is as follows. Specifically, the embodiment includes the following steps: Step S1: acquiring undersampled K-space data in each respiratory state by using a self-navigation imaging sequence: Step S1.1: The self-navigated gradient echo coronary imaging sequence acquires three-dimensional cardiac magnetic resonance imaging data during free breathing using the following method: The sequence employs ECG gating, applying a T2 preparation pulse after the R wave on the ECG to enhance myocardial blood contrast. The T2 preparation time is 32 ms, after which imaging data is acquired. The timing of the T2 preparation pulse application and the data acquisition time window are determined based on cardiac cine images. The sequence uses a pseudo-radial sampling trajectory to acquire K-space data. Within the acquisition window of each cardiac cycle, data at the K-space center point is acquired first as navigation data and also as imaging data. Then, other imaging data are acquired from the K-space center outwards.

[0030] Step S1.2: Extract the navigation data collected for each heartbeat, perform a one-dimensional inverse Fourier transform on it to obtain a one-dimensional projection image of the imaging area along the head-to-toe direction, then select the one-dimensional projection image of the last heartbeat at the end of expiration as the reference image, and perform correlation calculations between it and the one-dimensional projection images of all heartbeats to obtain the respiratory curve; Step S1.3: Divide the respiratory curve into four respiratory phases at equal intervals based on the respiratory amplitude, corresponding to the four respiratory states from the end of inspiration to the end of expiration from bottom to top; among them, the four respiratory states are divided into four respiratory states at equal intervals from the end of inspiration (the lowest point of the respiratory curve) to the end of expiration (the highest point of the respiratory curve) based on the respiratory amplitude.

[0031] Step S1.4: Divide the data collected for each heartbeat into the respiratory phase in which the heartbeat occurs, and merge all K-space data collected under the four respiratory phases to obtain the undersampled K-space data for the corresponding respiratory state.

[0032] Step S2: Construct a loss function as the optimization objective, which includes constraints on motion compensation data fidelity, motion field smoothness regularization, image smoothness regularization, and image edge sharpness regularization. (1) Motion compensation data fidelity item ( First, the images for each breathing state are calculated based on the motion-compensated reconstructed image x and the solved motion field U for each breathing state. Then, the reconstructed K-space data obtained after coil sensitivity encoding and Fourier transform of the images for each breathing state is calculated. Finally, the square of the 2-norm of the difference between the reconstructed K-space data and the actual acquired K-space data at the sampling location is calculated. in, Indicates the first K-space data collected from each respiratory state. Indicates the first K-space sampling mask corresponding to each breathing state Indicates Fourier transform, Indicates coil sensitivity. The first part represents the estimation of the motion solution network. The respiratory motion field corresponding to each breathing state This represents the motion-compensated reconstructed image estimated by the image reconstruction network.

[0033] (2) Regularization term for smoothness of sports field ( ): Calculate the smoothness measure of the motion field U for each breathing state, using the square of the 2-norm of the finite difference result: in, This represents the finite difference operator.

[0034] (3) Image smoothness regularization term ( ): Calculate the smoothness measure of the motion-compensated reconstructed image. The smoothness measure uses the 1-norm of the finite difference result, i.e., the total variation: (4) Image edge sharpness regularization term ( To improve the smoothness of motion-compensated reconstructed images while preserving high-frequency information such as image details and edges, this regularization term enhances edge sharpness in motion-compensated reconstructed images through image denoising and image sharpening algorithms. in, This represents the implicit image edge sharpness loss term.

[0035] The four losses above are weighted and summed using weighting coefficients. Then, the losses for the four breathing states are summed to form the complete loss function. The weighting coefficients are hyperparameters.

[0036] Step S3: Construct an iterative reconstruction algorithm based on implicit neural expression to minimize the loss function, and jointly solve for the respiratory motion field and motion-compensated reconstructed image. Specifically, a half-quadratic splitting algorithm is used to optimize the loss function. By introducing an auxiliary variable z, the optimization problem is split into two sub-problems: The loss function for subproblem 1 is: Where z represents an auxiliary variable, This indicates the respiratory phase, or the respiratory state. Indicates the first K-space data collected from each respiratory state. Indicates the first K-space sampling mask corresponding to each breathing state Indicates Fourier transform, Indicates coil sensitivity. a first respiratory state corresponding to the first respiratory state estimated by the motion estimation network, a respiratory motion field corresponding to the first respiratory state estimated by the motion estimation network, a motion-compensated reconstructed image estimated by the image reconstruction network, a finite difference operator, , , a weight coefficient corresponding to a different term.

[0037] The loss function of the sub-problem 2 optimization is: wherein, an image edge sharpness loss term, a weight coefficient corresponding to a different term.

[0038] wherein, it can be known from the principle of Half-Quadratic Splitting algorithm. The sub-problem 1 is to solve an optimization problem about the reconstructed image x, and the sub-problem 2 is to solve an optimization problem about the auxiliary variable z. The reconstructed image x in the sub-problem 2 is the reconstructed image x obtained after solving the sub-problem 1, and the auxiliary variable z is optimized on this basis, and the sub-problem 1 is solved with the current optimal auxiliary variable z in the next iteration process.

[0039] Step S3.1: constructing a joint motion estimation and motion-compensated reconstruction model based on implicit neural representation to solve the sub-problem 1 (first sub-problem) as shown in Figure 2 The model includes the following modules: (1) a respiratory motion estimation module: the input of the module is a 4D coordinate (x, y, z, b) composed of a 3D spatial coordinate and a 1D respiratory phase coordinate, and the output is the displacement amount of the voxel with a spatial coordinate (x, y, z) in the end-expiratory image along the left-right, up-down and front-back directions when the end-expiratory image is transformed to the bth respiratory phase. First, the input 4D coordinate is hashed to obtain a feature vector, and then the feature vector is input into the motion estimation network. The network architecture is a multi-layer perceptron with 3 layers, the number of hidden neurons is 32, 32, 3 respectively, and the activation function is Elu except for the last layer. The output length is 3, corresponding to the displacement amount along the left-right, up-down and front-back directions.

[0040] (2) Motion-compensated image reconstruction module: The input of this module is 3D spatial coordinates (x, y, z), and the output is the real and imaginary parts of the signal of the reconstructed three-dimensional cardiac magnetic resonance image voxel with spatial coordinates (x, y, z). First, the input 3D coordinates are hashed to obtain a feature vector, which is then input into the image reconstruction network. The network architecture is a multi-layer perceptron with 6 hidden layers, and the number of hidden neurons is 128, 128, 128, 128, 128, and 2 respectively. The activation function of all layers except the last layer is GELU. The output length is 2, corresponding to the real and imaginary parts of the reconstructed voxel signal.

[0041] Step S3.2: Define the motion field and the spatial coordinates of the magnetic resonance image voxel: The motion field and the magnetic resonance image are equally divided to obtain 3D spatial coordinates; the 4 breathing states are defined as 1D coordinates of (1), (2), (3), and (4), respectively; and the 4D coordinates of each breathing state motion field are composed of 3D spatial coordinates and 1D breathing state coordinates.

[0042] Step S3.3: Train the joint motion estimation and motion-compensated reconstruction model using the undersampled K-space data of each breathing state to be reconstructed: In each iteration, a motion field and a motion-compensated reconstruction image under a breathing state are randomly reconstructed. First, the 4D coordinates of the breathing state motion field selected in this iteration are input into the respiratory motion estimation module to obtain the motion field under the breathing state; the 3D spatial coordinates of the magnetic resonance image are input into the motion-compensated image reconstruction module to obtain the real and imaginary parts of the motion-compensated reconstruction image; then, the loss function corresponding to sub-problem 1 in step S3 is calculated according to the outputs of the relevant modules and the undersampled K-space data corresponding to the selected breathing state, and the learnable parameters of the motion solving network and the image reconstruction network and the hash table are updated by minimizing the loss function until the loss function no longer decreases, thereby obtaining the trained joint motion estimation and motion-compensated reconstruction model. In this embodiment, the Adam optimizer is used, the initial learning rate of the motion solving network is 0.001, the initial learning rate of the image reconstruction network is 0.01, the learning rate is exponentially decayed, the decay coefficient is 0.1, and a total of 500 iterations are trained.

[0043] Step S3.4: Input the 4D coordinates of each breathing state motion field voxel and the 3D spatial coordinates of the magnetic resonance image voxel into the trained motion solving network and image reconstruction network, respectively, to obtain the respiratory motion field estimation result of sub-problem 1 and the motion-compensated reconstruction image.

[0044] Step S3.5: Solve sub-problem 2 (second sub-problem): First, perform image denoising on the motion-compensated reconstruction image of sub-problem 1, and then perform image sharpening on the denoised image to obtain the solution of sub-problem 2. In this embodiment, the Non-Local Means algorithm is used to realize image denoising, and the Unsharp Masking algorithm is used to realize image sharpening.

[0045] Step S3.6: iteratively solve sub-problem 1 and sub-problem 2 according to the above steps to obtain the final motion estimation result and the motion compensated reconstructed image.

[0046] Figure 3 For the results of this embodiment, the PSNR is the peak signal-to-noise ratio, and the SSIM is the structural similarity index. Among them, the three rows correspond to the reconstruction results of three healthy volunteers; the first column is the image reconstructed by the traditional method NR-SENSE from 3-fold accelerated acquisition data as the reference; the second column is the image reconstructed by the method of the present application under 6-fold acceleration; the third and fourth columns are the reconstruction results of different control methods under 6-fold acceleration. Compared with the control methods, the image reconstructed by the method of the present application retains more complete details, has no obvious artifacts, and has obvious improvement in the peak signal-to-noise ratio and the structural similarity index relative to the reference image.

[0047] Based on the present application, the respiratory motion field can be estimated directly from the undersampled K-space data acquired under free breathing, and the motion compensated reconstructed image can be obtained. Compared with the existing methods, the present application does not need additional registration of magnetic resonance images containing artifacts, nor does it need additional training data, and is expected to improve the clinical usability of three-dimensional cardiac magnetic resonance imaging.

[0048] The above is only the preferred specific embodiment of the present application, but the protection scope of the present application is not limited thereto, any person skilled in the art can easily think of changes or replacements within the technical scope disclosed by the present application, which should be covered within the protection scope of the present application. Therefore, the protection scope of the present application should be subject to the protection scope of the claims.

Claims

1. A rapid three-dimensional cardiac magnetic resonance imaging method with free breathing, characterized in that, include: Undersampled K-space data of the heart under different respiratory states were collected, wherein the undersampled K-space data is undersampled Fourier space cardiac magnetic resonance data. Construct motion solving networks and image reconstruction networks based on implicit neural expression; Construct a loss function for joint motion estimation and motion-compensated image reconstruction; Based on the sampled K-space data under different respiratory states, the motion solving network and image reconstruction network are iteratively updated and reconstructed through a loss function to obtain the final respiratory motion field estimation result and the motion-compensated reconstruction image of cardiac magnetic resonance imaging. The loss function includes different loss terms for data fidelity, data smoothing, and image detail constraint.

2. The method according to claim 1, characterized in that, The process of acquiring the undersampled K-space data includes: Magnetic resonance imaging (MRI) data of the heart during free breathing was acquired using a self-navigation imaging sequence; for the MRI data, the center point data of the K-space was obtained from the MRI data within the acquisition window of each cardiac cycle as navigation data; A breathing curve is generated based on the navigation data. The breathing curve is divided according to the breathing amplitude of the breathing curve to obtain the breathing phases, where each breathing phase corresponds to a breathing state from the end of inspiration to the end of expiration. By integrating the magnetic resonance imaging data relative to the respiratory phase, undersampled K-space data under different respiratory states are obtained.

3. The method according to claim 1, characterized in that, The loss function includes a motion compensation data fidelity term, a motion field smoothness regularization term, an image smoothness regularization term, and an image edge sharpness regularization term; wherein the motion compensation data fidelity term is a loss term between the real undersampled K-space data and the K-space data synthesized from the motion compensation reconstructed image, the motion field smoothness regularization term is a smoothness measurement loss term of the breathing motion field, the image smoothness regularization term is a smoothness measurement loss term of the motion compensation reconstructed image, and the image edge sharpness regularization term is an edge sharpness loss term of the motion compensation reconstructed image.

4. The method according to claim 3, characterized in that, The loss function is a weighted sum of motion compensation data fidelity, motion field smoothness regularization, image smoothness regularization, and image edge sharpness regularization.

5. The method according to claim 1, characterized in that, The motion solving network employs a multilayer perceptron, which includes three hidden layers. Each hidden layer has 32, 32, and 3 hidden neurons. The first and second hidden layers use the Elu activation function. The output layer has a length of 3. The input is the encoded information of the respiratory motion field, which is obtained by hashing the spatial coordinates of the voxel points of the respiratory motion field and the respiratory state. The output is the estimated result of the respiratory motion field, which is the displacement of each voxel point in different directions.

6. The method according to claim 1, characterized in that, The image reconstruction network employs a multilayer perceptron, which includes six hidden layers. The input is the encoded information of the reconstructed image, which is obtained by hashing the spatial coordinates of the voxel points in the reconstructed image. The output consists of the real and imaginary parts of the voxel signals of the voxel points in the reconstructed image.

7. The method according to claim 1, characterized in that, The process of iteratively updating and reconstructing the motion solving network and the image reconstruction network includes: The loss function is optimized by a semi-quadratic splitting method, and by introducing auxiliary variables, the optimization steps of the loss function for the motion solving network and the image reconstruction network are split into a first sub-problem and a second sub-problem. For the first subproblem, the spatial coordinates of voxel points in the reconstructed image and motion field are obtained. The motion solving network and the image reconstruction network are used to reconstruct the motion field based on the spatial coordinates of voxel points in the reconstructed image and motion field, thus obtaining the breathing motion field estimation result and the motion-compensated reconstructed image for the first subproblem. The motion solving network and the image reconstruction network are then optimized based on the breathing motion field estimation result and the motion-compensated reconstructed image for the first subproblem. For the second subproblem, the motion-compensated reconstructed image of the first subproblem is denoised and sharpened to obtain the processed image of the second subproblem; the motion solving network and the image reconstruction network are optimized based on the processed image of the second subproblem. The optimization steps of the first and second subproblems are executed iteratively and alternately to obtain the final respiratory motion field estimation results and motion compensation reconstruction images.

8. The method according to claim 1, characterized in that, The loss function for optimizing the motion solving network and image reconstruction network for the first subproblem. for: Where z represents an auxiliary variable, This indicates the respiratory phase, or the respiratory state. Indicates the first K-space data collected from each respiratory state. Indicates the first K-space sampling mask corresponding to each breathing state Indicates Fourier transform, Indicates coil sensitivity. The first part represents the estimation of the motion solution network. The respiratory motion field corresponding to each breathing state This represents the motion-compensated reconstructed image estimated by the image reconstruction network. Represents the finite difference operator. , , This represents the weight coefficients corresponding to different items.

9. The method according to claim 8, characterized in that, The loss function for optimizing the motion solving network and image reconstruction network for the second subproblem. for: in, This represents the loss term for image edge sharpness. This represents the corresponding weighting coefficient.

10. A rapid three-dimensional cardiac magnetic resonance imaging system with free breathing capability, characterized in that, Used to perform the method described in any one of claims 1-9.

Images