Method and apparatus for generating a T1/T2 map

Abstract

A method and apparatus for generating a T1 or T2 map for a three-dimensional (3D) image volume of a subject. The method includes acquiring first, second, and third 3D images of the image volume of the subject. Signal evolutions of voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images. A simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations is obtained. The T1 or T2 map is generated by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

Claims

1. A method of generating a longitudinal magnetic relaxation time constant (T1) or transverse magnetic relaxation time constant (T2) map for a three-dimensional (3D) image volume of a subject, the method comprising: applying, by a gradient system, a magnetic field gradient; applying, by an excitation system, an excitation pulse to a subject and to receive signals from the subject; communicating with and controlling, by a computing system, the excitation system and the gradient system, to receive the signals from the excitation system; and executing, by the computing system, program code to control the gradient system and the excitation system to: acquire first, second, and third 3D images of the image volume of the subject at different times using a 3D magnetic resonance (MR) acquisition sequence with different MR acquisition parameters, wherein the 3D MR acquisition sequence used to acquire the first 3D image comprises a 3D data acquisition stage preceded by a preparation pulse, and the 3D MR acquisition sequence used to acquire the second 3D image comprises a 3D data acquisition stage which is not preceded by a preparation pulse; determine signal evolutions of voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images; obtain a simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations; and generate the T1 or T2 map by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

2. A method as claimed in claim 1, wherein the preparation pulse is an inversion recovery preparation pulse.

3. A method as claimed in claim 1, wherein the 3D MR acquisition sequence used to acquire the third 3D image comprises a 3D data acquisition stage preceded by a preparation pulse.

4. A method as claimed in claim 3, wherein the preparation pulse is a T2 preparation pulse.

5. A method as claimed in claim 1, wherein at least one of the 3D MR acquisition sequences uses variable-density Cartesian trajectories with spiral profile ordering.

6. A method as claimed in claim 1, wherein the 3D MR acquisition sequences are generated during successive consecutive cardiac cycles of the subject.

7. A method as claimed in claim 1, wherein at least one of the 3D MR acquisition sequences further comprises a motion estimation stage preceding a 3D data acquisition stage of the 3D MR acquisition sequence, wherein the motion estimation stage is used to acquire navigator data for correcting motion in the data acquired in the 3D data acquisition stage.

8. A method as claimed in claim 1, wherein obtaining at least one of the first, second and third 3D images comprises: obtaining a k-space dataset for each 3D MR acquisition sequence; performing a motion correction operation on each of the k-space datasets; and reconstructing the first, second, and third 3D images from the motion corrected k-space datasets.

9. A method as claimed in claim 1, wherein at least one of the first, second, and third 3D images are fat-suppressed 3D images in which signal contributions of fat tissue are suppressed.

10. A method as claimed in claim 9, wherein acquiring at least one of the first, second, and third 3D images comprises: acquiring, during each 3D MR acquisition sequence, a pair of echo 3D images; and using the pair of echo 3D images as an input to a water-fat separation algorithm which generates, as an output, the fat-suppressed 3D image.

11. A method as claimed in claim 1, wherein the simulation dictionary is obtained using extended phase graph (EPG) simulations.

12. A non-transitory computer readable medium having instructions recorded thereon, which, when executed by a computer, cause the computer to perform the method of claim 1.

13. A magnetic resonance (MR) apparatus comprising: a gradient system configured to apply a magnetic field gradient; an excitation system configured to apply an excitation pulse to a subject and to receive signals from the subject; and a computing system configured to communicate with and control the excitation system and the gradient system, to receive the signals from the excitation system, and execute program code to control the gradient system and the excitation system to: acquire first, second, and third three-dimensional (3D) images of the image volume of the subject at different times using a 3D MR acquisition sequence with different MR acquisition parameters, wherein the 3D MR acquisition sequence used to acquire the first 3D image comprises a 3D data acquisition stage preceded by a preparation pulse, and the 3D MR acquisition sequence used to acquire the second 3D image comprises a 3D data acquisition stage which is not preceded by a preparation pulse; determine signal evolutions of voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images; obtain a simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations; and generate a T1 or T2 map by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

