RESUME: Turning an SWI acquisition into a fast qMRI protocol

Susceptibility Weighted Imaging (SWI) is a common MRI technique that exploits the magnetic susceptibility differences between the tissues to provide valuable image contrasts, both in research and clinical contexts. However, despite its increased clinical use, SWI is not intrinsically suitable for quantitation purposes. Conversely, quantitative Magnetic Resonance Imaging (qMRI) provides a way to disentangle the sources of common MR image contrasts (e.g. proton density, T1, etc.) and to measure physical parameters intrinsically related to tissue microstructure. Unfortunately, the poor signal-to-noise ratio and resolution, coupled with the long imaging time of most qMRI strategies, have hindered the integration of quantitative imaging into clinical protocols. Here we present the RElaxometry and SUsceptibility Mapping Expedient (RESUME) to show that the standard acquisition leading to a clinical SWI dataset can be easily turned into a thorough qMRI protocol at the cost of a further 50% of the SWI scan time. The R1, R2*, proton density and magnetic susceptibility maps provided by the RESUME scheme alongside the SWI reconstruction exhibit high reproducibility and accuracy, and a submillimeter resolution is proven to be compatible with a total scan time of 7 minutes.

Numerous studies have demonstrated the clinical relevance of this technique across a wide range of pathological conditions [6], including neuroinflammation, neurodegeneration and head trauma. Moreover, it is valuable in the detection of cerebral micro-bleeds in aging population, in the follow-up of hemorrhagic infarctions, and in the assessment of cavernous malformations, venous thrombosis, and calcium or iron deposition, as thoroughly reviewed in [7][8][9].
Although several MRI Gradient Echo (GRE) signals (Steady-State Free-Precession-SSFP [10], Echo-Shifted-ES [11], etc.) acquired with different spatial encoding strategies PLOS  (Projection Acquisition-PA [12], Periodically Rotated Overlapping ParallEL Lines with Enhanced Reconstruction-PROPELLER [13], Echo Planar Imaging-EPI [14], etc.) can be processed to obtain an SWI dataset, it is common practice to rely on a standard 3D single spoiled GRE with a relatively long echo time (usually T E * 20 ms at 3 T [15]). Such choice, however, seems suboptimal, as at least one more echo at short T E could fit within an otherwise dead-time of the sequence. Indeed, even at moderately low receiver bandwidths (BW of the order of 100 Hz/pixel), the sampling time of an echo t s = 1/BW is of the order of 10 ms, thus allowing for the acquisition of one more signal at T E * 8 ms. The rationale behind the acquisition of more echoes should not be found in their availability without inducing major changes that may alter the original GRE signal at long T E (in terms of contrast, SNR, etc.). Indeed, multi-echo GRE allows for R Ã 2 quantification [16], along with a more accurate estimate of Quantitative Susceptibility Mapping (QSM) compared to what obtained with the single-echo acquisition [17]. These quantitative parameters provide a robust and comprehensive characterization of the tissue susceptibility, which may shed light in many areas of clinical interest (from mineral metabolism [18] to vein segmentation [19,20]).
On the other hand, a more thorough approach to quantitative MRI (qMRI) is desirable in a wide spectrum of clinical conditions [21], since R 1 and Proton Density (PD) give complementary and more detailed information on water tissue content. These parameters were shown to provide unique insights into the status of white matter in normal aging and demyelinating disease [22][23][24] or into tumor characterization [25,26]. In [27], PD was used to derive Macromolecular Tissue Volume Fraction (MTVF) maps, which provide a sensitive measure of the disease status in multiple sclerosis patients, and, together with R 1 maps, allow for an estimation of the apparent volume of interacting water protons (strictly related to the properties of local physico-chemical environment). Moreover, R 1 , coupled to R Ã 2 , is a key to determine the stage of haemorrhages [28] or to obtain absolute measures of contrast agent concentrations within tissues [29,30]. In addition, the joint availability of R 1 , R Ã 2 and susceptibility measures in Fabry's disease has been recently exploited to redefine the pathogenetic mechanisms and the incidence of the pulvinar sign in that metabolic disorder [31]. R 1 and PD maps are usually derived from two or more spoiled GRE sequences acquired at fixed T R and variable Flip Angle (FA) [32]. Therefore, one may consider to append one more GRE sequence (e.g. with a lower FA) to that used for SWI, QSM and R Ã 2 mapping [16]. However, an unpleasant and straightforward consequence of R 1 mapping by variable FAs consists in the increased total acquisition time, which clearly doubles the original duration of the SWI protocol, and certainly limits its clinical applicability.
Here we present the RElaxometry and SUsceptibility Mapping Expedient (RESUME) that halves the duration of the supplementary acquisition, still providing the complete set of qMRI maps (R 1 , R Ã 2 , PD and QSM) along with the SWI processing.

