Accelerated Imaging Methods#
MRI scans can be accelerated by applying methods that reduce the number of k-space samples required, also known as undersampling. These require a corresponding advanced image reconstruction method to allow for the reduced sampling. The techniques covered here include parallel imaging, compressed sensing, and deep learning methods.
Learning Goals#
Describe how images are formed
Describe how accelerated imaging methods work
Understand the constraints and tradeoffs in MRI
Identify accelerated imaging methods
Understand when accelerated imaging methods can be applied
Manipulate MRI sequence parameters to improve performance
Define which parameters are modified when accelerated imaging methods are applied
Manipulate and analyze MRI data
Reconstruct an image from undersampled raw data
General Formulation of MRI Reconstruction#
For advanced image reconstruction methods, it is helpful to reformulate the reconstruction problem as a linear system and using standard mathematical notation describing linear systems as:
where \(y\) is the acquired data, \(\mathbf{E}\) is the encoding matrix, \(x\) is the spatial distribution of the transverse magnetization (e.g. image), and \(n\) is noise. In this formulation, the image is vectorized. For example, 2D FT sampled data and the corresponding 2D image would be converted to:
where \(s_m(t_n)\) is the data from the \(m\)th TR at the \(n\)th sample.
where \(x_i\) and \(y_i\) are equally spaced locations in image space.
The encoding matrix, \(\mathbf{E}\), at a minimum includes a discrete Fourier Transform matrix, \(\mathbf{F}\), representing our data is the Fourier Transform of the transverse magnetization, evaluated at k-space locations. For this matrix, the entries are defined by k-space location being sampled, \(\vec{k}_i\) and the spatial locations, \(\vec{r}_j\):
In this case, the measurement is the discrete Fourier Transform of the image
and image reconstruction is performed by inverse discrete Fourier Transform, which for fully sampled 2D FT imaging is well defined by matrix inversion:
Partial Fourier Imaging#
Partial Fourier imaging, also known as fractional NEX or partial k-space, utilizes the conjugate symmetry property of the Fourier Transform. This property states that the Fourier Transform of a real-valued function has conjugate symmetry. If our image, \(m(\vec{r})\), is real-valued, then the k-space data satisfies
In this case, only half of k-space is required, and the other half can be filled in by this property.
In practice, there is typically low spatial frequency phase components in the image so the initial assumption of a real-valued image is violated. However, the data acquisition can still be accelerated by acquiring slightly more than half of k-space, and the fully sampled center of k-space can provide sufficient information about the low spatial frequency phase to still allow for accelerationg based on this conjugate symmetry property.
Parallel Imaging#
For parallel imaging (PI), we need to consider the coil sensitivity profiles, \(\mathbf{C}_q\), for each RF coil into encoding matrix along with a Fourier Transform encoding matrix, \(\mathbf{F}\), as well as a k-space sub-sampling operator, \(\mathbf{S}\), for the measurements from each RF coil, \(\mathbf{y}_q\):
The coil sensitivity profiles \(\mathbf{C}_q\) is a diagonal matrix with entries corresponding to the coil sensitivity profile at each location in \(\mathbf{x}\).
The k-space sub-sampling operator \(\mathbf{S}\), is a diagnoal matrix with entries of 1 or 0, describing whether an expect grid location was sampled. It is convenient to use when describing image reconstruction as filling in points in k-space.
Here we can return to our original formulation by creating augmented matrices, concatenating along the coil elements dimension for example as:
resulting in:
Image-space Methods#
Image-space parallel imaging methods (e.g. SENSE) can be formulated as the following optimization problem
This can be solved directly by a pseudo-inverse, but requires measurements of the coil sensitivity profiles to calculate \(\mathbf{E}\):
The k-space sampling patterns used for these methods typically use regular undersampling, meaning there is a consistent pattern of acquired and skipped k-space lines.
Auto-calibrated K-space Methods#
K-space parallel imaging methods (e.g. GRAPPA) utilize a calibration kernel, computed from the data itself and captured in the matrix \(\mathbf{G}\). These can also be generally formulated as the following optimization problem
Where \(\lambda\) is a tuning or regularization parameter and \(\mathbf{I}\) is the identity matrix. This simultaneous enforces data consistency (first term) as well as self-consistency of the calibration (sescond term).
The calibration kernel, \(\mathbf{G}\), can either be directly measured from the data as in GRAPPA, or iteratively estimated while solving the optimization, as in SPIRiT.
The k-space sampling patterns used for these methods typically use regular undersampling, meaning there is a consistent pattern of acquired and skipped k-space lines, combined with full sampling in the center of k-space.
SNR in Parallel Imaging#
There is an SNR penalty when using these methods that varies in severity depending on the conditioning of the undersampled reconstruction matrix. Typically, more RF coil elements in the direction of the undersampling leads to a more well-conditioned reconstruction and lower SNR penalty. Conversely, fewer RF coil elements in the direction of the undersampling leads to a more ill-conditioned reconstruction matrix and higher SNR penatly. Also, regions will little difference betweeen RF coils (typically in the center of the body), also tend to have larger SNR penalties.
This is characterized by the “g-factor”, where \(g(\vec{r}) \geq 1\) characterizes a spatially-varying SNR loss that is dependent on the coil loading and geometry, k-space sampling pattern, and parallel imaging reconstruction method. When appyling an acceleration factor of \(R\), the total readout time reduces the SNR by \(\sqrt{R}\) as well, leading to the SNR relationship:
Artifacts with Parallel Imaging#
Parallel Imaging artifacts typically appear as aliasing. The undersampled data will have aliasing when not using a parallel imaging reconstruction method, and sometimes the aliasing is not completely removed. There are also aliasing like artifacts when there is a mismatch between the data and measured coil sensitivity profiles.
Compressed Sensing Methods#
Compressed Sensing (CS) theory says that an image that is compressible in some domain can then be reconstructed from a subset of data samples. In other words, we can further accelerate our data acquisition.
Compressed Sensing is formulated as the following optimization problem, specifically using the \(\ell_1\)-norm (\(\|x\|_1 = \sum_{i=1}^n |x_i|\)) that promotes sparsity in the solution:
which includes a data consistency term where the data multiplied by the encoding matrix must match the reconstructed image, and a regularization term that enforces that the image is sparse using the \(\ell_1\) norm in some domain through the sparsifying transform, \(\mathbf{W}\). There is a regularization factor, \(\lambda_{CS}\), that must be chosen to balance the data consistency and sparsity terms.
These methods specifically require k-space sampling patterns with pseudo-random undersampling, which is illustratic by the diagram below of the undersampled reconstruction procedure:
Caption (Figure from https://doi.org/10.1109/MSP.2007.914728 ): Heuristic procedure for reconstruction from undersampled data. A sparse signal (a) is 8-fold undersampled in its 1-D k -space domain (b). Equispaced undersampling results in signal aliasing (d) preventing recovery. Pseudo-random undersampling results in incoherent interference (c). Some strong signal components stick above the interference level, are detected and recovered by thresholding (e) and (f). The interference of these components is computed (g) and subtracted (h), thus lowering the total interference level and enabling recovery of weaker components.
Popular sparsifying transforms include total variation (TV), total generalized variation (TGV), and wavelets, where MRI data is typically more sparse. The pseudo-random undersampling typically uses a variable density that preferentially has more samples near the center of k-space.
This regularization term, \(\| \mathbf{Wx} \|_1\) can also be combined with iterative parallel imaging reconstruction methods as well, for example in the \(\ell_1\)-ESPIRiT method.
SNR in Compressed Sensing#
It is difficult to define SNR when using compressed sensing reconstruction methods, as they inherently perform some denoising when constraining the reconstruction based on sparsity. The apparent SNR can also vary significantly depending on the choice of the regularization factor, \(\lambda_{CS}\).
Artifacts with Compressed Sensing#
The most common artifacts from compressed sensing come from overfitting to the sparsity penalty. This result in either an over-smoothed, unnatural appearance or, in the case of a wavelet sparsifying transform, blocking artifacts.
Deep Learning Reconstructions#
Deep learning (DL) can be used to reconstructed undersampled k-space into image data from a subset of data samples as well. Conceptually, these methods can be trained to learn how to incorporate coil sensitivity information like parallel imaging and typical image sparsity patterns like compressed sensing. Since they learn from real-world data, they can learn information that is the most relevant to MRI data and thus have been shown to support higher acceleration factors. For example, they can learn what patterns within the image and/or k-space domain and most relevant and most common in MRI.
These methods can often be considered general solutions to a regularized least-squares objective
where \(\mathcal{R}(\cdot)\) can be a parallel imaging or compressed sensing regularizer as above, but can also be implicitly implemented via machine learning techniques. This means that the machine learning reconstruction is attempting to constrain or regularize the reconstruction.
The most straightforward approach using deep learning is to train a CNN to invert the encoding operation in image space, e.g. after applying the inverse Fourier Transform to the data:
The learned deeep neural network, \(g_{NN}()\), includes many weights that must be trained prior to being able to perform the reconstruction.
However, so-called un-rolled networks have become more popular for MRI reconstruction, as they can better incorporate parallel imaging information, data consistency, and also are often more compact. They are designed to mimic iterations of classical algorithms to solve optimization problems, and each iteration has been “un-rolled” into a series of cascaded neural networks. Some popular methods include MoDL and Variational Networks.
Variational Networks (VarNets) for undersampled MRI reconstruction are formulated as the above regularlized least-squares objective, but models the architecture after a single gradient update in iterative gradient descent (used to solve optimization problems) and uses a CNN (e.g. UNet) within the iterations
Caption (from https://arxiv.org/abs/2004.06688): The E2E-VN model takes under-sampled k-space as input and applies several cascades, followed by IFT and RSS. Bottom: The DC module brings intermediate k-space closer to measured values, the Refinement module maps multi-coil k-space to one image, applies a U-Net and maps back to k-space, and the SME estimates sensitivity maps used in the refinement step.
The Model-based Deep Learning (MoDL) method uses the formulation where the regularizer is a CNN that estimates the noise and aliasing patterns
and the estimate depends on a set of learned parameters, \(w\), in the neural network, \(\mathcal{N}_\mathbf{w}(\cdot)\). It’s structure is overall similar to VarNets, as both are unrolled architectures, but with slightly different formulations.
Deep Learning Requirements#
Data
Generally a lot of data (\(10^4\) in many radiology applications, but also depends a lot on the question being answered and data quality)
Fully sampled data typically required for supervised learning of undersampled MRI reconstruction
Data augmentation – realistic modifications of existing data to increase training data size
For reconstruction, fastMRI provides thousands of high quality training data: http://fastmri.med.nyu.edu, https://doi.org/10.48550/arXiv.1811.08839
Architecture
U-net, ResNet popular for working in image-domain
Un-rolled networks popular for image reconstruction.
Generative models possible
Transformers
Computation
GPUs are essential for training of deep networks
Training
Send input data forward through network, compare to the expected output (e.g. fully sampled image or data) update the network weights through back-propagation
Choose loss function, such as l1 or l2
Choose batch sizes, Learning rate, etc
SNR in DL Reconstructions#
It is difficult to define SNR when using machine learning reconstruction methods, as they inherently perform some denoising when constraining the reconstruction.
Artifacts with DL Reconstructions#
A challenge with machine learning reconstructions is that artifacts can be difficult to identify. They can result in an over-smoothed appearance, but often a well-trained network will retain realistic texture and noise appearance features. These methods can also overfit to the learned anatomy and data, meaning features might be erased or hallucinated, although using physics-based techniques provide some assurances of data-consistency.
Reference#
For a recent, comprehensive reference on thes methods I reccommend:
Hammernik K, Küstner T, Yaman B, Huang Z, Rueckert D, Knoll F, Akçakaya M. Physics-Driven Deep Learning for Computational Magnetic Resonance Imaging: Combining physics and machine learning for improved medical imaging. IEEE Signal Process Mag. 2023 Jan;40(1):98-114. doi: 10.1109/msp.2022.3215288. Epub 2023 Jan 2. PMID: 37304755; PMCID: PMC10249732.