14. A magnetic resonance (MR) apparatus comprising: a gradient system configured to apply a magnetic field gradient; an excitation system configured to apply an excitation pulse to a subject and to receive signals from the subject; and a computing system configured to communicate with and control the excitation system and the gradient system, to receive the signals from the excitation system, and execute program code to control the gradient system and the excitation system to: acquire first, second, and third three-dimensional (3D) images of the image volume of the subject at different times using a 3D MR acquisition sequence with different MR acquisition parameters, wherein at least one of the 3D MR acquisition sequences uses variable-density Cartesian trajectories with spiral profile ordering; determine signal evolutions of voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images; obtain a simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations; and generate a T1 or T2 map by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

15. A magnetic resonance (MR) apparatus comprising: a gradient system configured to apply a magnetic field gradient; an excitation system configured to apply an excitation pulse to a subject and to receive signals from the subject; and a computing system configured to communicate with and control the excitation system and the gradient system, to receive the signals from the excitation system, and execute program code to control the gradient system and the excitation system to: acquire first, second, and third three-dimensional (3D) images of the image volume of the subject at different times using a 3D MR acquisition sequence with different MR acquisition parameters, wherein at least one of the 3D MR acquisition sequences comprises a motion estimation stage preceding a 3D data acquisition stage of the 3D MR acquisition sequence, wherein the motion estimation stage is used to acquire navigator data for correcting motion in the data acquired in the 3D data acquisition stage; determine signal evolutions of voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images; obtain a simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations; and generate a T1 or T2 map by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) Examples of the present disclosure will now be described with reference to the accompanying drawings, in which:

(2) FIG. 1 shows a schematic of an example 3D MR acquisition sequences in accordance with the present disclosure;

(3) FIG. 2 is shows a detailed view of a first of the three 3D MR acquisition sequences shown in FIG. 1;

(4) FIG. 3 is shows a detailed view of a second of the three 3D MR acquisition sequences shown in FIG. 1;

(5) FIG. 4 is shows a detailed view of a third of the three 3D MR acquisition sequences shown in FIG. 1;

(6) FIG. 5 shows a schematic outlining an example of how the echo images obtained during the three MR acquisition sequences of FIG. 1 are used to generate water images for the three MR acquisition sequences;

(7) FIG. 6 shows a schematic outlining how the obtained water images are used to determine the signal evolutions for the voxels;

(8) FIG. 7 shows a schematic outlining how the determining signal evolutions are compared with a simulation dictionary so as to generate the T1 and T2 maps;

(9) FIG. 8 is a process diagram of an example method according to aspects of the present disclosure;

(10) FIG. 9 is a process diagram of a more detailed example method according to aspects of the present disclosure; and

(11) FIG. 10 is a simplified schematic diagram of an example MR apparatus according to aspects of the present disclosure.

DETAILED DESCRIPTION

(12) Referring to FIG. 1, there is shown a schematic of an MR pulse sequence chart according to aspects of the present disclosure. The MR pulse sequence chart shows example 3D MR acquisition sequences used to acquire first, second, and third 3D images of the image volume of the subject, wherein each of the first, second, and third 3D images are acquired at a different time using a 3D MR acquisition sequence with different MR acquisition parameters. Three 3D MR acquisition sequences are shown occurring during three consecutive cardiac cycles of the subject.