Materials and methods
RESUME is a qMRI scheme consisting of a 3D acquisition protocol coupled with the analytical solution of the Bloch equations associated to the acquired MRI signals. The protocol is derived from a standard SWI acquisition, modified to sample at least two echoes, and by the addition of one more spoiled GRE sequence with halved repetition time.
First, the analytical solution of the qMRI problem is given for a single-echo spoiled GRE (with repetition time T R,0 , echo time T E,0 and FA θ 1 ) and a multi-echo spoiled GRE (with repetition time 2 Á T R,0 , echo times {T E, i } and FA θ 2 ). Finally, more details on the actual acquisition setup will be given.

Analytical solution of the qMRI equations
Spoiled GRE signal. The complex signal of a spoiled GRE sequence is where K is the complex coil sensitivity, M 0 is the equilibrium magnetization, E 1 exp(−T R Á R 1 ), γ n is the gyromagnetic ratio of the imaged nucleus, ΔB is the local magnetic field variation and ϕ 0 is the phase-shift induced by the RF-pulse [16]. R Ã 2 map. If s i represents the magnitude of the i-th echo of the multi-echo spoiled GRE, from Eq 1 it follows that for Therefore, a simple regression via weighted least squares leads to the following estimates of and where Eq 4 reduces to in case of double-echo acquisitions. R 1 map. Eqs 1, 4 and 5 allow the estimation of the signal magnitude at T E = 0 of the single-echo (S 1 ) and multi-echo (S 2 ) spoiled GREs as (s 0 being the magnitude of the single-echo GRE) and (see Eq 3).
In order to determine which of the solutions in Eq 15 is physically relevant, E þ 1;0 is recognized as Looking for the parametric domain in which the condition is satisfied, after some algebraic manipulation, it results that: Of note, the condition in Eq 20 explicitly depends on the unknown E 1,0 (see Fig 1b)). However, the graphics in Fig 1c) and 1d) show that for E 1,0 > e −1 (corresponding to the realistic condition in which T R,0 < T 1 ), for any 0 < θ 2 < π and 0 < θ 1 < 0.47 Á θ 2 , E þ 1;0 provides the 1;0 satisfy the physical constrain 0 < E 1 < 1. b) Region in the parameter space (cosθ 1 , cosθ 2 , E 1,0 ) in which solution E þ 1;0 applies. c) Region in the parameter space (θ 1 /θ 2 , θ 2 , E 1,0 ) in which solution E þ 1;0 applies. d) Family of curves parametrized by E 1,0 (ranging from 0 [black] to 1 [light gray] in step of 0.1) that pose an upper limit to θ 1 /θ 2 for E þ 1;0 to be acceptable. The red curve corresponds to E 1,0 = e −1 and entirely lies in the half-plane θ 1 /θ 2 > 0. 47 correct solution to Eq 14. Therefore, in the following it will be assumed that A thorough analysis of the FAs θ 1 and θ 2 that optimize the R 1 map SNR in the general case of a multi-echo GRE with multiple echoes and variable BW i is quite complex and, probably, not really insightful. However, in the specific case of a double-echo GRE acquired with T E,1 = T E,0 and BW 1 = BW 0 (which is the most realistic condition-see the section dedicated to acquisition), a useful numerical approach can be adopted. From BW 1 = BW 0 , it follows that S 1 and S 2 share the same noise power s 2 S . Therefore the variance of the R 1 estimate can be expressed as The optimal FAs " y 1 and " y 2 are thus obtained by minimizing s 2 R 1 for a given R 1 ; the numerical results are shown for a wide range of practical PD map. From Eq 11, PD can be estimated as QSM. There are plenty of techniques that can be used to derive the QSM from the complex signal in Eq 1 (for a recent overview, see [33]), each more or less recommended depending on the specific application. A highly accurate reconstruction scheme is given by the algorithm introduced in [34], which is here adopted to obtain the map shown in the Results from the phase images of the multi-echo GRE sequences.

