Patch-Based Super-Resolution of MR Spectroscopic Images: Application to Multiple Sclerosis

Purpose: Magnetic resonance spectroscopic imaging (MRSI) provides complementary information to conventional magnetic resonance imaging. Acquiring high resolution MRSI is time consuming and requires complex reconstruction techniques. Methods: In this paper, a patch-based super-resolution method is presented to increase the spatial resolution of metabolite maps computed from MRSI. The proposed method uses high resolution anatomical MR images (T1-weighted and Fluid-attenuated inversion recovery) to regularize the super-resolution process. The accuracy of the method is validated against conventional interpolation techniques using a phantom, as well as simulated and in vivo acquired human brain images of multiple sclerosis subjects. Results: The method preserves tissue contrast and structural information, and matches well with the trend of acquired high resolution MRSI. Conclusions: These results suggest that the method has potential for clinically relevant neuroimaging applications.


INTRODUCTION
Magnetic resonance spectroscopy imaging (MRSI) of the brain provides information on the chemical composition of tissues and may reveal underlying metabolic changes that are not visible on conventional magnetic resonance imaging (MRI) such as T1-weighted MRI. For example in multiple sclerosis (MS), MRSI can show the extent of damage in normal appearing white matter and therefore is more sensitive in detecting pathological changes compared to conventional MRI (Filippi and Agosta, 2010;Rovira et al., 2013). Additionally, metabolite concentrations in normal appearing brain tissue are better correlated with clinical scores than MRI lesion load (Sajja et al., 2009).
Despite the fact that MRSI provides attractive complementary information with respect to conventional MRI, it is still not widely used in practice. Possible reasons include (a) low signal to noise ratio, (b) poor spatial resolution, (c) prolonged acquisition time, (d) lack of standardized protocols, (e) and limited quality control. Since the metabolite concentrations are very low, a large voxel size is required in order to capture sufficient signal, resulting in a poor spatial resolution (typically around 1 cm 3 ). Additionally, this results in a broad point spread function which further causes spectral signal leakage in the neighboring voxels. This limits the ability to differentiate among spatial features of interest, e.g.,between different tissue types. High spatial resolution MRSI cannot be acquired in the clinical routine with standard acquisition sequences as it requires a long scanning time, e.g.,in the range of tens of minutes.
Therefore, there is a need to design alternative, nonconventional approaches to create high resolution MRSIbased metabolic images within a limited acquisition time. The acquisition time can be reduced by optimising the conventional phase encoded MRSI protocol (Posse et al., 1995;Pohmann et al., 1997;Adalsteinsson et al., 1998) or by using parallel imaging reconstruction based techniques (Dydak et al., 2001;Lin et al., 2007;Otazo et al., 2007;Banerjee et al., 2009) or combination of both Dreher et al. (2011), and Posse et al. (2013) for a detailed review. Conventional MRSI protocols contain a series of periodic spatial and spectral encoding schemes. These could be optimized by reading out the spatial information during spectral acquisition. This optimized technique forms the basis of the echo planar spectroscopic imaging (EPSI) (Posse et al., 1995) which provides an improved spectral and spatial resolution compared to conventional MRSI protocols for fields up to 3T. Derived from the EPSI technique, spiral MRSI (Adalsteinsson et al., 1998;Andronesi et al., 2012;Bogner et al., 2014) uses a spiral trajectory in k-space, which allows faster acquisition than EPSI however at the cost of increased complexity in the data reconstruction. Parallel imaging reconstruction techniques are based on the principle of acquiring multiple undersampled k-space data using phased array coils. This undersampled data is then reconstructed into a single image using techniques like SENSE (Pruessmann et al., 1999) or GRAPPA (Griswold et al., 2002). Since optimized MRSI protocols like EPSI and parallel imaging are independent in nature, they can be combined to further reduce the acquisition time (Dreher et al., 2011). This reduction in the acquisition time generally comes at the cost of decreased signal to noise ratio (SNR) in the data. In spite of advances in speeding up the MRSI acquisition these techniques are yet not clinically feasible as the time required to attain the spatial resolution of conventional MRI sequences such as T1-weighted imaging still remains very high, although approaches that combine phase-encoding and EPSI (Ma et al., 2017) appear quite promising.
An ultra-high magnetic field strength such as 7T (used for non-routine research) offers more SNR, but the inherent limitations (e.g., related to specific absorption (SAR) or B 0 /B 1 homogeneity) have complicated high-resolution MRSI acquisition until the recent years. Free induction decay (FID) MRSI (Bogner et al., 2012;Považan et al., 2015;Strasser et al., 2013) with short acquisition delays was proposed to avoid SARintensive localization schemes and SNR-loss due to T2 * decay. The application of parallel imaging like SENSE (Zhu et al., 2013;Kirchner et al., 2015), GRAPPA , or CAIPIRINHA  allowed acceleration factors up to 9. Combining parallel imaging with a short TR of 200 ms, the acquisition of high-resolution MRSI within a 128 × 128 mm 2 field of view with a voxel volume of 23 µL and full slice coverage in around 10 min was successfully demonstrated in Hangel et al. (2016).
An alternative is to increase the MRSI resolution by superresolution (SR) techniques. Super resolution methods include k-space based reconstruction methods that improve the spectral quality of reconstructed high resolution MRSI data using high resolution spatial features, such as edges, from other imaging modalities, in particular T1-weighted MRI (for a detailed review see Kasten et al., 2016). Super-resolution methods can be categorized into two sub-groups: methods based on a linear regression framework, or on a Bayesian framework. Linear regression based SR methods (LRSR) assume that the acquired data can be explained by a linear combination of a set of independent variables. Therefore, they aim at optimising the coefficients of these independent variables such that the error between the predicted and the actual measurements is minimized. One such method is spectral localization by imaging (SLIM) (Hu et al., 1988) which assumes that the high resolution MRSI consists of a linear combination of L anatomical compartments that are spectrally homogeneous. Although easy to implement, the assumption of spectral homogeneity may not hold true. To cope with this, two extensions have been proposed: generalized SLIM (GSLIM) (Liang and Lauterbur, 1991) and SLIM with explicit B 0 field inhomogeneity compensation (BSLIM) (Khalidov et al., 2007). GSLIM uses spatial Fourier harmonics to absorb any spatially dependent spectral variations and BSLIM assumes that the spectral variations are solely due to local (static) field inhomogeneity, therefore, requires an additional B 0 map for its correction. It has been observed that simple LRSR methods not always result in physically plausible solutions, therefore, an additional term commonly known as "regularizer" is added to the optimization problem so that the reconstruction remains well-behaved. For example, in Jacob et al. (2007), a local B-spline basis function is added as a regularizer to the GSLIM model to allow local intensity variations within each compartment. In Haldar et al. (2008), a smoothness regularization term (controlled by pre-computed anatomically derived weight factors) has been added to penalize the local spectral variation between neighboring voxels. The second subgroup of methods are based on Bayesian theory and model the reconstruction in the k-space as a likelihood function where the anatomical information acts as prior knowledge to estimate the optimized model parameters via the expectationmaximization algorithm. For example in Bao and Maudsley (2007), the likelihood function consists of a combined spectralspatial model where the tissue segmentations acts as prior information in estimating additional high frequencies in k-space. As an extension of this work, the likelihood model in Kornak et al. (2010) also addresses the spectral fitting problems and additionally uses prior information on the relationship between tissue segmentation and spatial metabolite distribution. These Bayesian based models are very complex in nature and have many variables to be optimized that often result in locally optimal solutions.
Contrary to the k-space based super-resolution approaches presented above, patch-based super-resolution (PBSR) methods have been used to reconstruct high resolution T2 and diffusion weighted MR images (Manjón et al., 2010a,b;Rousseau and The Alzheimer's Disease Neuroimaging Initiative, 2010;Coupé et al., 2013) in the spatial domain using a high resolution T1-weighted image as reference. PBSR methods are based on the principle of image redundancy (ability to find similar patches in an image) and aim at finding similar voxels [defined using self similarity or using the high resolution image (e.g., T1-weighted image)] in the neighborhood of the central voxel that could be used to guide the reconstruction of the central voxel. As opposed to classical interpolation methods (nearest neighbor, linear interpolation and B-splines), PBSR method provides better tissue contrast and thus yields better image quality reconstructions.
In view of good performance of PBSR methods in upsampling T2 and diffusion weighted MR images, we propose a multi parametric patch-based super-resolution method that acts directly on quantified MRSI metabolite maps from lowresolution MRSI data. In addition to the conventionally used high resolution T1-weighted image, a FLAIR image is also used for lesion segmentation in multiple sclerosis subjects. The reconstruction process is guided by the intensities of T1weighted and FLAIR images along with the brain segmentation, thus providing better estimation of brain tissue boundaries. To the best of our knowledge, a PBSR technique has never been applied before for upsampling MRSI-based metabolic images. We compare our method against the classical upsampling methods (nearest neighbor, linear interpolation and B-splines) on (a) an acetate image from a phantom, which contains different acetate proportions in different spatial locations, (b) simulated MS brain and quantified N-acetylaspartate images and (c) real images of MS patients with quantified N-acetylaspartate and myo-inositol images.

METHODS
The patch-based super-resolution pipeline upsamples the low resolution quantified MRSI metabolite map using the high resolution T1-weighted and FLAIR MR images. In the next sections, we describe the pipeline that has three steps: (1) a preprocessing step, including (a) MRI brain tissue segmentation, (b) metabolite quantification, (c) super-resolution initialization, (2) a reconstruction step that estimates the metabolite concentration at each higher resolution MRSI voxel using the tissue segmentations along with image intensities of bias corrected T1-weighted and FLAIR images, (3) a mean correction step that rectifies the estimated metabolite concentration in every high resolution voxel by taking into account a point spread function (PSF) and the error toward the corresponding lower level voxel's metabolite value. The method iterates between step 2 and step 3 such that the corrected metabolite concentration map from previous iteration is used to initialize the metabolite concentration prior map for the current iteration. The convergence of our method is detected when the relative metabolite concentration difference between the current and previous iteration is negligible. It takes generally five iterations for the algorithm to converge. An overview of the pipeline is shown in Figure 1.

MRI Brain Tissue Segmentation
The anatomical brain MR images are segmented into gray matter (GM), white matter (WM), cerebrospinal fluid (CSF) and lesions using MSmetrix (Jain et al., 2015). The method iteratively segments the T1-weighted image into GM, WM, and CSF, segments the WM lesions on the FLAIR image as an outlier to normal brain using Mahalanobis distance, and performs lesion filling in the T1-weighted image to improve tissue segmentation at next iteration. After convergence, soft segmentations of GM, WM, and CSF are created together with binary lesion segmentation. In addition, bias corrected T1weighted and FLAIR images are also produced.

Metabolite Quantification
The acquired MRSI data, after Fourier transformation from kspace to the spatial domain, can be fitted in a voxel-by-voxel fashion using a simulated metabolite basis set, matching the acquisition parameters (e.g., echo time) of the given signals. The fitting can be based either on time-domain MRSI signals, or after Fourier transformation to the frequency domain. In this work we used the SPID software 1 (Van Cauter et al., 2013) and LCModel 2 for metabolite quantification. The output consists of a set of metabolite images ("metabolite maps"), as many as metabolites in the basis set, and having the same spatial resolution as the input MRSI grid.

Super-Resolution Initialization
The super-resolution pipeline is initialized by (1) upsampling the low resolution metabolite map using linear interpolation (scale factor of 2). This forms the starting value of the high resolution MRSI metabolite map, (2) by aligning the MRI and MRSI images. As the MRI and MRSI images are acquired on the same scanner and in a sequential manner, their affine transformations are used to project the MRI image along with the tissue segmentations into the space of the upsampled MRSI image. The total affine transformation that relates these two image spaces is defined as:

Background
Let y be a noise-free low resolution metabolite map and x be the unknown corresponding high resolution metabolite image defined over a high resolution voxel space . Theoretically, where H is a known blur and downsampling operator. In order to reconstruct x given y, the reconstruction model is typically formulated as an Euler Lagrange energy minimization problem: where the parameter λ controls the amount of regularization and the R(x) term preserves geometry and regularity and is defined as: is the conditional reconstruction of x i based on voxel intensity values in the "search volume" ζ i . We assume that the term E(x i |ζ i ) is consistent under strictly stationary conditions, i.e., this stationary condition amounts to saying that, as the resolution of the image increases, there are many similar patches for all the details of the image. Also, when the regions move far from each other, their correlation decreases (Buades et al., 2005). Both assumptions are approximately met in the case of MRSI. If E(x i |ζ i ) is modeled under these assumptions, then x i E(x i |ζ i ) which means R(x) 0 (Roussas, 1990). Therefore, the two terms of the energy functional can be decoupled and asymptotically approximated as the reconstruction and mean correction steps of the iterative super-resolution algorithm. Algorithm: SR(y = low res image; x = starting value for high res image) where p is the index of the corresponding voxel at low resolution of which x i is a part.

Reconstruction
The reconstruction term, E(x i |ζ i ) in its general form is defined as: where the weight w(x i , x j ) defines the contribution of neighborhood voxel x j in the reconstruction of x i . Mathematically, the general form of w(x i , x j ) is defined as Coupé et al. (2013): where Z i is a normalization constant such that N i is the intensity vector of the local neighborhood N i of length N, h i is the standard deviation of neighborhood voxel intensities in N i . N i is empirically defined as a cube of size 8 × 8 × 8 voxels around the centre voxel x i . From Equation (6), we observe that if N j is similar toN i , more weight is given to x j in the reconstruction of x i . Since we use both high resolution MRI (T1-weighted and FLAIR) and low resolution MRSI, the weights are refined for both modalities taking into account the prior information: Refining the weights for high resolution MR images w MRI (x i , x j ) are defined such that the weights in Equation (6) accommodate the prior knowledge on brain tissue segmentations that were down-sampled to match the MRSI resolution: Frontiers in Neuroscience | www.frontiersin.org where k ∈ K = {GM, WM, CSF}, p .,k denotes the probability that the voxel belongs to a particular tissue class k.

Refining the weights for low resolution MRSI
The w MRSI (x i , x j ) weights have the form (6), and the only unknown term is h i . Assuming white noise Gaussian distribution, h i can be calculated using the pseudo residual ǫ i : where N * i is the 6-neighborhood (4 in-plane and 2 out-of-plane neighbors) of x i . Now, where M is the number of voxels in ζ i that have non-zero weight in the reconstruction of x i . After computing weights for both modalities, E(x i |ζ i ) is defined similar to Rousseau and The Alzheimer's Disease Neuroimaging Initiative (2010): where α(x i ) is a weighing term between the MRI and MRSI weights. The choice of α(x i ) is application specific. For example, in case of non-pathological studies, the reconstruction process could be driven by MRI only, and therefore, α(x i ) can be set to zero. In case of pathological cases like MS, MRSI is more sensitive to pathological changes such as lesions compared to MRI. As we focus on MS subjects in this paper, α(x i ) is defined using binary lesion segmentation which was obtained in the preprocessing step (see Section 2.1.1). In particular, the weight α is defined for each voxel x i such that the α(x i ) = 1 implies that it is a lesion and thus the reconstruction process is driven by MRSI only. Otherwise, it is driven by MRI (see Equation 10).

Mean Correction
In the mean correction step of algorithm, the reconstructed values are first convolved with the PSF, which is typically a sinc function sinc(πx), where x = 1 nominal voxel size = no of phase encoding steps field of view (De Graaf, 2013;Posse et al., 2013). Then, the average of the PSF corrected reconstruction values that compose the low resolution voxel y p (of which x i is part) must be close to the corresponding original value of the low resolution image. This corresponds to a sinc and boxcar operator for H in Equation (2), but it could be easily replaced by a general smoothing and downsampling operator H.

Phantom Data
A cylindrical phantom (180 mm diameter, 168 mm height) is used, consisting of 7 small cylinders each having a diameter of 40 mm and depth of 150 mm, going through the thickness of the phantom as shown in Figure 2. The first cylinder was kept empty for orientation purposes, followed by the second cylinder having acetate (Ace) concentration of 6 mM. Then, the Ace concentration was increased by 2 mM for each next cylinder resulting in a concentration of 16 mM for the last cylinder. The rest of the empty space in the phantom was filled with water. MR imaging was performed on a 3T whole body scanner (Philips Achieva, Best, the Netherlands). The protocol contained two sequences: 3D T1-weighted fast field echo (FFE) sequence (TR = 6.6 ms, TE = 3.12 ms, FA = 9 • 180 × 180 mm 2 FOV, 3 axial slices, 0.8 × 0.8 × 10.0 mm 3 voxel resolution) and the 1 H MRSI point resolved spectroscopy (PRESS) sequence (TR = 1500 ms, TE = 120 ms, 300 × 300 mm 2 FOV, 3 axial slice, 10 × 10 × 10 mm 3 voxel resolution). MRI individual cylinder segmentation was performed manually using the Slicer tool (version 4.3.1) 3 and the super-resolution method was then performed at three different higher resolution voxel sizes: 5 × 5 × 5 mm 3 , 2.5 × 2.5 × 2.5 mm 3 and 1.25 × 1.25 × 1.25 mm 3 . As the phantom is homogeneous, we present in this paper only the phantom data results at 1.25 × 1.25 × 1.25 mm 3 voxel resolution. MRSI data was quantified using the SPID software 1 (Van Cauter et al., 2013) to obtain metabolite images for Ace.

Evaluation Criteria
We compare our method's accuracy against three conventional interpolation techniques: nearest neighbor, linear interpolation and B-splines. The nearest neighborhood (NN) function assigns the value of the new point with the closest old neighbor value. Linear interpolation (LIN) function linearly interpolates the new point between the old points. B-splines interpolation (BS) function uses piecewise cubic polynomials of degree 3 to interpolate the new point using four old points (two on each sides of the new point).
In order to avoid the need of absolute quantification for the phantom data, we compare ratios between the estimated metabolite values obtained by each method in cylinders 5 and 2, and 7 and 3, respectively, as we know that the ground truth ratio is equal to 2. Table 1 presents the quantitative results of Ace concentration ratio between cylinders 5 and 2, and 7 and 3 respectively, both at low resolution (acquired at 10 × 10 × 10 mm 3 ) and at high resolution (reconstructed at 1.25 × 1.25 × 1.25 mm 3 ) for all methods. On average, PBSR and LIN seem to preserve better median metabolite ratio compared to other methods. However, the visual assessment for all methods (Figure 2) shows that contrary to conventional interpolation methods, PBSR reduced the partial volume effects considerably by incorporating tissue segmentation information. Therefore, overall PBSR provides better results.

Evaluation Criteria
We compare our method's accuracy against nearest neighbor, linear interpolation and B-splines interpolation techniques. For validating the accuracy on the simulated brain dataset, the global image similarity is measured with the structural similarity index (SSIM) (Wang et al., 2004) which is more compatible with the human visual assessment and is defined as: where (µ x , µ y ) and (σ x ,σ y ) are the respective means and standard deviations of image x and y. In our case, the metabolite concentrations are not zero or close to zero within brain, therefore, the stability constants c 1 and c 2 were set to zero. Additionally, statistical differences in the metabolite concentrations between the lesions and in white matter surrounding lesions are tested using Welsh's t-test (Welch, 1947) and the magnitude of this difference i.e., the effect size is measured using Cohen's d (Cohen, 1988). Figures 3-5 show the respective visual results for all methods on mild, moderate and severe MS subjects. On close inspection, it can be seen that PBSR maintains better tissue contrast compared to other methods. For PBSR, this contrast is also more   evident between lesions and their surrounding white matter, both visually and quantitatively (see Table 2). Moreover, the median value in white matter surrounding lesions for PBSR is closest to the actual value considered for WM, 25, compared to other methods. For lesions, the median value for all methods is close to 20 with slightly better effect size for PBSR. SSIM score for all methods is very high with PBSR being marginally higher.

Evaluation Criteria
We compare our method's accuracy against nearest neighbor, linear interpolation and B-splines interpolation techniques. For validating the accuracy on the human dataset acquired in vivo, SSIM score is used. Also, statistical differences in the metabolite concentrations between the lesions and their neighboring WM tissue are tested using Welsh's t-test and the effect size is measured using Cohen's d.

Accuracy on Human Dataset Acquired In vivo
Figures 6, 7 show the high resolution NAA and myo-Ins maps respectively, for all methods on a representative example. In contrast to the results on the simulated brain dataset, visually it is very difficult to see a correlation between any method and the acquired high resolution maps. However, quantitatively, a contrast between lesions and surrounding white matter can be computed ( Table 3). The acquired high resolution shows a significant contrast between lesions and their surrounding WM, which is not seen at acquired low resolution. For NAA, this significant contrast can only be seen in PBSR. For myo-Ins, this significant contrast is highest for LIN followed by BS, NN and PBSR. The average SSIM score is highest for PBSR for both NAA and myo-Ins compared to the other methods, which have equal average SSIM scores.

DISCUSSION AND CONCLUSION
In this paper, we presented a patch-based super-resolution method for upsampling the low resolution quantified metabolite maps. The method incorporates the high resolution spatial information from T1-weighted and FLAIR MR images in guiding the reconstruction process. The method iteratively estimates and corrects for the metabolite concentration at a high resolution.
In contrast to the k-space based techniques where anatomically derived information is used to improve the spectral quality of high resolution MRSI data, in our case, the method reconstructs each central voxel x i using a weighted average of voxels that have similar tissue composition as the central voxel in the search volume. Moreover, as MRI and MRSI have complementary information, separate weights are defined for each modality which are regulated with a parameter α(x i ). Depending on the application, α(x i ) gives the flexibility to control the influence of each modality in the reconstruction process. In this paper, when α(x i ) = 1, the reconstruction is guided by MRSI only resulting in a smooth image and when α(x i ) = 0, the reconstruction is guided by MRI alone resulting in a better tissue contrast compared to conventional interpolation techniques (see Tables 2, 3). The reconstruction quality is affected by the potential registration errors and the reliability of the input quantified metabolite  map. Additionally, the method could not be expected to recover small scale features that were not sufficiently picked up by the low resolution measurement. On a similar note, the use of tissue segmentations to guide the reconstruction might affect metabolite values at tissue interfaces. In particular, through such an effect, some hyper intense voxels are introduced near the brain boundary, requiring further investigation. Finally, if the difference between the image resolution of MRI and MRSI images is big (scale factor more than 4), then the downsampled tissue segmentations for super-resolution process have high partial volume effect and thus do not carry much information on tissue type. That is why in this study, we limited the in-plane resolution of acquired low resolution MRSI to ∼ 3.4 × 3.4 mm 2 .
Our method can also be implemented in a multi-scale fashion. However, validating such high resolution reconstructed image would be a problem in case of real datasets (e.g., Experiment 3) where acquiring such high resolution image is not possible due to practical limitations. Indeed, multi-scale super-resolution method can easily be tested on a simulated data, and, as our results show better tissue contrast at the current high resolution, we expect similar or better results at even higher resolution. Our main aim in this paper was to validate the method on real datasets and thus we opted not to use the multi-scale super-resolution. Among the methods proposed in the literature for superresolution, our approach has some similarities to Rousseau and The Alzheimer's Disease Neuroimaging Initiative (2010), which is also based on patch-based framework. In contrast with Rousseau and The Alzheimer's Disease Neuroimaging Initiative (2010), which uses a neighborhood averaging strategy to define MRI weights, we used brain tissue segmentations. In Rousseau and The Alzheimer's Disease Neuroimaging Initiative (2010), α(x i ) defines the correlation between high resolution and low resolution images and is adapted in each iteration. In our case, α(x i ) is defined using binarized lesion segmentation and is fixed throughout the process.
One of the biggest problems that we encountered in this study was having ideal data for method validation. Challenges include (a) absence of high resolution ground truth data, (b) lack of a publicly available datasets, and (c) difference between low and high resolution acquired MRSI data. Acquiring high resolution ground truth data in reality is not feasible because low metabolite concentration requires large voxel size to capture sufficient signal, and prolonged acquisition time. Although several group studies have reported metabolite concentration in normal and pathological cases like MS (Sajja et al., 2009;Bogner et al., 2012), conclusions are not easily generalizable when it comes to individual patients. Lack of a standardized acquisition protocol for high spatial resolution MRSI hampers the creation of publicly available datasets against which state-of-the-art methods' performances can be compared. Moreover, most of the state-of-the-art super-resolution methods for MRSI improve MRSI data using sophisticated reconstruction techniques (kspace based), which make their reproducibility for comparison very difficult. However, from methodological point of view, our method requires fewer variables to be optimized simultaneously which simplify the underlying optimization problem. For example, in k-space based methods, the MRI spatial parameters (e.g., bias field inhomogeneity, tissue class information) and the MRSI spectral parameters (e.g.,metabolite amplitudes, B 0 shifts, zero-order phase and lineshape parameters) have to be optimized simultaneously, which may result in a sub-optimal solution. In our case, we deal with the MRI spatial parameters and MRSI spectral parameters separately. Bias field correction and tissue class segmentation are performed by MSmetrix (2.1.1) and the interaction between neighborhood voxels is exploited in the proposed super-resolution method. The spectral parameters are estimated by state-of-the-art metabolite quantification methods such as LCmodel or SPID.
In human MRSI data acquired in vivo, the median metabolite concentrations acquired at low resolution are greater than those at higher resolution (see Table 3), although the same quantification method [LCModel] 2 and parameters have been used to quantify metabolic maps at both resolutions. This generates a degree of uncertainty in the interpretation of the performance analysis. This issue was not present in the simulated data, which explains a significant decrease in the SSIM score from simulated brain dataset's results (see Table 2) to human dataset acquired in vivo (see Table 3). This decrease can also be explained by the fact that simulated low resolution MRSI data has twice the resolution of low resolution MRSI human dataset acquired in vivo, resulting in a loss of structural information. Finally, it may seem like the conventional interpolation methods provide better results than PBSR for myo-Ins (see Table 3), however, this is not completely true. Four out of five MS subjects in the human dataset acquired in vivo, have low lesion load in the acquired MRSI plane and the remaining subject has high lesion load. If this subject is removed, none of the methods show any statistical difference in the metabolite concentrations between the lesions and their neighboring WM tissue for myo-Ins. However, this statistical difference still holds true for NAA. Tuning the parameter α(x i ) may address this issue and will be explored.
In conclusion, we presented a patch-based super-resolution approach for upsampling the low resolution quantified metabolites maps using T1-weighted and FLAIR MR images. The proposed method preserves tissue contrast and structural information compared to conventional interpolation methods, and matches well with the trend of acquired high resolution MRSI. These results suggest that the method has potential for clinically relevant neuroimaging applications.