(13) The first cardiac cycle is shown as beginning with ECG pulse 101. The first 3D MR acquisition sequence is ECG-gated, meaning that the generation of the first 3D MR acquisition sequence is triggered by the ECG pulse 101. The first 3D MR acquisition sequence comprises an inversion recovery (IR) pulse 103. A motion estimation sequence 105 is performed following the IR pulse 103 so as to acquire an image navigator (iNav) for the subject. The iNav in this example is a low-resolution 2D iNav. A data acquisition sequence 107 follows the motion estimation sequence 105. The data acquisition sequence 107 is performed an inversion time TI 108 after the generation of the inversion recovery pulse 105. The data acquisition sequence 107 obtains k-space data for the whole-heart region of the subject.

(14) The second cardiac cycle begins with the ECG pulse 109. The second 3D MR acquisition sequence is ECG-gated such that the second 3D MR acquisition sequence is triggered by the ECG pulse 109. The second 3D MR acquisition sequence comprises a motion estimation sequence 111. The motion estimation sequence 111 is performed so as to acquire an iNav for the subject, and in this particular example a low-resolution 2D iNav. A data acquisition sequence 113 follows the motion estimation sequence 111. Unlike the first 3D MR acquisition sequence, the second 3D MR acquisition sequence does not comprise a preparation pulse before the data acquisition sequence 113. This means that the second 3D MR acquisition sequence does not comprise an IR pulse or another preparation pulse such as a T2 preparation pulse (T2 prep). The data acquisition sequence 113 obtains k-space data for the whole-heart region of the subject.

(15) The third cardiac cycle begins with the ECG pulse 115. The third 3D MR acquisition sequence is ECG-gated such that the third 3D MR acquisition sequence is triggered by the ECG pulse 115. The third 3D MR acquisition sequence comprises a T2 preparation pulse (T2prep) 117 followed by a motion estimation sequence 119. The motion estimation sequence 119 is performed so as to acquire an iNav for the subject, and in this particular example a low-resolution 2D iNav. A data acquisition sequence 121 follows the motion estimation sequence 119. The motion estimation sequence 119 is generated shortly after the T2prep pulse 117. The data acquisition sequence 121 is performed shortly after the motion estimation sequence 119. There is no IR pulse following the T2prep pulse 117 in the third imaging sequence. The data acquisition sequence 121 obtains k-space data for the whole-heart region of the subject.

(16) Referring to FIG. 2, there is shown a detailed overview of an example first 3D MR acquisition sequence according to aspects of the present disclosure and as shown in FIG. 1. The data acquisition sequence 107 is a 3D gradient echo acquisition sequence that is under-sampled (i.e. accelerated). The acquisition sequence 107 follows a 3D Cartesian trajectory 123 with spiral profile ordering (CASPR). This means that the phase-encoding plane is sampled along approximate spiral interleaves on the Cartesian grid (Kz-Ky). The angular step between two consecutive successive spiral interleaves in the phase-encoding plane is, in this example, given by the golden ratio. The 3D acquisition scheme 123 shown in FIG. 2 also uses a variable-density approach. This means that the sampling density changes along the spiral trajectory, and in particular the centre-region of k-space in the phase-encoding plane is sampled more densely than the outer region of k-space in the phase-encoding plane. This data acquisition sequence 107 is known as 3D Cartesian variable-density (VD-CASPR).

(17) Two echo readouts are performed during the data acquisition sequence 107 so as to obtain two image datasets acquired with different echo times. In other words, a dual-echo Dixon imaging procedure is performed so as to obtain a two-point Dixon read-out. A first of the image datasets may be acquired with the fat and water signals in phase at the centre of the echo. This is known as an in-phase image dataset. The second of the image datasets may be acquired with the echo time (TE) adjusted by a few milliseconds so that the fat and water signals are out-of-phase. This is known as an out-of-phase image dataset. The first and second image datasets which are k-space data are reconstructed to generate 3D MR image volumes 125, 127. The 3D MR image volumes 125, 127 are bright blood image volumes 125, 127. The reconstruction of the MR image volumes 125, 127 from the image datasets is discussed in more detail below.