MRI acquisition
MRI was performed using a 3 T scanner (Trio, Siemens Medical Systems, Erlangen, Germany) with a volume transmitter coil and an 8-channel head receiver coil. Sequences were acquired after having obtained an informed consent for experimentation with human subjects. The "Carlo Romano" ethics committee for biomedical activities of "Federico II" University of Naples (Italy) specifically approved the study and the written informed consent form, which was signed by the subject undergoing the MR scan.
The 3D RESUME protocol actually acquired for this study consists of: Both sequences shared the same geometry (pseudo-axial orientation along the anterior commissure-posterior commissure; Field of View (FOV) = 230 × 194 × 166 mm 3 ; voxel size ΔV = 0.65 × 0.65 × 1.3 mm 3 ; GRAPPA acceleration factor = 2) and were acquired with a BW of 190 Hz/pixel in a total acquisition time of 7' 9". For each echo, a magnitude/phase reconstruction was enabled, thus obtaining a total of 3 complex volume datasets.
According to the chosen FAs, the SNR of the resulting R 1 map can be estimated from Eq 22. This quality index can be made independent from the noise power within the acquired images by normalizing it to the SNR of the input GRE images (possibly after some denoising) computed in an average brain parenchyma (with a nominal T 1,par = 1 s). The result is shown in Fig 3 as a function of the expected T 1 .

Denoising
Before computing the R Ã 2 , R 1 and PD maps according to Eqs 4, 21 and 23, the acquired datasets were denoised following a revised version of the Multi-spectral Non-Local Means (MNLM) approach first described in [16].
Briefly, in this context, M different datasets are considered as a multi-component image The MNLM is, then, a class of endomorphisms of the image space, each identified by 3 parameters (κ, ρ and B), that act as follows: In Eq 24 B r ½x is the ball centered inx, whose radius r is defined such that Z κ * 10 2 Á ΔV provides a criterion to dynamically allocate a search window, r $  The denoising of the magnitude datasets (s 0 , s 1 and s 2 ) required by R Ã 2 , R 1 and PD calculation was performed by exploiting also the information contained within the QSM derived from the phase images. Therefore, a total of M = 4 real datasets were jointly processed in Eq 24.
Due to the high computational complexity of the above scheme, a multi-Graphics Processing Unit (GPU) implementation of the MNLM algorithm [36] was adapted to Eq 24 and exploited for fast image processing.

Field inhomogeneity correction
B AE 1 non-ideal profiles critically influence the accuracy of the R 1 and PD estimates. In particular, while R 1 just depends only on B þ 1 through the presence of θ 1,2 in Eq 15 (the receiver coil sensitivity B À 1 is eliminated in the ratio of signal intensities in Eq 12), PD is also affected by B À 1 via S 2 in Eq 23.
If a protocol for a fast measurement of B 1 profiles of slab-selective RF-pulses is available on the scanner, the actual B 1 maps can be directly used in Eqs 15 and 23 to obtain bias-free R 1 and PD maps.
In the results shown below, to get rid of such biases, the effective approach based on the information content of the derived maps proposed in [16] was adopted.

