Linear projection-based chemical exchange saturation transfer parameter estimation

Isolated evaluation of multiparametric in vivo chemical exchange saturation transfer


Funding information
information. In this work, a linear projection-based multiparametric CEST evaluation method is introduced that offers fast B 0 and B 1 inhomogeneity correction, contrast generation and feature selection for CEST data, enabling reduction of the overall measurement time. To that end, CEST data acquired at 7 T in six healthy subjects and in one brain tumor patient were conventionally evaluated by interpolation-based inhomogeneity correction and Lorentzian curve fitting. Linear regression was used to obtain coefficient vectors that directly map uncorrected data to corrected Lorentzian target parameters. L1-regularization was applied to find subsets of the originally acquired CEST measurements that still allow for such a linear projection mapping.
The linear projection method allows fast and interpretable mapping from acquired raw data to contrast parameters of interest, generalizing from healthy subject training data to unseen healthy test data and to the tumor patient dataset. The L1-regularization method shows that a fraction of the acquired CEST measurements is sufficient to preserve tissue contrasts, offering up to a 2.8-fold reduction of scan time. Similar observations as for the 7-T data can be made for data from a clinical 3-T scanner. Being a fast and interpretable computation step, the proposed method is complementary to neural networks that have recently been employed for similar purposes. The scan time acceleration offered by the L1-regularization ("CEST-LASSO") constitutes a step towards better applicability of multiparametric CEST protocols in a clinical context. studies focus on the CEST effects of amide, amine, and guanidine protons related to peptides and proteins. 1 Furthermore, CEST was shown to give insights into pH 2 and the metabolite content of creatine 3 or glutamate. 4 Additionally, multiparametric CEST protocols yield contrasts related to the semisolid compartment provided by the semisolid magnetization transfer (ssMT) effect, as well as relayed nuclear Overhauser (NOE) effects that are known to correlate with protein content and conformation. 5,6 In the context of brain cancer, CEST is of clinical interest; for example, amide CEST has been shown to correlate with gadolinium enhancement, 7 and changes in NOE have been reported to correlate with histology 8 and also to be a measure for tumor therapy response. 9 However, extraction of the CEST contrast parameters of interest often requires complex mathematical modeling, for example, by nonlinear curve fitting of Bloch-McConnell 10 or Lorentzian models. [11][12][13][14] These are time-consuming, depend on initial and boundary conditions, and thus remain difficult. Presumably, this is why simple metrics like asymmetry, 15 ratios or linear interpolations of certain points in the Z-spectra 16,17 are often preferred for CEST contrast generation. Common to all of these metrics is that they describe the target contrasts as linear expressions of acquired CEST measurements at certain offset frequencies.
Linear transforms play a tremendous role in all branches of science, and MRI is no exception. Most prominently, the Fourier transform describes the decomposition of an image into spatial frequencies, which can be acquired as gradient-encoded MR signals and reconstructed using the inverse Fourier transform. In the case of a 1D Fourier transform, the contribution of a certain harmonic frequency to a signal is expressed as linear projection of the signal onto the respective harmonic signal. The harmonic signals of different frequencies thus act as basis vectors spanning a linear space of representable signals. An example of this situation is displayed in the left column of Figure 1. Such linear transforms have the advantage of being stable, fast to calculate, and insightful for theoretical analysis.
To utilize the full potential of linear transforms for CEST data evaluation, in this work we aim to find the best linear combination of acquired points in the Z-spectrum to generate a contrast as close as possible to a desired target. Such optimal linear combination weights can be found by linear regression applied to conventionally evaluated training data, in the present case using a Lorentzian fit model. Contrast generation can then be expressed analogously to the discrete Fourier transform (Figure 1, left column) as linear projection of the raw acquired data onto the respective weight vectors (Figure 1, right column). In the context of current exploration of neural networks for such tasks, 18-21 the linear transform forms the simplest learning-based approach and by its linearity allows for direct interpretation and thus also guidance for more sophisticated approaches.
Furthermore, the linear projection approach can be extended to address the scan time issue of multiparametric CEST protocols. Isolated evaluation of in vivo CEST effects usually requires sufficiently densely sampled Z-spectra, to allow for separation of concomitant exchange effects as well as correction of field inhomogeneity. Correction of B 1 inhomogeneity, which is increasingly severe at high and ultrahigh field scanners, even requires acquisition of spectral data at multiple saturation amplitude levels. 14 The number of acquired offsets and saturation amplitudes largely determines the total CEST acquisition time, as the entire sequence needs to be repeated for each amplitude and frequency offset. However, the acquired CEST data are known to be partially redundant, which can be exploited for denoising. 22 In view of such redundancies, the question arises if the desired target contrasts could be generated from only a subset of the originally acquired data. To find such subsets, we extend the linear projection approach by automatic feature selection using the well-established least absolute shrinkage and selection operator (LASSO) technique. 23 This type of L1-regularization is known from compressed sensing MRI, 24 where it is used to enforce sparsity of reconstructed MR images in certain transform domains. In the present work, the proposed "CEST-LASSO" is set up to promote sparsity of required frequency offsets and B 1 amplitude levels simultaneously, and by doing so offers a direct reduction of acquisition time.