(18) Referring to FIG. 3, there is shown a detailed overview of an example second 3D MR acquisition sequence according to aspects of the present disclosure and as shown in FIG. 1. The data acquisition sequence 113 is a 3D undersampled gradient echo acquisition sequence. Similar to the data acquisition sequence 107 of FIG. 2, the data acquisitions sequence follows a 3D Cartesian trajectory 129 with spiral profile ordering and variable density sampling (VD-CASPR). Again, two Dixon image readouts are performed during the data acquisition sequence 113 so as to obtain two image datasets acquired with different echo times. The first and second image datasets which are k-space data are reconstructed to generate 3D MR image volumes 131, 133. The 3D MR image volumes 131, 133 are bright blood image volumes 131, 133. The reconstruction of the MR image volumes 131, 133 from the image datasets is discussed in more detail below.

(19) Referring to FIG. 4, there is shown a detailed overview of an example third 3D MR acquisition sequence according to aspects of the present disclosure and as shown in FIG. 1. The data acquisition sequence 121 is a 3D undersampled gradient echo acquisition sequence. Similar to the data acquisition sequence 107 of FIG. 2, the data acquisitions sequence follows a 3D Cartesian trajectory 135 with spiral profile ordering and variable density sampling (VD-CASPR). Again, two Dixon image readouts are performed during the data acquisition sequence 121 so as to obtain two image datasets acquired with different echo times. The first and second image datasets which are k-space data are reconstructed to generate 3D MR image volumes 137, 139. The 3D MR image volumes 137, 139 are bright blood image volumes 131, 133. The reconstruction of the MR image volumes 137, 139 from the image datasets is discussed in more detail below.

(20) In summary, the MR pulse sequence chart shown in FIGS. 1 to 4 results in the generation of three ECG-triggered interleaved 3D bright blood volume datasets. These three interleaved 3D bright blood volume datasets are acquired with data acquisitions sequences 107, 113, 121 that use gradient echo two points Dixon readouts and a variable density Cartesian trajectory with spiral profile order (VD-CASPR). The data acquisition sequences 107, 113, 121 are undersampled (i.e. accelerated). The first dataset is acquired with an IR preparation pulse 103. The second dataset is acquired without preparation pulses, and the third dataset is acquired with a T2prep pulse 117.

(21) Image Reconstruction

(22) The interleaved 3D volume data sets acquired using the first through third 3D MR acquisition sequences are reconstructed so as to generate the MR image volumes 125, 127, 131, 133, 137, and 139 shown in FIGS. 2 through 4.

(23) Prior to the image reconstruction, the acquired 2D image navigators are used to estimate the superior-inferior (SI) and left-right (LR) translational respirator motion of the subject. Translational motion correction is then performed on each echo and each volume is reconstructed independently. This means that data does not need to be rejected due to subject motion, and a 100% scan efficiency can be obtained using a free-breathing sequence.

(24) The motion corrected k-space data for each echo readout is then reconstructed so as to generate the corresponding 3D image 125, 127, 131, 133, 137, and 139.

(25) In general terms, the objective of a reconstruction process is to find an MR image X for which the following inverse problem is minimised:
{tilde over (X)}=argmin1/2∥EX-Y∥2/2

(26) Here, X denotes the reconstructed MR image X ∈.sup.M.sup.z.sup.×M.sup.y.sup.×M.sup.z.

(27) Here, argmin.sub.xf(X) is the argument of the minimum operation and finds the value of the variable X at which the function f(X) is minimised.

(28) Here, Y refers to the k-space data, which can be denoted as Y ∈.sup.M.sup.z.sup.×M.sup.y.sup.×M.sup.z.sup.×N.sup.c, M.sub.x denotes the number of voxels in the x spatial direction. M.sub.y denotes the number of voxels in the y spatial direction. M.sub.z denotes the number of voxels in the z spatial direction. N.sub.c denotes the number of channels, i.e. the number of receiver coils of the MR apparatus. In practice, the k-space data is undersampled and thus may be defined as Y ∈.sup.z×L×N.sup.c, where Z«M.sub.m×M.sub.y×M.sub.x.