RESUME assessment
The reproducibility and the accuracy of the RESUME maps were estimated on a pool of 10 subjects (5 patients and 5 healthy controls-HCs-were recruited in order to test the method across a broad range of tissue parameters) by applying the following procedure on each of them (hereafter identified by the index j = 1, . . ., 10).
First, to obtain the confidence interval of the RESUME maps, each sequence of the protocol was acquired twice during the same session. Given the redundant set of complex datasets, an ensemble of 2 3 (cardinality of the Cartesian product of magnitude datasets s 0 , s 1 and s 2 ) × 2 2 (cardinality of the Cartesian product of phase datasets from the dual-echo GRE-used for QSM) choices of different complete RESUME protocols was produced. Thirty-two different estimates of the relaxometry and susceptibility maps were thus derived and used to estimate the overall reproducibility (accounting for image uncorrelated noise, variations in signal amplification of the MR scanner, tissue temperature fluctuations, etc.) via voxelwise mean (μ j ðxÞ) and standard deviation (σ j ðxÞ) maps.
Moreover, the RESUME R 1 , R Ã 2 and PD maps were compared with the corresponding maps derived on the same subject using another qMRI approach established in [16]. This approach is based on two or more multi-echo spoiled 3D GRE sequences, in which all acquisition parameters are kept constant, except the FA, which is varied in order to provide different T 1 -weights. The acquired signals are then conveniently combined in order to derive high SNR maps of R 1 , R Ã 2 , PD and magnetic susceptibility. On the other hand, the RESUME susceptibility map was compared to the QSM derived from an 8-echo GRE sequence (Repetition time: 28 ms; Echo Times: {T E, i = [3.38 + (i−1) Á 2.82] ms}; FA: 20˚; BW: 400 Hz/pixel), acquired with bipolar readout gradients (no flow-compensation) and the same geometry of the RESUME sequences.
The accuracy was finally evaluated via the difference magnitude (δ j ðxÞ) between RESUME maps and the maps used for comparison.
For evaluation purposes, a healthy control (HC1) was chosen as the common spatial reference for the pool of subjects. The reproducibility (fσ j ðxÞ j j ¼ 2; . . . ; 10g) and accuracy (fδ j ðxÞ j j ¼ 2; . . . ; 10g) maps obtained for each other subject were, then, spatially normalized to HC1 according to the deformation field given by the elastic registration of the dual GRE image acquired at T E,1 on the corresponding HC1 image. Furthermore, for a numerical assessment of the maps, the median values within the brain mask were extracted for each patient from fσ j ðxÞg and fδ j ðxÞg.
To compare the RESUME measures with the values found in the literature, an expert neuroradiologist manually drew a set of bilateral Regions of Interest (ROIs) on the R 1 map of HCs in the white matter, in the cortical grey matter, in the head of the caudate nucleus and in the putamen; the mean values of R 1 , R Ã 2 , PD and QSM were then extracted for further comparisons.
Full access to the obtained results is provided at https://figshare.com/s/6d467faa11579c7c0 c02.

Image processing
All data processing, except the QSM (which was obtained by the free tool described in [34]), was performed using an in-house developed library for Matlab (MATLAB 1 Release 2012b, The MathWorks, Inc., Natick, MA, USA), partly described in previous works [16, 35-37], on a commercial workstation (Intel 1 Core™ i7-3820 CPU @ 3.6 GHz; RAM 16 GB) equipped with 2 GPU boards (NVIDIA GeForce 1 GTX 690). The demonstrative application of qMRI for vein segmentation, and Oxygen Extraction Fraction (OEF) and MTVF mapping was respectively provided by the MAVEN algorithm [20] and the in-house implementations of the methods described in [27,38].