| Linear projection
Target CEST parameters in each voxel are expressed as linear projections, that is, dot product of the acquired raw data x ! with N measurements, onto a vector of regression coefficients β ! , as In the case of multiple target parameters, the vectors of regression coefficients for each of these M targets can be assembled into the coefficient matrix B (N Â M), such that a vector of target parameters y ! can be obtained from the vector of input data x ! in a given voxel as a matrix-vector product Evaluating this expression for K different input vectors, for example, from multiple voxels, assembled to an input data matrix X (K Â NÞ, leads to the formulation as a general linear model F I G U R E 1 Analogy of discrete Fourier transform and the proposed linear projection approach for chemical exchange saturation transfer (CEST) evaluation. For the Fourier transform, a signal vector S t ð Þ is projected onto basis vectors β 1 , …, β n consisting of the respective harmonics to yield the Fourier coefficients A 1 , …,A n for different frequencies. In the case of linear CEST evaluation, acquired raw data are projected onto coefficient vectors to yield desired target contrasts like amide proton transfer (APT), nuclear Overhauser effect (NOE), and semisolid magnetization transfer (ssMT) amplitudes. The optimal coefficients can be found by linear regression from conventionally evaluated training data with the target data matrix Y K Â M ð Þ : Due to the structure of rows and columns in the defined matrices and the right-hand side multiplication of X by, the forward model Equation (3) only mixes spectral but no spatial information, such that the linear projections operate on each voxel independently. For a collection of known input-target pairs X and Y, referred to as training data from now on, the matrix of optimal regression coefficients can be obtained by solving the ordinary least squares problem, which has the global analytical solution with the Frobenius norm Á j j j j F and the Moore-Penrose pseudoinverse X þ ¼ X T X À1 X T : 25

| LASSO
The LASSO regression objective for a single scalar target parameter y can be formulated as the L1-regularized linear least squares problem, 23 In general, there is no analytical solution for b β ! LASSO . However, the optimization problem is still convex, which means that there are globally optimal solutions that can be found iteratively. 26 Depending on the choice of the regularization parameter λ, the L1-regularization leads to a sparse solution, where a certain number of coefficients are zero. As the corresponding input components do not contribute to the linear projection (Equation (1)), they can be removed entirely from the model. In the case of CEST data, this means that the corresponding measurements (at particular frequency offsets and saturation amplitudes) do not need to be acquired at all to generate the desired target contrast parameter y by linear projection, once the regression coefficients β ! are obtained.
In the case of multiple target parameters, solutions to the standard LASSO problem (Equation (5) can be calculated for each target parameter y individually. In this case, the inputs that can be removed from the model will be different for each target parameter. Consequently, the result will be differently subsampled sets of the original inputs that are all suited for just one particular target parameter, but not a common reduction scheme that simultaneously optimizes the linear prediction of multiple target parameters. To overcome this limitation, the LASSO regularization can be modified to an instance of the multivariate group LASSO, 27 namely, the row-sparsity-enforcing L2-L1 LASSO (rowLASSO). 28 A row of the coefficient matrix B corresponds to the contribution of a particular input component to all target parameters in Equation (3). The L2-L1 norm used as regularization for the rowLASSO forces entire rows of B to become zero, which means that the corresponding input components do not contribute to any of the included target parameters and can therefore be removed. Solving the rowLASSO for CEST data thus provides a single reduced list of measurements that still allows for simultaneous generation of multiple target parameters. Denoting the number of retained inputs by N red and the original number of inputs by K, the reduction factor R ¼ K=N red can be defined.
The regularization term added to the ordinary least squares objective (Equation (4)) in the case of LASSO (Equation (5)) and rowLASSO (Equation (6)) is known to introduce a bias to the obtained nonzero regression coefficients in b B rowLASSO , 29 which can deteriorate the performance of the obtained linear model when applied to new data. This bias can be removed by refitting: a reduced input data matrix X red is formed by removing all columns from X, for which the corresponding row in b B rowLASSO is zero. With that, the unregularized ordinary least squares regression According to this procedure, the regularized problem is used only for variable selection, whereas the unregularized regression provides the final coefficient estimates. Having only a reduced input vector x ! red at hand-for example, from an accelerated scan, for which the measurements rendered irrelevant by LASSO have been omitted-the contrast parameters can then be obtained as It was found throughout the present work that the refitted coefficients b B 0 rowLASSO provided better results than the original b B rowLASSO .

| Data acquisition
Data were acquired from six healthy subjects and one patient with a brain tumor (glioblastoma WHO grade IV) after written informed consent at a MAGNETOM Terra 7T scanner (Siemens Healthineers AG, Erlangen, Germany) with a 32Rx/8Tx-channel head coil (Nova Medical, Wilmington, MA). All in vivo examinations were approved by the local ethics committee.
Homogeneous saturation was realized as in 30 using the MIMOSA scheme (120 Gaussian pulses, pulse duration t p = 15 ms, interpulse delay t d = 10 ms, duty cycle DC sat = 60.56%, recovery time between previous readout and start of saturation block t rec = 1 s). Two B 1 maps were acquired for circularly polarized (CP) and 90 mode, which were combined according to 31

| Conventional evaluation
A schematic of the employed postprocessing pipeline is shown in Figure 2. All acquired 3D volumes were coregistered onto the chosen reference volume at Δω = 3.5 ppm and B 1 = 0.72 μT, to correct for subject motion, using the SPM toolbox. 34 B 0 inhomogeneity correction was applied by fitting the water peak of the Z-spectra with a smoothing spline and shifting the spectra according to the spline's minimum on the frequency axis.
The obtained frequency shift in each voxel provides a relative B 0 inhomogeneity map. During postprocessing, an interpolated baseline correction 35 with the offsets acquired at ±100 ppm turned out to yield more stable results than M 0 normalization with -300 ppm, which is why the two -300 ppm scans were not considered for further evaluations. Spectra were denoised using principal component analysis (PCA), 22 retaining the first 11 principal components. After that, two-point Z-B 1 -correction 14 was applied using the acquired spectra at both B 1 amplitudes and the relative MIMOSA B 1 map.
F I G U R E 2 Schematic of the conventional multiparametric chemical exchange saturation transfer (CEST) data evaluation pipeline (bottom), for which the proposed linear projection-based evaluation (top) forms a shortcut. In both cases, acquired CEST data are corrected for subject motion and spectrally normalized. Conventional evaluation consists of interpolation-based B 0 and B 1 inhomogeneity correction, followed by spectral principal component analysis (PCA) denoising and nonlinear least-squares fitting of a five-pool Lorentzian model. In the case of the proposed linear approach, all of these steps are accomplished by projection onto regression coefficients previously obtained from conventionally evaluated data Based on the models, [11][12][13][14]35 the resulting spectra were fitted using a five-pool Lorentzian model (water, amide, rNOE, amine, and ssMT, resulting in 16 free fit parameters), using the same initial and boundary conditions for all datasets. The output of the complete conventional evaluation including B 0 and B 1 correction and Lorentzian fitting was used to generate the training datasets for linear regression. All computations were carried out in MATLAB (MathWorks, Natick, MA).

| Linear regression and LASSO
Uncorrected but normalized Z-spectra of both B 1 levels were assembled to form the input data matrix X. Additionally, the B 1 -MIMOSA and B 1 -CP field map values in each voxel were provided as additional inputs by concatenating them to the columns of X, resulting in a total of 110 (54 + 54 + 2) input features. Lorentzian parameters from each voxel obtained by conventional fitting of B 0 -corrected data and the B 0 inhomogeneity in each voxel calculated from the non-B 0 -corrected data were likewise assembled in the target data matrix Y, resulting in 17 (16 + 1) target parameters.
The introduced formulation of the general linear model (Equation (3)) and LASSO methods (Equations (5) and (6)) assumes that both input data matrix X and target data matrix Y are centered, meaning they have column-wise mean zero. With that, there is no need for an additional constant intercept term. 29 Because of that, and to avoid any scaling issues, input and target matrices were standardized to column-wise mean zero and variance one before performing the least-squares and LASSO fits.  (5)) and rowLASSO (Equation (6)) problems were solved using the FISTA algorithm. 26,36 For the typical data matrix sizes occurring in this work (number of inputs N ¼ 110, number of targets M ¼ 17, number of training voxels K ≈ 300 000), calculation of a solution takes approximately 0.25 s. This allows calculating many solutions for increasing values of the regularization parameter λ, such that any possible number of retained inputs, meaning arbitrary acceleration factors, can be found within minutes.
The approach of generating CEST contrast maps by linear projections from the unprocessed raw Z-spectra was first validated in a healthy subject test dataset. The data of five healthy subjects were used as training set to obtain regression coefficients b B according to Equation (4), which were then applied to a sixth healthy test dataset according to Equation (3).
Predictions of the linear and LASSO method for amplitude parameters were assessed in terms of the normalized root-mean-squared error between predicted values b y and the reference parameters y with mean value y. For ΔB 0 as a target parameter, non-normalized RMSE was used, because of the mean value being close to 0.

| RESULTS
The resulting maps in Figure 3 show that for APT, NOE, and ssMT contrasts, as well as field inhomogeneity ΔB 0 , the linear projection results ( Figure 3B) preserve the general contrast of the reference maps ( Figure 3A) obtained by conventional evaluation, with NRMSE of 11%, 4.8%, and 3.9% (for APT, NOE, and ssMT amplitudes, respectively), and RMSE of 0.035 ppm for ΔB 0 . The NRMSE between the reference and projection result for the individual CEST contrasts coincides with the observed effect strengths: the strongly pronounced ssMT effect, manifested as the highest of the Lorentzian peak amplitudes, can be best predicted by the linear projection, followed by the smaller NOE and APT effects. For the amine contrast, which is the least pronounced CEST effect in the acquired data, the linear projection result is the least accurate (NRMSE = 15%). The difference maps in Figure 3C  Results for the remaining parameters of the five-pool Lorentzian model are displayed in Figure S1. From the figure it can be seen that the linear projection works better for amplitude parameters than for peak widths and positions.
The linear regression coefficients B obtained by Equation (4) allow a directly interpretable insight into how the target contrast parameters emerge from the raw input data by linear projection. In Figure 4, the regression coefficients used to generate the contrast maps in Figure  More insight into the regression coefficient patterns can be gained by simulations with artificial data, which are shown in Figures S2-S7. From those figures, it can be seen that the oscillatory sign switches of the regression coefficients effectively form weighted sums of Z-spectral points around each CEST resonance that are suitable for isolating the respective amplitude parameters from concomitant effects.
As a next step, the LASSO procedure was applied to reduce the number of inputs to the linear model. As described in the Methods section, values of the regularization parameter λ were found such that the number of retained inputs (i.e., the number of nonzero rows in B) decreases gradually from 110 (all original inputs) to one. The models were trained again on the same five healthy subject datasets and tested on the sixth healthy subject dataset, as shown in Figure 3. The rowLASSO objective was set up here to fit the targets APT, NOE, ssMT, and amine amplitudes, and ΔB 0 simultaneously, which means that the obtained reduced input lists are a compromise to yield predictions of all of these targets at the same time. A comparison of alternative rowLASSO objectives, including one reduction set for all Lorentzian parameters (amplitudes, peak widths, and positions), as well as individual reduction sets for all parameters separately (standard LASSO), is provided in Figure S8. As a general observation, including more target parameters to be predicted simultaneously by a single reduced input list yields less accurate results than individual reduced input lists for each target parameter separately in most of the cases.
An overview of the input frequency offsets that were removed by the LASSO procedure in each step is given in Figure S9. Figure 5 shows the influence of offset reduction on the NRMSE between the linear projection result from the rowLASSO-reduced inputs and reference data (red curve). Retaining 39 of the originally 110 inputs (corresponding to a reduction factor R = 2.8) still allows linear projections with NRMSE = 15%, 6.7%, 5.9%, and 21% for APT, NOE, ssMT, and amine amplitudes, respectively, and RMSE = 0.06 ppm for ΔB 0 over the test dataset.
For comparison, randomly subsampled lists of inputs were generated 100 times for each respective reduction factor, and linear regressions (Equation (4) (4) with training data generated from five healthy subject measurements. In blue, an example of the corresponding voxel-wise input is given, consisting of Z-spectra at two different saturation amplitudes and the values of the multiple interleaved mode saturation (MIMOSA) and circularly polarized (CP) B1+ field maps (not plotted). Contrast parameters in each voxel are then obtained by a simple dot product between input and coefficient vectors random reductions that yield better predictions than the rowLASSO solution. In these cases, the subsampled input lists may be incidentally well suited for the prediction of just one single target, whereas the rowLASSO objective is always to find a compromise that works as best as possible for all target parameters simultaneously. This can be confirmed by calculating the NRMSE for all standardized target parameters simultaneously ( Figure 5F), showing the superiority of the rowLASSO solution over the random reductions, especially for very high reduction factors (number of retained offsets ≤ 20).
With regard to these results, the question arises if the conventional Lorentzian fit could be used directly to generate accurate contrast maps from a reduced number of acquired frequency offsets. To investigate this, conventionally corrected Z-spectra were retrospectively undersampled (acceleration factors: R = 3/2, i.e., removing every third offset; R = 2, i.e., removing every second offset; and R = 3, i.e., retaining every third offset) and processed by the nonlinear least squares fit. The fitted results from such reduced spectra were then compared with the linear projection results from the original uncorrected spectra of both B 1 values for the same acceleration factors. As regular undersampling in the case of the Lorentzian fit permits different undersampling patterns for each acceleration factor (e.g., for R = 2 retaining offsets 1,3,5 … or 2,4,6 …), the respective pattern that yielded the best results in terms of NRMSE was chosen to make the analysis less biased towards the optimized rowLASSO.
The obtained results shown in Figure 6 illustrate that, for moderate acceleration up to R = 2, the equidistantly reduced Lorentzian fit provides lower NRMSE with respect to the reference data compared with rowLASSO. However, for an acceleration of R = 3 in the case of APT, NOE, and MT amplitudes, rowLASSO performs better, which is especially visible in the degraded NOE contrast in Figure 6E and the difference maps in Figure 6F. For amine amplitudes, the reduced Lorentzian fit performs better than rowLASSO for all the considered accelerations. Overall, it is noticeable that the reduced Lorentzian fit contrast maps fluctuate for different undersampling schemes, while rowLASSO contrasts appear stable over different acceleration factors.
The maintained low NRMSE in the case of rowLASSO reduction (NRMSE ≤ 15% in APT, NOE, and MT amplitudes for R ≤ 2.8) indicates that a linear projection reconstruction of these CEST contrast maps is still possible with only a fraction of the originally acquired frequency offsets at both B 1 amplitude levels. Indeed, the obtained parameter maps for acceleration factors R = 2 ( Figure 7C) and R = 2.8 ( Figure 7D) show no major  (6)) was set up with the target parameters APT, NOE, MT, and amine amplitudes, and ΔB 0 simultaneously deviation from the full projection result (R = 1, Figure 7B) and therefore are in agreement with the reference contrasts ( Figure 7A) as well. Only in the case of a very strong reduction ( Figure 7E, R ≈ 37, only 3 inputs retained) were the obtained APT and amine contrasts severely corrupted.
Remarkably, even for such a strong reduction, the NOE and ssMT amplitudes and ΔB 0 still correlate with the reference data to a certain extent, and the coarse spatial structure and even some anatomical contrast between gray and white matter (in the case of NOE and ssMT amplitudes) are roughly preserved. The rowLASSO-reduced regression coefficients used for generating the contrast maps in Figure 7B-E are displayed in Figure S10.
Having established the linear projection and rowLASSO method in healthy subject training and test datasets, as a next step, generalization of the method to pathology was investigated. Linear regression coefficients from the five healthy subject training datasets were obtained and applied to the glioblastoma patient dataset. Figure 8 shows the CEST results next to clinical contrasts for this patient. Interestingly, this tumor did not show typical gadolinium uptake ( Figure 8A) as expected for glioblastoma, although amide CEST ( Figure 8D, first row) showed the same hyperintensity as reported previously. 37 Despite the regression coefficients being obtained from only healthy subject data, the linear projection approach appears to generalize to tumor data, and the resulting maps ( Figure 8E) still match the general contrast of the reference Lorentzian fit maps ( Figure 8D), albeit with higher NRMSE ( Figure 8I) than in the healthy test case ( Figure 3D). In particular, the linear projections preserve the amide hyperintensity in the tumor. The same holds for the rowLASSO result with a reduction factor of R = 2.8 ( Figure 8F), although with slightly reduced contrast.
Finally, the linear projection and rowLASSO approach was applied to 3T data to assess the applicability of the method in a broader clinical context. For that, data from 18 were retrospectively re-evaluated and compared with the recent deepCEST 3T approach introduced therein. Three healthy subject datasets were used for training and a glioblastoma patient (WHO grade IV) dataset for testing. The linear projections from all 55 acquired frequency offsets ( Figure 9E), as well as a rowLASSO-reduced projection from 18 retained offsets (reduction factor R ≈ 3, Figure 9F amplitudes, and ΔB 0 , respectively. In terms of NRMSE, the deepCEST 3T neural network ( Figure 9D) performs best for all target parameters. Still, the linear and rowLASSO predictions preserve the principal contrast in the tumor region, especially the ring-shaped APT hyperintensity, which coincides with the gadolinium contrast-enhanced hyperintensity ( Figure 9A).

| DISCUSSION
In this work, it was shown that CEST parameter maps, which would conventionally have been obtained by iterative nonlinear least squares fitting of a multipool Lorentzian model, can be obtained as simple linear projections of the acquired raw Z-spectra onto regression coefficients generated from conventionally corrected and evaluated training data. The linear projection thus integrates B 0 correction, B 1 correction, and contrast generation into a single computation step. Applying the coefficient vectors is fast (fractions of seconds) compared with conventional evaluation (>10 minutes for all required steps). Compared with the training of neural networks (several hours), obtaining solutions for the supervised learning problem (i.e., calculating the regression coefficients from a training set) is also faster by several orders of magnitude. The approach was shown to translate from healthy subject training data to a tumor patient test dataset, which was also observed in recent works on neural networks for CEST evaluation 18,19 and appears plausible as long as tumor tissue spectra can be approximated by linear combinations of healthy tissue spectra. 19 The introduced linear projection approach has similarities with PCA methods, which were recently introduced in the field of CEST for denoising 22,38 and data-driven feature extraction. 39 For the latter method, test Z-spectra were linearly projected onto principal component vectors generated from a training dataset of healthy subject spectra. This was shown, similar to the finding in this work, to yield interesting tissue con- this work (i.e., the linear regression coefficients) describe the correlations between spectra and tailored target contrast parameters defined by the conventional evaluation, making the present approach a more supervised one.
Extending the proposed linear projection method by LASSO regularization provided subsets of the original measured spectral offsets at both saturation amplitude levels, which still allowed linear mapping to the target contrast parameters. By that, potential acceleration factors for the CEST acquisition of up to R = 2.8 could be achieved with only minor impact on prediction quality compared with the full linear solution (Figures 7 and 8). The LASSO reduction method is data-driven and thus specific to the chosen CEST protocol. Consequently, the offset reduction scheme depends on B 0 , B 1 , and other sequence parameters, which is why new offset lists should be generated if the acquisition protocol has changed.
As shown in Figure 6, conventional Lorentzian fitting on equidistantly undersampled spectra can also be used for acceleration. However, for acceleration R = 3 in the case of APT, NOE, and MT, rowLASSO performed better than Lorentzian fitting. Even more conservatively, Goerke et al. 40 manually reduced offsets for Lorentzian fitting and reported a maximum reduction of 19% (corresponding to acceleration R ≈ 1.23) to still yield acceptable results. Better performance of undersampled conventional Lorentzian fitting might be achievable by dedicated offset list optimization instead of equispaced undersampling, but in the case of nonlinear least squares Lorentzian fitting, the required combinatorial optimization would suffer from computational complexity. The observed dependence on the sampling pattern (every first or every second removed, etc.) also indicates that detailed analysis of B 0 shifts is needed here, as these shifts have a similar influence on the effective sampling pattern.
The performance of CEST-LASSO for strong spectral undersampling is surprising, as to properly define the Lorentzian fit curve, many spectral points are necessary. However, (i) several Lorentzian parameters (e.g., the peak width and height) are correlated and the reduced training data effectively defines a subspace of the original data, so less information is needed than the full Lorentzian parameter space; and (ii) we also map solely on amplitude parameters, which again is just a part of the Z-spectral information, thus redundant data can be removed. If all parameters are reconstructed simultaneously, the performance decreases ( Figure S8). approach, by contrast, does not apply sparsity constraints to the image encoding, as the k-spaces of each acquisition are conventionally sampled and reconstructed, but to the number of CEST acquisitions. With that, the potential for CEST acceleration can be even greater than in the case of methods that focus only on acceleration of image encoding, as the saturation block of a CEST sequence can account for more than 80% of the total acquisition time, 41 depending on the protocol.
Most in vivo CEST studies use evenly distributed or manually tailored frequency sampling schedules that are not necessarily optimal.
Addressing this issue, Tee et al. 42 apply an optimal design approach that yields optimized sampling schedules for amine quantification. In contrast to the present work, their objective is to find optimal sampling points for a predefined fixed number of acquisitions, whereas CEST-LASSO finds subsets of arbitrary length from a fixed initial sampling schedule. Consequently, the optimal design approach optimizes parameter quantification for fixed acquisition time, whereas CEST-LASSO finds compromises between parameter quantification and shorter acquisition time. The optimal design method relies on differential sensitivity analysis performed on numerically simulated spectra. Consequently, transferring this approach from phantom to in vivo CEST imaging would require an accurate Bloch-McConnell model for in vivo spectra, which is currently still not well established, as reported quantification results vary strongly among research groups. 10,[43][44][45] By contrast, CEST-LASSO, as a more data-driven approach operating directly on measured spectra and evaluation results, is readily applicable to in vivo data.
The initial example of Figure 1 employing the analogy to Fourier transform hints at previous work performed by Yadav et al.,46 where an actual Fourier transformation on Z-spectra was performed. This should not be confused with our approach, which finds a new data-driven linear basis, whereas Yadav et al. used the Fourier basis directly to remove low and high frequency components in the time domain.
The Lorentzian fitting method employed for generating reference data is only one of the many approaches for analyzing CEST data. Being a nonlinear least squares fit with 16 free parameters, it suffers from typical drawbacks, such as depending on initial and boundary conditions and being susceptible to noise and fluctuations. In particular, the amine contribution at +2 ppm can only be poorly extracted by both the Lorentzian fit and the proposed linear method, which also manifests in a low correspondence between the amine signals from CEST-LASSO and the reference data. The reason is that the corresponding features in the spectra are only weakly pronounced due to the relatively low applied B 1 amplitudes (0.72 μT and 1.00 μT) compared with the relatively high reported exchange rates of 700 Hz-10 kHz. 47 The resulting weak CEST labeling, 48 together with the low reported amine pool size, 49 leads to an overall small effect strength. Furthermore, the amine resonance at $2 ppm is close to the water peak where spillover dilution effects occur. The choice of five pools is of course an oversimplification of the complex in vivo situation, which comprises many concomitant exchange effects (e.g., of -OH groups), which are not taken into account here. Nevertheless, similar five-pool Lorentzian models have been successfully applied as estimators of in vivo CEST effects before. 11,14 Similar to the observations made for the tumor patient dataset in this work, hyperintensities in tumor areas have been observed before for APT and amine contrasts obtained by Lorentzian fitting. 11,50 Although in the present work the MIMOSA saturation scheme was employed for achieving more homogeneous saturation at 7 T, the proposed linear projection and LASSO method does not rely on this scheme and can be used for other CEST saturation schemes as well, as demonstrated for a 3T protocol without MIMOSA in Figure 3.

| Explainable AI
Recently, deep neural networks have been suggested for CEST data evaluation. [18][19][20][21] Like most deep learning methods, these are ascribed to be "black boxes" as their nested nonlinear structures do not allow interpretable insight into how inputs are mapped to targets. This makes careful assessment of such methods difficult. Addressing this issue, the recent trend of "explainable AI" aims at making machine learning systems more accessible for human interpretation. 51 In this sense, the presented linear projection approach could be considered as a simple step towards "explainable AI" applied to CEST MRI, as linear models are intrinsically interpretable, 51 which was demonstrated in Figure 4 and In general, the relationship between CEST contrast parameters and measured Z-spectra is nonlinear and can be described by solutions of the Bloch-McConnell differential equations or simplified Lorentzian models, as employed in this work. Given that, the proposed linear model can be interpreted as a first order approximation in the sense of a Taylor expansion of the generally unknown nonlinear function solving the inverse problem of inferring the parameters of interest from the measured data.
Resorting to a linear model instead of a nonlinear neural network thus enables interpretability (and much faster computation), however, at the expense of lower prediction performance, as shown in comparison with the deepCEST method 18 in Figure 9. This represents a general trade-off between model complexity and capacity on the one hand and simplicity and interpretability on the other. 51 Given the advantages of the linear approach in terms of speed and insight, first applying a linear model before advancing to more complex nonlinear models like neural networks might be a reasonable choice and is recommended for other learning-based CEST evaluations.
Linear models are more robust to errors in the input data than nonlinear models like neural networks. 52 Small fluctuations, which could be caused, for example, by motion in the case of CEST acquisitions, can in the worst case translate to arbitrarily high errors in the parameters predicted by a neural network. By contrast, the prediction error of a linear model in such cases is always bounded and can be estimated from the "slope" of the model (i.e., the regression coefficients). Furthermore, as there are no cross terms between different spectral input points (e.g., products like Z Δω 1 ð ÞÁZ Δω 2 ð Þ), a single corrupted frequency offset can only have limited impact on the result, which is not the case for highly entangled models like neural networks.

| CONCLUSION
Multiparametric CEST contrasts including field inhomogeneity correction can be well approximated by a simple linear projection of the acquired uncorrected Z-spectra onto regression coefficients fitted from conventionally evaluated data. The method translates from healthy to tumor patient datasets and is fast and interpretable, the latter being in contrast to neural networks employed for similar purposes.
Extending the approach by L1-regularization yields reduced frequency offset acquisition schedules offering a potential reduction of total scan time by factors of up to 2.8 with only moderate quality losses, compared with directly fitting subsampled Z-spectra. This could help make multiparametric CEST protocols more viable for clinical application.
APPENDIX A.