(29) Here, E is an encoding operator that transforms the MR images, X, into k-space. E may be defined as being equal to AFS.sub.c. A is the undersampling operator that acquires k-space data for the image. F is the Fourier transform operator. S.sub.c is the complex receive-coil sensitivity maps for the N.sub.c channels and may be denoted as {S.sub.c}.sub.c=1.sup.N.sup.c ∈.sup.M.sup.x.sup.×M.sup.y.sup.×M.sup.z.sup.×N.sup.c. These coil-sensitivity maps S.sub.c may be determined from MR acquisitions of a fully-sampled k-space centre using a known eigenvalue-based approach commonly referred to as ESPIRIT.

(30) Here, ∥.∥ 2/2 denotes the I2-norm.

(31) The above equation is an inverse problem because it is attempting to determine the source (the MRI images, X, of the image region) from a limited set of observations (the undersampled k-space data, Y of the image region). This inverse problem is generally ill-posed because there may either be several possible solutions (MRI images) to the problem and no way to determine which solution is the correct one. Further, a solution to the inverse problem may not exist at all. Therefore, generally image reconstruction techniques use prior assumptions of the unknown solution, X, when solving the inverse problem.

(32) One example image reconstruction technique that may be used with the present disclosure is parallel imaging. Parallel imaging involves using the known placement and sensitivities of the receiver coils of the MR apparatus to assist in the spatial localization of the MR signal. Parallel imaging techniques process the received signals in parallel along separate channels, one for each receiver coil. This contrasts with non-parallel techniques which first generate a combined signal from the receiver coils before performing the image reconstruction. The parallel imaging techniques consider the relative intensities sensed by the receiver coils and use this information to estimate the location of the tissue corresponding to the sensed signal. The prior information about the receiver coils is therefore exploited for the image reconstruction, which allows for an under-sampled k-space scheme to be used (e.g. by reducing the number of phase encoding steps during image acquisition).

(33) Another example image reconstruction technique that may be used with the present disclosure is compressed sensing. Compressed sensing techniques can involve reconstructing the undersampled k-space data by solving an optimisation problem of the form:

(34) $\hat{x}=\arg \underset{x}{\min }\left\{\frac{1}{2}{.Math.\mathrm{Ex}-k.Math.}_2^2+{\mathrm{\gamma \Psi }}_t\left(x\right)\right\}$

(35) x is the reconstructed image series, k is the undersampled k-space data. The encoding operator E=AFS incorporates the sampling matrix A for the k-space data, Fourier transform Fund coil sensitivities S.

(36) ∥.∥ 2/2 denotes the I.sub.2-norm defined as ∥u∥ 2/2=√(Σ.sub.i|u.sub.i|.sup.2), for a given vector u defined by vector components u.sub.i.

(37) Ψ.sub.c is a 1D temporal total variation (TV) function or a wavelet-based regularization. γ is a regularization parameter.

(38) The component ½∥Ex-k∥ 2/2 of the optimisation problem is known as a data consistency component of the compressed sensing procedure.

(39) The component Ψ.sub.t(x) of the optimisation problem is known as a transform sparsity component of the compressed sensing procedure.

(40) The parameter γ is a regularization parameter that determines the trade-off between the data consistency and the sparsity.

(41) Another example image reconstruction technique that may be used with the present disclosure exploits the anatomical similarity within patches of an image and is known as ‘3D-Patch-based low-rank ReconstructiOn with Self-similariTy learning (3D-PROST)’. 3D-PROST integrates anatomical structure information from 3D patch neighbourhoods through sparse representation, exploiting the redundancy of 3D patches in the acquired data itself. The optimization problem iterates between a data consistency step, which reconstructs a high-resolution isotropic MR volume x, and a low-complexity 3D patch-based denoising step, which provides a reconstructed volume as prior information for the next step. A detailed explanation of the 3D-PROST procedure is disclosed in A. Bustin et al., “Five-minute whole-heart coronary MRA with sub-millimeter isotropic resolution, 100% respiratory scan efficiency, and 3D-PROST reconstruction,” Magn. Reson. Med., 2019; 81:102-115, the disclosures of which are hereby incorporated by reference.