Results
Each RESUME processing took 10 minutes on the above described workstation.
The derived maps exhibit a uniform quality throughout the entire FOV, and provide adequate spatial and contrast resolution for clear identification of the brain structures at different axial levels (see Fig 4). Moreover, the resolution of the acquired datasets allows for satisfactory multi-planar reconstructions (see Fig 5).
Reproducibility and accuracy measures of RESUME were available both for visual inspection (see Fig 6) and in terms of summarizing numerical indices. In particular, in the pool of the analyzed subjects, the median values of fσ j ðxÞg within the brain mask were normalized to the range of the associated fμ j ðxÞg maps, and their mean values were (1. 38  The RESUME values measured in the ROIs manually drawn in different brain structures of the HCs largely overlap with the corresponding ranges derived from the reference method [16] in the same cohort of subjects or found in the literature for different groups of HCs (see Table 1).
The output of an extended quantitative processing entirely based on the RESUME protocol is shown in Fig 7, where R 1 , R Ã 2 , PD and QSM come along with the MAVEN vein segmentation, the OEF and MTVF maps.

Discussion
The search of quantitative and reproducible measures in MRI is of unquestionable interest for the scientific community. Indeed, qMRI methods directly focus on physical and objective  parameters that could both increase the sensitivity of clinical MR scans, thus providing unique information about brain pathophysiology, and simplify the data harmonization in longitudinal or multicentric studies.
However, despite the growing interest of the neuroimaging community for qMRI (which provides hints on longitudinal and transverse relaxation processes, spin density, magnetic susceptibility, etc.), a certain skepticism usually creeps up at the moment of designing an acquisition protocol that could include it, for three main reasons.
First, the total acquisition time of a qMRI protocol that yields adequate resolution and SNR maps may be disheartening: in [47,48], R 1 and R 2 maps (neglecting susceptibility phenomena and PD) were obtained in about 15 minutes, while the method described in [16] squeezed in the same duration a complete high-resolution qMRI scheme; conversely, other authors succeeded in lowering the acquisition time (in the order of 5 minutes), but they greatly sacrificed the resolution ( [49,50] stooped to a slice thickness of 5 mm).
Second, the sequences required by the great majority of the qMRI approaches are usually of poor radiological relevance as source of traditional MR image contrasts: in a clinical context, for instance, a T 2 -weighted Fast Spin Echo (FSE) is commonly preferred to the TrueFISP [16, 51, 52]; a Magnetization-Prepared RApid Gradient-Echo (MPRAGE) or T 1 -weighted SE are usually acquired instead of short T E -GRE; a PD-weighted FSE, if necessary, substitutes the low FA-GRE [42,48]; etc.
The derived maps share the same resolution of the SWI, thus being suitable for multi-planar reconstructions (see Fig 5) and for the study of thin structures such as the cortex, small gray matter nuclei and medullary veins (see Figs 4 and 7).
Despite the short acquisition time (overall protocol duration of * 7 minutes), the derived qMRI datasets exhibit a considerably high SNR even at a simple visual inspection. This is further confirmed by the confidence interval maps (see Fig 6), whose extent encompasses general reproducibility aspects besides the phenomena strictly connected to noise propagation.
A comparison with the maps obtained by an established state-of-the-art qMRI protocol [16] highlights an excellent accuracy of the RESUME results: indeed, the magnitude of the difference between the values estimated by the two methods falls within the limits fixed by the reproducibility issues. Moreover, the ranges of RESUME measures in different brain structures of HCs are consistent with those found in the literature for a variety of different qMRI approaches (see Table 1).
In this context, it is worth stressing the 3D nature of the acquired sequences that lead to the estimation of the relaxation parameters. Indeed, as most clearly pointed out in [32], relaxometry schemes based on 2D imaging sequences are prone to severe and barely rectifiable inaccuracies, due to the intrinsic intra-voxel FA variations of the slice-selective RF pulses. Furthermore, a comparison of Inversion Recovery (IR)-based T 1 -mapping schemes [41] showed that 2D multislice acquisitions, such as IR-FSE, lead to an underestimation of T 1 , due to several source of errors (such as flow, perfusion, through-plane motion, magnetization transfer caused by offresonance excitation, non-ideal slice profiles, etc.) that are much less relevant on 3D sequences. On the other hand, relaxometry from 3D sequences can easily include B 1 -inhomogeneity correction schemes, e.g. based on the acquisition of the B 1 -map, and this seems consistent with the substantial agreement of the longitudinal relaxation rates measured in several brain structures by RESUME and by multiple MPRAGEs acquired at different inversion times [41].
From a computational point of view, the analytical RESUME solutions for the relaxometry equations are particularly beneficial for the processing execution time: the negligible burden associated with the computation of Eqs 4, 21 and 23 allows for an inexpensive application of the iterative B 1 -inhomogeneity correction adopted, so that the denoising step remains the bottleneck of an acceptable 10-minutes pipeline.
Finally, the interest in acquiring a thorough qMRI protocol like RESUME may also go beyond the extraction of classical relaxometry or susceptibility maps. In fact, the reconstructed parameters allow for additional quantitative processing steps that rely on the primary RESUME maps, such as the vessel-based analyses leading to vein segmentation and OEF maps, or the assessment of macromolecule density and local physico-chemical environment (see Fig 7).
In conclusion, the proposed acquisition and processing scheme allows for an accurate, reproducible and time-affordable strategy for obtaining different quantitative measures from a brain MRI scan. These features make RESUME feasible not only for research aims, but also for clinical practice, in light of an always increasing incorporation of qMRI protocols in the neuroradiological routine.