|Home | About | Journals | Submit | Contact Us | Français|
Parallel imaging reconstruction has been successfully applied to magnetic resonance spectroscopic imaging (MRSI) to reduce scan times. For undersampled k-space data on a Cartesian grid, the reconstruction can be achieved in image domain using a sensitivity encoding (SENSE) algorithm for each spectral data point. Alternative methods for reconstruction with undersampled Cartesian k-space data are the SMASH and GRAPPA algorithms that do the reconstruction in the k-space domain. To reconstruct undersampled MRSI data with arbitrary k-space trajectories, image-domain-based iterative SENSE algorithm has been applied at the cost of long computing times. In this paper, a new k-space domain-based parallel spectroscopic imaging reconstruction with arbitrary k-space trajectories using k-space sparse matrices is applied to MRSI with spiral k-space trajectories. The algorithm achieves MRSI reconstruction with reduced memory requirements and computing times. The results are demonstrated in both phantom and in vivo studies. Spectroscopic images very similar to that reconstructed with fully sampled spiral k-space data are obtained at different reduction factors.
Magnetic resonance spectroscopic imaging has been increasingly used in diagnosis and treatment planning to provide information of metabolite distribution. The major limitation of magnetic resonance spectroscopic imaging (MRSI) is the long scan times required to acquire data in both spatial and spectral domains. In particular, when 3D MRSI scans are conducted, data in 4 dimensions, kx, ky, kz, and kf, have to be collected, resulting in impractical scan times with standard phase encoding (1,2).
MRSI scan times are reduced using fast imaging techniques. By using time-varying gradients, spectral and spatial k-space data can be more efficiently sampled during each repetition time. Among them, echo-planar and spiral k-space trajectories have been widely used (3,4).
Another approach, independent to MRSI with time-varying gradients, to reduce MRSI scan times is to apply parallel imaging techniques to reconstruct MRSI using partial k-space data with the knowledge of each coil’s sensitivity profile. For undersampled k-space data on a Cartesian grid, the sensitivity encoding (SENSE) algorithm (5) can be applied to each spectral point in image domain after performing the Fourier transform to obtain the unaliazed MRSI spectra. With k-space data on a Cartesian grid, MRSI reconstruction can also be achieved in the k-space domain using the simultaneous acquisition of spatial harmonics (SMASH) (6) or generalized autocalibrating partially parallel acquisitions (GRAPPA) (7) algorithms. To date, these parallel imaging methods have been successfully applied to MRSI with conventional phase encoding (8,9) and echo planar k-space trajectories (10,11).
To further speed up MRSI with non-Cartesian k-space sampling such as spiral k-space trajectories, parallel imaging reconstruction with arbitrary location k-space data can be applied (12–16). Similar to parallel imaging with Cartesian k-space data, the reconstruction with non-Cartesian k-space data can be achieved in both image and k-space domains. The typical image domain-based algorithm, iterative SENSE (12), solves the reconstruction problem by expressing k-space data as a linear combination of the spatially encoded functions and then solving it using the conjugate gradient (CG) method. The k-space domain-based algorithms, such as generalized GRAPPA and parallel magnetic resonance imaging with adaptive radius in k-space (PARS) algorithms (14,16), on the other hand, achieves the reconstruction by estimating the k-space data on a Cartesian grid using its neighboring data points on arbitrary k-space trajectories. Although effective, these reconstruction methods using undersampled arbitrary k-space data suffer from long computing times with large MRSI data sets.
In this work, we propose a rapid parallel imaging reconstruction method with arbitrary trajectories using k-space sparse matrices (KSPA) (17). The image reconstruction problem is formulated in the k-space as a linear algebra problem. Taking advantage of the compactness of the convolution kernel defined by the coil sensitivity, the reconstruction matrix can be approximated as a sparse matrix. The reconstruction using this algorithm is demonstrated using undersampled spiral k-space data for both phantom and in vivo studies with different reduction factors.
From the Fourier transform theorem, the k-space data at arbitrary locations acquired from each coil are the convolution of the nonsensitivity weighted k-space data with the Fourier transform of the coil’s sensitivity. Suppose data are acquired at nk sampling locations in k-space with nc receiving coils and the fully sampled Cartesian grid is N × N in dimension. The KSPA parallel reconstruction algorithm forms the reconstruction problem as a linear system:
where d is a column vector stacked with the k-space data acquired by all coils and m is a column vector with Cartesian k-space value to be estimated. G coefficient matrix composed of a coil sensitivity spectrum interpolated at an arbitrary k-space location. The problem can be solved by finding a pseudoinversion matrix G+ of the coefficient matrix G.
Due to the large size of G, the direct computation of G+ requires prohibitively large storage and a long computing time. However, because the coil sensitivity is a smooth function and contains only low spatial frequency components, the coefficient matrix G is a sparse matrix. Exploiting the sparsity of the coefficient matrix G, KSPA algorithm solves G+ by finding a sparse approximate inverse (GHG–1).
From the fact that most entries in GHG are zero and the assumption that a k-space sample on a Cartesian grid is only affected by its its neighboring k-space samples, its approximate inverse M+ can be obtained by solving equation:
where I is the identity matrix, row by row using the least-square method with significantly less complexities than a direct inversion (18). After M+ is obtained, the G+ and the k-space data on a fully sampled Cartesian grid can be readily calculated with matrix multiplication (17).
The KSPA parallel reconstruction algorithm’s ability of estimating k-space data points on a fully sampled Cartesian grid with undersampled arbitrary k-space trajectory data can be readily applied to MRSI reconstruction to reduce scan times. To encode chemical shift information, fast MRSI using time-varying gradients samples data points on repeated arbitrary k-space trajectories. Therefore, large number of fully sampled k-space data on repeated Cartesian grids need to be estimated with the same arbitrary k-space trajectories. The reconstruction time is thus substantially reduced by using one identical reconstruction matrix to estimate the fully sampled k-space data on repeated Cartesian grids obtained along the kf dimension.
Using multiple receiving coils for data acquisition, the KSPA reconstruction is applied to high-resolution MRSI with spiral k-space trajectories. The spiral gradient waveforms are generated in real time according to the prescription (19). Reconstruction with different reduction factors are performed using the KSPA algorithm with partial k-space data. For example, for a reduction factor of 2, k-space data on every the other spiral interleaf are used for the MRSI reconstruction and for a reduction factor of 3, k-space data on every the third spiral interleaf are used for the MRSI reconstruction.
The coil sensitivity maps required for the KSPA reconstruction are estimated with two steps. First, the raw sensitivity map of each coil is estimated by dividing the coil water image with the sum-of-squares water image. Second, a 2D thin-plate spline fitting is conducted to both the real and imaginary part of the sensitivity map independently to obtain a smoothed sensitivity map.
Using the KSPA algorithm, the reconstruction matrix is obtained with data on the first spiral k-space trajectories of each repetition time. K-space data on fully sampled Cartesian grids for each repeated spiral trajectory are estimated by applying the reconstruction matrix to each repeated spiral k-space trajectory data. Inverse Fast Fourier Transform (FFT) and 4-Hz Gaussian apodization are then performed to obtain the reconstructed spectrum. The spectrum of each voxel is phased with zero and first-order phasing to obtain the spectrum in the absorption mode.
The KSPA reconstructions with different reduction factors are compared to a reference standard reconstructed with fully sampled k-space data using a conventional gridding algorithm. The fully sampled k-space data for each receiving coil are gridded on a Cartesian grid to perform the inverse FFT before apodization with a 4-Hz Gaussian filter. After correcting for coil sensitivity, the individual images are summed to produce the combined spectrum.
To evaluate the noise enhancement of the KSPA recontruction algorithm, g-factor maps are calculated for different reduction factors for regions with metabolite signals. For voxels within the excited region, g-factors are determined as:
where SNAA,1 is the NAA peak reconstructed without reduction, σ1 is the standard deviation of the noise without reduction, SNAA,R is the NAA peak reconstructed with a reduction factor of R, and σR is the standard deviation of the noise with reduction factor of R. The standard deviation of noise is obtained from the signal free region, 5.5–8 ppm, in the spectrum.
All data were collected at 1.5T using a General Electric (G.E. Health Care, Milwaukee, WI, USA) scanner with an eight-channel phased array coil and a gradient system with maximum gradient of 4 Gauss/cm and maximum slew rate of 15 Gauss/cm/ms. A phantom study was performed using a standard GE spherical MRS phantom. A 2D single slice was localized with a PRESS module. The PRESS box was prescribed to cover the whole slice of the MRS sphere with the slice thickness of 2 cm, echo time (TE)/pulse repetition time(TR) = 144/1500 ms, 32 cm field of view (FOV), 64 × 64 matrix size, and 350-Hz spectral bandwidth. This resulted in a voxel size of 0.5 cc and the spiral trajectories for this prescription required 16 spatial interleaves. With 64 NEX, the total scan time was 26 min. Before the MRS measurement, high-order shimming was performed over the region of interest with a resulting full width at half maximum (FWHM) line width of 5 Hz.
The in vivo study was conducted on a healthy young male subject. Informed consent approved by the institutional review board was obtained from the subject before the study. Similar to the phantom study, a 2D single slice was localized with a PRESS module. However, the PRESS box was prescribed inside the brain across the ventricle region to avoid lipid signal contamination. The sequence was prescribed with slice thickness of 2 cm, TE/TR = 144/1500 ms, 32 cm FOV, 64 × 64 matrix size, 0.5 cc voxel, 16 spiral interleaves, 64 NEX, and 26-min acquisition. Once again, high-order shimming was performed before the MRS measurement over the region of interest, resulting in an FWHM line width of 9 Hz.
For the phantom study, reconstruction with reduction factors 1, 2, 3, and 4 were obtained over the whole FOV using both conventional gridding and KSPA algorithms. MRS images reconstructed using conventional gridding with fully sampled k-space data were used as the reference standard for obtaining difference maps. The reconstructed water images and the difference maps in magnitude are shown Figure 1, with the difference maps scaled up by 10 for visualization purposes. For reduction factors of 2, 3, and 4, the water images reconstructed using conventional gridding showed significant aliasing artifacts, whereas the water images reconstructed using the KSPA algorithm showed no visible aliasing artifacts.
Spectra reconstructed from representative voxels in the phantom for reduction factors of 1, 2, 3, and 4 are shown in Figure 2. The results showed similar spectra for different reduction factors with an expected enhanced noise level. Metabolite maps of Cho, Cre, and NAA maps for different reduction factors were obtained by integrating the corresponding peaks. The metabolite maps generated were similar and only NAA maps, and their difference maps with the reference standard are presented in Figure 3. Figure 4 shows the g-factor maps for reduction factors of 2, 3, and 4. Only regions in the volume of interest are shown for the g-factor maps as an accurate estimation of g-factors require referencing to the NAA signal. The mean and standard deviation of g-factors are 1.03 and 0.17 for reduction factors of 2, 1.16, and 0.22 for a reduction factor of 3, and 1.28 and 0.25 for a reduction factor of 4. As shown in the g-factor maps, there was no significant spatial variability in noise enhancement for reduction factors of 2, 3, and 4.
For the in vivo study, reconstruction with reduction factors of 1, 2, 3, and 4 were also obtained over the whole FOV using the KSPA algorithm. Similar to the phantom study, MRS images reconstructed using conventional gridding with fully sampled k-space data were used as the reference standard for obtaining difference maps. The reconstructed water images and the difference maps with the reference standard reconstructed with fully sampled k-space data using conventional gridding algorithm are shown in Figure 5 with the difference maps scaled up by 10 for visualization purposes.
Spectra reconstructed from representative voxels in the brain for reduction factors 1, 2, 3, and 4 are shown in Figure 6. The results showed similar spectra for different reduction factors with an expected enhanced noise level. As in the phantom study, similar metabolite maps of Cho, Cre, and NAA maps for different reduction factors were obtained by integrating the corresponding peaks, and only the NAA maps and their difference maps with the reference standard are presented in Figure 7. Figure 8 shows the g-factor maps for reduction factors of 2, 3, and 4. Once again, only regions in the volume of interest are shown for the g-factor maps. The mean and standard deviation of the g-factors are 1.05 and 0.15 for reduction factors of 2, 1.21, and 0.25 for a reduction factor of 3, and 1.31 and 0.28 for a reduction factor of 4. Similar to g-factor maps from the phantom study, there was no significant spatial variability in noise enhancement for reduction factors of 2, 3, and 4.
This study proves the feasibility of using the KSPA algorithm for parallel reconstruction of MRSI with spiral k-space trajectories. Using the KSPA algorithm, minimum scan times can be readily reduced by omitting part of the spiral interleaves at the cost of increased computation time and decreased signal-to noise ratios (SNRs). As demonstrated in the mean and standard deviation of g-factor values for both phantom and in vivo studies, the noise amplification is more significant as the reconstruction matrix is more ill-conditioned with high reduction factors. By exploiting the sparsity of the coefficient matrix, the KSPA algorithm efficiently calculates the reconstruction matrix, saving computing time and storage requirement. Table 1 shows the comparison of memory requirements using both direct inverse and KSPA reconstruction algorithm for an eight-channel phased array acquisition with different matrix sizes and reduction factors. The saving on the memory requirements is significant, and more memory can be saved with a higher matrix size and more receiving coils. Unlike parallel MRSI reconstruction with iterative SENSE algorithm (12) where repeated iterative procedure is used for reconstruction of each spectral point, the KSPA algorithm estimates the fully sampled Cartesian k-space data by applying one identical reconstruction matrix to the rest of repeated spiral trajectories. As a result, MRSI reconstruction time using KSPA is significantly less than that using iterative SENSE algorithm. Table 2 shows the computing times using both iterative SENSE and KSPA reconstruction algorithms for an eight-channel phased-array acquisition with different matrix sizes and reduction factors on a 2-GHz Pentium PC. The number of iterations of iterative SENSE algorithm is chosen to achieve less than 4% mean square error between estimated and acquired data. In addition, the parallel MRSI reconstruction times using the iterative SENSE reconstruction algorithm are approximately linear to the number of sampled k-space data points and the number of receiving coils. As higher spatial and spectral resolution MRSI data are acquired with more receiving coils, the saving of reconstruction time is more significant using the KSPA algorithm.
The KSPA MRSI reconstruction algorithm sees more application at high fields. With intrinsic high SNRs at high fields, MRSI is driven to reveal metabolite spectra with high spatial and spectral resolution for a large volume of interest. In these cases, the bottleneck of MRSI falls to the minimum scan time and the reconstruction time. Another MRSI application that can use the KSPA reconstruction algorithm is the multi-TE MRSI, such as 2D J-resolved spectroscopy (20–22) and CT-PRESS experiments (23,24), with arbitrary k-space trajectories. These applications, even combined with fast MRSI techniques with time varying gradients, are limited by minimum scan time.
Recently, hyperpolarized C13 has opened a new field in spectroscopic imaging (25,26). With SNRs over 10,000, this technique is capable of acquiring MRSI images that reveal metabolism in real time. The major technical challenge is to perform multiple MRSI within the T1 of hyperpolarized C13. To meet this requirement, both fast and parallel MRSI methods have to be combined, and the KSPA MRSI reconstruction algorithm can certainly be utilized in this effort.
Parallel spectroscopic imaging reconstruction with arbitrary trajectories using k-space sparse matrices is implemented and tested on MRSI with spiral k-space trajectories in both phantom and in vivo studies. Spectroscopic images very similar to that reconstructed using fully sampled k-space data are obtained using KSPA MRSI reconstruction with reduction factors up to 4. The algorithm demonstrates its flexibility of reconstruction using undersampled data on arbitrary k-space trajectories and its efficiency with significantly reduced computing times. The KSPA MRSI reconstruction can be readily applied to MRSI applications where the minimum scan time is the limitation.
Grant sponsor: Lucas foundation. Grant sponsor: NIH; grant numbers: AG18942, CA 48269, and RR 09784mtbn.