(42) Fat Suppression

(43) Referring to FIG. 5, there is shown an overview of a process by which the MR image volumes 125, 127, 131, 133, 137, and 139 for each dataset are input to a water fat separation algorithm indicated by the process blocks 139, 141, 143 so as to generate a water image 145, 147, 149 for each dataset.

(44) In each data set, one of the MR image volumes (e.g. the MR image volume 125 in the first dataset) was acquired when the fat and water signals were in phase at the centre of the echo. The MR image volume 125 can be considered as an in-phase image 125. In simple terms, an in-phase (IP) image can be defined as IP=W+F, where W denotes the signal contributions from water and F denotes the signal contributions from fat. The other of the MR image volumes (e.g. the MR image volume 127 in the first dataset) was acquired when the fat and water signals are out-of-phase. The MR image volume 127 can be considered as an out-of-phase image 127. In simple terms, an out-of-phase (OP) image can be defined as OP=−F, where W denotes the signal contributions from water and F denotes the signal contributions from fat.

(45) In a simple example, a water image (i.e. an image with the signal contributions from fat supressed) can be obtained by performing a pixel-wise addition of the in-phase and out-of-phase images and averaging each of the resultant summed pixel values. Specifically:
½┌IP+OP┐=½┌(W+F)+(W−F)┐=½2W┐=W

(46) In reality, a simple addition of pixel values does not lead to a reliable estimate of the water image due to effects such as B0 inhomogeneity.

(47) Some existing two-point Dixon techniques explicitly estimate the B0 map by performing phase unwrapping. A particular approach is disclosed in: J. Liu, D. C. Peters, and M. Drangova, “Method of B0 mapping with magnitude-based correction for bipolar two-point Dixon cardiac MRI,” Magn. Reson. Med., vol. 78, no. 5, pp. 1862-1869, December 2016., the disclosure of which is hereby incorporated by reference. This approach seeks to improve the reliability of B0 mapping in the water-fat separation algorithm by performing unwrapping error correction. In summary, this approach determines a magnitude-based fat/water mask from the in-phase, out-of-phase image pair (the pair of echo images) and uses this mask as a reference to identify pixels being mismatched by a phase-based mask. The phase-based mask was derived from the B0-corrected phase term of the Hermitian product between in-phase, out-of-phase image pair. Then, the phase error is corrected on a region-by-region basis. Of course, the skilled person will appreciate that other water-fat separation algorithms may be used to generate the water image. In addition to the water image, the water-fat separation algorithm also generates a complementary fat image (water-supressed image) which can also be used to provide useful tissue information regarding the tissue characteristics.

(48) Referring to FIG. 6, there is shown a schematic diagram of another stage of the process in accordance with the present disclosure. In this stage of the process the change in voxel intensity value for each corresponding voxel in the three image datasets is considered so as to obtain the signal evolution.

(49) It will be appreciated that a voxel in the first 3D water image 145 at a voxel location (x1, y1, z1) may have a different signal intensity value to a corresponding voxel located at position (x1, y1, z1) in the second 3D water image 147 due to the difference in the preparation stage of the first and second 3D MR acquisition sequences. This change in signal intensity value from the first 3D water image 145 to the second 3D water image 147 is known as the signal evolution of the voxel. The signal evolutions are between the first 3D water image 145 and the second 3D water image 147 and between the second 3D water image 147 and the third 3D water image 149 are obtained for each of the voxels. It will be appreciated that the voxel values may be normalized in time (i.e. across the three water images 145, 147, 149) before the signal evolutions are determined. That is the signal evolutions may be normalized to the unit norm.

(50) FIG. 6 shows a plot 151 of the acquired signal evolutions. The y-axis denotes the signal intensity values and the x-axis denotes time. The plot 151 shows how the signal intensity values change from the first image volume 145 to the second image volume 147 and from the second image volume 147 to the third image volume 149.

(51) Referring to FIG. 7, there is shown a schematic diagram of another stage of the process in accordance with the present disclosure. In this stage of the process the determined signal evolutions 151 of the voxels are compared to a simulation dictionary 153 representing the signal evolutions for a number of different tissue parameter combinations. Like the determined signal evolutions, the signal evolutions in the simulation dictionary may be normalized to the unit norm.

(52) The simulation dictionary may be obtained using a computational algorithm that simulates the spin behaviour for a number of different tissue parameters during 3D MR acquisition sequences using the same MR acquisition parameters as used to acquire the 3D images in accordance with the present disclosure. That is the simulations use the same acquisition parameters as the actual parameters of the 3D MR acquisition sequences.

(53) In one example, Bloch equations may be used to generate the simulation dictionary. Bloch equations may be time consuming, so a preferred method in accordance with the present disclosure uses extended phase graph (EPG) simulations to as to obtain the simulation dictionary. EPG simulations work by representing a spin system affected by the sequence as a discrete set of phase states. An existing approach for using EPG simulations is provided in: M. Weigel, “Extended phase graphs: Dephasing, RF pulses, and echoes—pure and simple,” J. Magn. Reson. Imaging, vol. 41, no. 2, pp. 266-295, February 2015, the discloses of which are hereby incorporated by reference.

(54) Further, the T1 map 155 and the T2 map 157 for the subject are generated by comparing the determined signal evolutions 151 to entries in the dictionary 153 and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution. In other words, a pattern recognition operation is performed to find the dictionary entry that is the best representation of the acquired signal evolution for each voxel.

(55) In one example, the pattern recognition operation takes the inner product between the determined signal evolution and each of the entries in the dictionary. The entry that returns the highest inner product value can then be considered as the dictionary entry that best represents the tissue properties, and the respective T1 and T2 values are then assigned to the voxel.

(56) The matching operation may be accelerated by compressing the dictionary either in the time dimension or in the parameter combinations dimension, thus reducing the total number of comparisons that need to be performed. For example, the singular value decomposition (SVD) can be applied to compress the dictionary in the time dimension and thus reduce the matching time. Another approach to reduce the computational time by matching is to reduce the parameter combination dimension.

(57) For example, dictionary entries that have strong correlations may be grouped together and a new signal that best represents the group may be generated.

(58) Referring to FIG. 8, there is shown an example process diagram for a method according to the first aspect of the present disclosure.

(59) In step 201, the method comprises acquiring first, second, and third 3D images of the image volume of the subject, wherein each of the first, second, and third 3D images are acquired at a different time using a 3D MR acquisition sequence with different MR acquisition parameters.

(60) In step 202, the method comprises determining signal evolutions of the voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images.

(61) In step 203, the method comprises obtaining a simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations.

(62) In step 204, the method comprises generating the T1 and/or T2 map by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

(63) Referring to FIG. 9, there is shown a detailed overview of a method according to the first aspect of the disclosure which combines the individual steps shown in FIGS. 1 to 7. In this example implementation three ECG-triggered interleaved 3D bright-blood volumes are acquired with gradient echo two points Dixon read out and a variable density Cartesian trajectory with spiral profile ordering (VD-CASPR) 123, 129, 135 and an acceleration factor of greater than 1. In other words, the k-space is undersampled. The acceleration factor used in this particular example is 4.

(64) The first dataset is acquired using an 3D MR acquisition sequence that comprises an IR pulse 103 that precedes a motion estimation sequence 105 used to acquire a 2D image navigator and the 3D data acquisition sequence. The second dataset is acquired using a 3D MR acquisition sequence that comprises a motion estimation sequence 111 and a 3D data acquisition sequence 113 but does not comprise a preceding preparation pulse. The third dataset is acquired using a 3D MR acquisition sequence that comprises a T2 prep pulse 117, followed by a motion estimation sequence 119 and a 3D data acquisition sequence 121. The

(65) The 2D image navigators acquired prior to the 3D data acquisition sequences 107, 113, 121 are used to estimate superior-inferior (SI) and left-right (LR) translational respiratory motion. Translation correction is then performed on each echo and each volume is independently reconstructed, preferably by using a 3D patch-based undersampled reconstruction (3D-PROST). The resultant 3D echo images 125, 127, 131, 133, 137, 139 are used as inputs to a water fat separation algorithm so as to generate fat images and water images 145, 147, 149 for each acquired dataset. The 3D motion corrected water images 145, 147, 149 are normalized in time in order to obtain the signal evolution 151, through the three acquired volumes, of each voxel. Extended phase graph (EPG) simulations, matching the acquisition parameters, are carried out to generate a subject-specific simulation dictionary 153. Quantitative T1 and T2 maps 155, 157 are generated by matching each measured signal evolution 151 to a specific dictionary entry 153, corresponding to a unique T1/T2 pair.

(66) Referring to FIG. 10, there is shown an example MR apparatus 400 according the second aspect of the disclosure. The MR apparatus 400 comprises a computing system 401, a gradient system 403, and an excitation system 405.

(67) The gradient system 403 is configured to apply a magnetic field gradient. The gradient system 403 may be configured to apply magnetic field gradients along three spatial axes.

(68) The excitation system 405 may comprise a transmitter (not shown) and a receiver (not shown). The excitation system 405 can be an RF system with one or more RF coils (not shown). The excitation system 405 is configured to apply an excitation pulse to the subject and to receive signals from the subject.

(69) The MR apparatus 400 includes a magnet (not shown) for establishing a stationary magnetic field. The magnet can include a permanent magnet, a superconducting magnet or other type of magnet.

(70) The computing system 401 is in communication with the excitation system 305, and the gradient system 403 for controlling these components. The computing system 401 is configured to receive the signals from the excitation system 405.

(71) The computing system 401 is further configured to execute program code to control the gradient system 403 and the excitation system 405 to acquire first, second, and third 3D images of the image volume of the subject, wherein each of the first, second, and third 3D images are acquired at a different time using a 3D MR acquisition sequence with different MR acquisition parameters. It will be appreciated that the gradient system 403 and the excitation system 405 will acquire k-space datasets that may then be converted into image data by the computing system.

(72) The computing system 401 is further configured to determine signal evolutions of the voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images. The computing system 401 is further configured to obtain a simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations. The computing system 401 is further configured to generate the T1 and/or T2 map by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution.

(73) In summary, there is provided a method and apparatus for generating a T1 and/or T2, map for a three-dimensional, 3D, image volume of a subject. The method comprises acquiring first, second, and third 3D images of the image volume of the subject (S201). Signal evolutions of the voxels through the first to third 3D images by comparing voxel intensity levels of corresponding voxel locations in the first, second, and third 3D images (S202). A simulation dictionary representing the signal evolutions for a number of different tissue parameter combinations is obtained (S203). The T1 and/or T2 map is generated by comparing the determined signal evolutions to entries in the dictionary and by finding, for each of the determined signal evolutions, the entry in the dictionary that best matches the determined signal evolution (S204). The method and apparatus can be used to generate high resolution joint T1/T2 maps.

(74) Attention is directed to all papers and documents which are filed concurrently with or previous to this specification in connection with this application and which are open to public inspection with this specification, and the contents of all such papers and documents are incorporated herein by reference.

(75) All of the features disclosed in this specification (including any accompanying claims, abstract and drawings), and/or all of the steps of any method or process so disclosed, may be combined in any combination, except combinations where at least some of such features and/or steps are mutually exclusive.