Noise dependency in vascular parameters from combined gradient-echo and spin-echo DSC MRI

Dynamic susceptibility contrast (DSC) imaging is a widely used technique for assessment of cerebral blood volume (CBV). With combined gradient-echo and spin-echo DSC techniques, measures of the underlying vessel size and vessel architecture can be obtained from the vessel size index (VSI) and vortex area, respectively. However, how noise, and specifically the contrast-to-noise ratio (CNR), affect the estimations of these parameters has largely been overlooked. In order to address this issue, we have performed simulations to generate DSC signals with varying levels of CNR, defined by the peak of relaxation rate curve divided by the standard deviation of the baseline. Moreover, DSC data from 59 brain cancer patients were acquired at two different 3 T-scanners (N = 29 and N = 30, respectively), where CNR and relative parameter maps were obtained. Our simulations showed that the measured parameters were affected by CNR in different ways, where low CNR led to overestimations of CBV and underestimations of VSI and vortex area. In addition, a higher noise-sensitivity was found in vortex area than in CBV and VSI. Results from clinical data were consistent with simulations, and indicated that CNR < 4 gives highly unreliable measurements. Moreover, we have shown that the distribution of values in the tumour regions could change considerably when voxels with CNR below a given cut off are excluded when generating the relative parameter maps. The widespread use of CBV and attractive potential of VSI and vortex area, makes the noise-sensitivity of these parameters found in our study relevant for further use and development of the DSC imaging technique. Our results suggest that the CNR has considerable impact on the measured parameters, with the potential to affect the clinical interpretation of DSC-MRI, and should therefore be taken into account in the clinical decision-making process.


Introduction
Dynamic susceptibility contrast (DSC) magnetic resonance imaging (MRI) is widely used to assess cerebral perfusion in both clinical settings and research studies, and has proven clinically useful for various brain pathologies, including brain tumors, stroke and neurodegenerative disorders [1][2][3].
DSC MRI is usually performed using either gradient-echo (GE) or spin-echo (SE) acquisitions, and measures of cerebral blood volume (CBV), cerebral blood flow (CBF) and mean transit time (MTT) are typically obtained from the relaxation rate curves. Simultaneous acquisitions of GE and SE enable calculation of the relaxation rate curves R2 GE (t) and R2 SE (t) from a single contrast agent injection. The different sensitivity to magnetic susceptibility between GE and SE read-outs leads to differences in R2 GE (t) and R2 SE (t), enabling the estimation of vascular parameters beyond conventional DSC. This includes the vessel size index (VSI), a measure of the weighted mean of the vessel sizes [4,5], and the vortex area, a parameter that has been associated with the underlying vessel architecture of the tissue [6][7][8]. Several studies suggest that these parameters may hold great clinical value regarding stroke [9], tumor grading [10] and treatment response [7,[11][12][13].
However, the information-and thus the clinical value-provided by these parameters is limited by the inherent noise in the signal. The presence of noise in DSC signals has previously been investigated to determine optimal experimental parameters and the effect on conventional DSC analysis and deconvolution [14][15][16][17]. However, how noise affects the estimated VSI and vortex area have not yet been examined. Moreover, a simple approach to determine the quality of clinical DSC data based on image noise is still lacking.
Evaluating noise in the DSC data can be done in several ways. The signal-to-noise ratio (SNR) of the baseline signal have been used by several studies, defined by the mean divided by the standard deviation of the baseline signal [17][18][19]. However, this definition of SNR does not account for the bolus peak and can lead to a high SNR value if the baseline signal is high, although the corresponding R2-curve is noisy. Other and more advanced definitions of SNR have been proposed [14,15], but these definitions are adapted to the specific aims of the studies and become impractical when it comes to a general noise evaluation of clinical DSC data. Here, we define the contrast-to-noise ratio (CNR) as the peak of R2-curve divided by the standard deviation of the R2-baseline, a simple parameter to extract from clinical DSC without any complex analysis or any knowledge about the experimental setup.
In this study, we have evaluated the CNR in clinical GE and SE DSC from 59 patients with brain metastases (N = 29) and glioblastoma (N = 30) scanned at two different 3 T scanners, and performed simulations to investigate the effect of typical CNR levels found in clinical data. Moreover, relative vascular parameters are often considered in clinical DSC, where the parameters are normalized to a reference value in normal-appearing brain tissue. We have therefore also examined the effect of CNR on the relative parameters in tumor regions of the brain cancer patients.

Simulations
Monte Carlo simulations were performed to obtain signal responses from GE and SE DSC acquisitions as previously described [8,20]. DSC signal responses were calculated based on the dephasing of proton spins in a vascular network consisting of arterioles, capillaries and venules with specified vascular characteristics. The vascular network was modelled by a vessel tree structure where each vessel branches into two smaller vessels on the arterial side and two vessels merge into one on the venous side. The radius of the vessel generations decreases with a factor of 2 1/3 , starting with an arteriole with R = 10.0 µm, via the smallest capillary with R = 3.2 µm, to the last venule with R = 12.6 µm, resulting in a vessel tree with 4, 3 and 5 generations of arterioles, capillaries and venules, respectively [6,21]. For each vessel generation, a large number of protons (N = 10 000) with a given water diffusion constant, D (= 800 µm 2 s −1 ) diffuse in a simulation space filled with randomly distributed vessels with a given oxygen saturation level, SO 2 (100% for arterioles, 50% for capillaries and venules), blood volume fraction, vf, (2%) and radius, R. The accumulated phase dispersion of each proton is calculated based on the magnetic field perturbation at the proton's location and the resulting signal intensity is given by: where k indicates the kth vessel generation, φ n (t) is the accumulated phase of the nth proton and N is the number of protons. The phase dispersion is calculated for both GE and SE acquisitions and for varying concentrations of a Gadolinium-based contrast agent. The contrast agent was assumed to reside within the vessels, i.e. leakage effects were not considered. The kth vessel generation is associated with a concentration-time curve, C k (t), given by: where * denotes convolution, C k = 0 (t) is the arterial input function and h(t) is the transport function of transit times within the vessels. The MTT of the whole vessel network is given by the sum of MTTs of the individual vessel generations and was set to 4.5 s. The MTT was assumed constant across the vessel generations. The signal response from the vessel network is then obtained by a weighted sum of the signal response from each vessel generation: where w k is the fraction of the total blood volume made up from vessel generation k.

Noise and signal processing
Gaussian noise was added to the complex simulated signal (SI(t)) and sampled with a temporal resolution of 1.5 s (typical temporal resolution of DSC). The relaxation rate curve, R2(t), was calculated from the modulus of the signal with added noise: where TE is the echo time (100 ms for SE and 30 ms for GE), and SI 0 is the mean signal intensity of the baseline (pre-contrast) signal. Here and in the following, R2(t) denotes the relaxation rate curves from both GE and SE acquisition, if not explicitly stated otherwise. The CNR was calculated as: where R2 max denotes the maximum value of the R2(t) and σ baseline denotes the standard deviation of the R2 baseline. The simulated CNR ranged from 1 to 50, and for VSI and vortex area, the ratio of CNR in SE and GE data were set to match the median ratio found in clinical data, where GE CNR ranged from 13 to 106 when SE CNR ranged from 1 to 50.
The CBV values were estimated from the area under the curve (figure 1(a)): where the subscript XE denotes either GE or SE. The VSI was estimated by [4,5]: where CBV GE is the measured CBV from GE acquisition and D is the water diffusion coefficient given by the input to the simulation model. The ratio between GE and SE ∆R2 was found by the slope of a linear regression of the ∆R2(t) GE vs. ∆R2 (t) 3/2 SE data-points ( figure 1(b)).
The vortex area was given by the area encircled by the ascending and descending branch of the ∆R2(t) GE vs. ∆R2 (t) 3/2 SE data-points (figure 1(b)) and divided by the measured CBV GE to correct for the blood-volume dependency [6,7].
Each CNR level were repeated N rep = 1000 times, and the mean and standard deviation of the parameters was calculated relative to the measurement from the corresponding noise-free curves with high temporal resolution. The variability (%) was defined as the ratio standard deviation/mean of the measured values.

MRI acquisitions
The study was approved by our Institutional Review Board and the Regional Committee for Medical and Health Research Ethics (ref: 2013/1033). Informed consent was obtained from all patients. MRI data was obtained from 29 patients with brain metastasis from non-small cell lung cancer (N = 24) or malignant melanoma (N = 5) (17 females, 12 males, median age = 63 years, range 42-78 years; clinicaltrials.gov identifier: NCT03458455) [22]. In addition, clinical MRI data was acquired from pre-surgical exams of 30 patients later found to have glioblastomas by histopathological analysis (6 females, 24 males, median age = 59 years, range = 37-77 years) [8].

Image analysis
For brain metastasis patients, probability maps of white matter (WM) were obtained from T1 or FLAIR-images using the segmentation tool in SPM12 (Statistical Parametric Mapping (SPM) toolbox version 12, University College London, England). Binary masks of WM tissue were obtained by including voxels with probability >0.85 (maximum probability = 1.0) and were co-registered to the DSC space using normalized mutual information co-registration using nordicICE (NordicNeuroLab AS, Norway). The tumor region of brain metastases was delineated on contrast enhanced T1-weighted images by a radiologist and subsequently co-registered to DSC space. For glioblastoma patients, probability maps of normal-appearing WM were obtained from the Look-Locker acquisition as previously described [8,24]. The enhancing tumor region was automatically identified from pre-and post-contrast T1-weighted and FLAIR images, using a deep learning technique based on Convolutional Neural Network classifiers [25]. The deep learning technique was only available for the data from the glioblastoma patients, and was therefore not used on the brain metastases data.
Voxel-wise GE and SE CNR, CBV GE , CBV SE , VSI and vortex area were calculated as described above, using the motion-corrected and leakage-corrected DSC data. Leakage correction was performed using the method proposed by Boxerman et al, where both T1 and T2 effects are taken into account [26]. Voxel-wise ADC from diffusion images and CBV GE were used to calculated the VSI, and the vortex area were divided by CBV GE to correct for the CBV-dependency [7]. The vascular parameters were subsequently normalized (denoted rCBV, rVSI, rVortex area) to a reference value given by the median value in normal-appearing WM. The reference tissue was limited to the WM region in order to obtain a homogenous reference tissue value and to adhere to the recent consensus recommendations [23].
In the clinical data, the variability of the relative parameters was calculated by the standard deviation divided by the mean of the subject-wise median, where the subject-wise median for varying levels of CNR was calculated by selecting voxels with CNR > n and CNR ⩽ n + 1, where n ranged from 0 to 30.
To investigate the effect of CNR on the relative values in tumour regions, the reference value was first calculated by the median value including all voxels in WM (as above). Then, the reference value was calculated based on the WM, but excluding voxels with CNR below given cut-off value ranging from 2 to 8. The relative change vs. cut-off values was then calculated for the resulting median value in tumor region: where X cut-off and X all voxels denote the resulting value when using a cut-off value or all voxels, respectively.

Results
In figure 2, typical CNR maps and corresponding histograms of the CNR values from a GE acquisition ( figure 2(a)) and a SE acquisition ( figure 2(b)) are shown for a patient with glioblastoma, where the GE data generally display higher CNR values compared to the SE data. Moreover, voxel-wise curves of ∆R2(t) are displayed for varying levels of CNR, ranging from CNR = 2 to CNR = 50 ( figure 2(c)). Across all patients, median (interquartile range) CNR in WM were 21.3 (14.1-27.3) for GE acquisitions, and 6.4 (4.6-6.9) for SE acquisitions. For tumour regions, median CNR were 40.8 (20.7-62.8) and 8.9 (4.4-11.6) for GE and SE acquisitions, respectively. Figure 3 shows boxplots of median CNR in WM and tumour for the two patient groups, where tumour regions display a considerably larger interquartile range of CNR compared to the WM. No significant difference in CNR was found in between the patient groups, except for a higher CNR in WM for brain metastasis patients compared to glioblastoma patients (27.2 vs. 17.4, p < 0.05, Mann-Whitney U test).
The simulations showed that, for all parameters, the estimations were highly unreliable for CNR < 4, with large deviations from the ground truth (i.e. mean values deviating from 1) and large standard deviations ( figure 4). For the CBV measurements from both GE and SE acquisitions, a slight overestimation (<15%) could be seen for CNR ≈ 4-10, (figures 4(a) and (b)). For vessel size measurements, a distinct underestimation of rVSI was seen for CNR < 10, with an underestimation of 60% for CNR = 4 (figure 4(c)). A general underestimation was also seen in the measured vortex area, with an underestimation of approx. 20% for CNR ≈ 5-10 ( figure 4 (d)).
The variability declined for increasing CNR for all parameters in the simulated data ( figure 5(a)). The variability of the rCBV measurements was around 80% for CNR = 4 and decreased to under 20% for CNR > 15. The variability of rVSI was lower than for the other parameters, with a variability around 30% at CNR = 4 and under 10% for CNR > 15. The vortex area displayed a higher sensitivity to noise, with large variability (>40%) for all CNR levels up to CNR = 15.
The variability of the measured parameters in clinical data displayed a similar trend as in simulated data, with a general decrease in variability for increasing CNR ( figure 5(b)). In the clinical data, a rapid decrease in the variability could be seen for up to CNR = 4, before a stabilization for higher CNR. The variability of rCBV GE and rVortex area stabilized around 50%, while the variability of rCBV SE and rVSI stabilized around 30%. Larger fluctuations in the variability curves can be seen for higher CNR due to a decrease in the number for voxels for increasing CNR, which is the underlying cause for the spike seen in the variability for rVSI at CNR close to 30. In clinical data, we saw that excluding voxels with low CNR affected the reference value (median value in WM) and the corresponding relative values in the tumour regions depending on the applied cut-off. In figure 6 the relative change in median tumour values vs. cut-off values are shown for all four parameters. The relative change is a result of the combined effects of the change in the reference value and the change in the distribution of values within the tumor region, but where the former had a larger impact on our data. When increasing the cut-off value, we saw a general decline in the median rCBV from both GE and SE acquisitions. The decline was generally larger for SE acquisitions compared with GE (figures 6(a) and (b)). For rVSI, a general increase in the relative values was seen for higher cut-off values (figure 6(c)), while rVortex area showed a decline for higher cut-off values ( figure 6(d)). For a cut-off value at CNR = 4, the relative change ranged within [−5%, 4%] for rCBV GE , [−18%, 8%] for rCBV SE , [−25%, 20%] for rVSI, and [−45%, 47%] for rVortex area. Moreover, the range of the relative change (error bars in the figure), generally increased for higher cut-off values. This is because the relative change is calculated based on the value for cut-off = 0 (no voxels are excluded). As the cut-off increases, more voxels are excluded, and the relative change (included the range across subjects) increases. The median number of voxels excluded for each cut-off value is displayed in table 1. Figure 2. CNR map from a GE acquisition (a) and a SE acquisition (b), and the corresponding CNR values displayed in normalized histograms (c), (d). CNR from GE acquisitions are generally higher compared to CNR from SE acquisitions. Moreover, high CNR signal typically coincides with major vessels with higher CBV (a), (b). Typical relaxation rate curves are shown for CNR ranging from 2 to 50 (e). Abbreviations: CNR-contrast-to-noise ratio; GE-gradient-echo; SE-spin-echo.
In figure 7, an example of how excluding voxels with low CNR (< 4) can affect the relative values in the tumour region is shown for rCBV SE and rVSI. For this glioblastoma patient, applying the cut-off value at CNR = 4 led to a decrease in the rCBV SE and an increase in the rVSI in the tumour region.

Discussion
In this study, we shed light on how noise affects the estimation of the vascular parameters CBV, VSI and vortex area from DSC imaging. Our results show that the parameters display different sensitivity to CNR, and that low CNR can potentially introduce large errors in the interpretation of clinical DSC data.
In our clinical data, the CNR was about a threefold higher in GE acquisition than in SE acquisitions. This is due to a higher response to the contrast agent bolus in GE acquisition [27]. Moreover, we saw that tumour regions generally displayed a higher CNR than the WM region, which is likely caused by the inherent CBV-dependency of the definition of CNR. High CNR may stem from high R2 peak and/or low standard deviation of the baseline. Tumours often have abnormally high CBV, causing a high R2 peak and thus a high measured CNR [28]. No difference in CNR between the patient groups was found, except for a higher CNR Figure 3. Boxplots of median CNR in white matter (yellow boxes) and tumour regions (red boxes) for GE acquisitions (a) and SE acquisitions (b). No significant difference in CNR was found in between the patient groups, except for a higher CNR in white matter for brain metastasis patients compared to glioblastoma patients (p < 0.05, Mann-Whitney U test). Abbreviations: CNR-contrast-to-noise ratio; GE-gradient-echo; SE-spin-echo; GBM-glioblastoma; MET-metastasis. in WM for brain metastases compared to glioblastomas in GE acquisitions. This may be due to the different scanners or contrast agent used.
In the simulated data, a general decrease in the measured VSI was seen for low CNR. This is likely due to the lower CNR in SE data compared with GE data, which results in a more pronounced overestimation of the SE data, and thus an underestimation of the slope value used in the VSI measurements. This is in line with a previous study that reported underestimated VSI in enlarged vessels [29]. For the vortex area, we saw a larger sensitivity to simulated noise compared with the CBV and VSI measurements. This is probably owing to how Figure 5. The variability of the measured parameters from simulated data (a) and clinical data (b). High variability was seen for CNR < 4 across all parameters, indicating unreliable measures of the parameters for low CNR. Abbreviations: CNR-contrast-to-noise ratio; GE-gradient-echo; SE-spin-echo; CBV-cerebral blood volume; VSI-vessel size index. Figure 6. Boxplots of the relative change in the median value in the tumour region of brain metastases (dark red) and glioblastomas (light red) for rCBVGE (a), rCBVSE (b), rVSI (c) and rVortex area (d). The cut-off value indicates the threshold for CNR applied on the maps (i.e. excluding voxels with CNR < cut-off value) before normalizing the parameters to the white matter reference value. The change is calculated relative to the median value for no cut-off, and the number sign next to the box indicates p < 0.05 (Wilcoxon signed rank test). Abbreviations: CNR-contrast-to-noise ratio; GE-gradient-echo; SE-spin-echo; CBV-cerebral blood volume; VSI-vessel size index; MET-brain metastases; GBM-glioblastomas. the vortex area is estimated, where each time-point on the vessel vortex can affect the resulted area considerably and is thus very sensitive to outliers caused by noise. In contrast, VSI is determined by the slope of a linear regression of the data-points and is less sensitive to outliers. In addition, VSI is affected by the noise in GE CBV and ADC as well, whereas the vortex area is affected by the noise in GE CBV. However, our simulations show that the variability in VSI and vortex area mainly originate from the variability in the Table 1. The median number of voxels excluded in the WM and tumour regions for each CNR cut-off. The corresponding percentage of the total ROI is given in the parenthesis.

Number of voxels excluded in WM regions
Number of voxels excluded in tumour regions estimation of the slope value and area, respectively. For CBV measurements, the variability in our simulated data was of the same order as results from previous studies [14,17,18], but a direct comparison could not be made due to the different study designs and definitions. The variability in both simulated and clinical data was high (>50%) for CNR < 4, which indicates that measurements from DSC data with CNR < 4 will be highly unreliable. In clinical data, the variability displayed a steeper decline compared with simulated data before a stabilization occurred. This may be due to the different way variability is measured in clinical vs. simulated data. In simulated data, the variability is given directly by the standard deviation/mean of the estimated parameters, whereas in clinical data, the variability is given by the standard deviation/mean of the median values across the patients. A direct comparison between the variability in simulated and clinical data should therefore be performed with care, yet, the result gives an indication of the general CNR dependency of the variability.
The CNR of the clinical data is measured from the leakage corrected curves. In the normal-appearing WM, the leakage correction is not likely to affect the result due to the absence of leaky vessels. In the tumor regions however, the leakage correction may influence both the measured CNR and the measured parameters. In addition, other leakage correction techniques may cause different effects on the resulting CNR. This should be investigated further, accompanied by simulations where different leakage effects and leakage corrections approaches are considered. This also includes assessing scenarios with use of a preload contrast agent administration.
In clinical data, we have demonstrated how excluding voxels with low CNR can cause changes in the relative median values in tumour regions. In our data, a higher CNR cut-off resulted in a general decrease in the median rCBV and rVortex area, and an increase in the rVSI of the tumour regions. These trends are likely caused by the CBV-dependency of the CNR. In clinical data, low CBV values are associated with low CNR values, and thus excluding low CNR values will exclude the low CBV values as well. The reference value in WM will therefore increase, causing a general decrease in the relative CBV values in the tumour regions. Correspondingly, the CBV-dependency of VSI measurements and vortex area, will cause an increase and a decrease in the relative values, respectively, when excluding low CNR voxels. Although a measure of CNR that is not dependent on CBV (or other underlying vascular properties) would be desirable, this is difficult to obtain in practice. To obtain information from the DSC data, a bolus passage, and thus a CBV > 0, is necessary.
In addition, in our clinical data, the relative change of the median tumour values was mainly due to a change in the reference value and not due to exclusion of voxels within the tumour region. However, this varied across subjects and the subject-wise distribution of CNR within the reference tissue and the tumour region.
Our study findings indicate that the voxel-wise CNR may have pronounced impact on the estimated vascular parameters, and the influence of CNR should therefore be taken into consideration when performing DSC analyses. By applying a reasonable cut-off level for CNR, the noisiest image voxels will be excluded and is thus likely to give more robust parameter values. However, the chosen cut-off must be weighed against the probability of discarding relevant information from the images. For example, excluding low CNR values will also exclude low CBV values, where the latter may hold clinically valuable information. Moreover, the cut-off level must also be weighed against the number of voxels present in the tissue region, as the number of voxels will decrease with increasing cut-off value. This is especially relevant for regions of interest with a low number of voxels, e.g. in small tumours. In our data, the combined results suggest that voxels should have at least CNR = 4, and the higher CNR, the more reliable the estimated parameter will be.
There is a large body of work in the literature addressing influence of image noise on the performance of different deconvolution techniques for estimating CBF and MTT [17,19,30,31]. Because estimation of VSI metrics is not dependent on CBF/MTT estimations, a more careful analysis on topic is assumed outside the scope of the current work. However, although a noise threshold usually is applied at some point in the signal processing of DSC imaging in both research and clinical settings, the noise level itself is seldom reported. We believe that measuring and reporting the noise, e.g. in terms of CNR, should be a standard practice when using these techniques to give more insight to how the noise varies across centres, scanners, protocols and subjects, and also to give an indication of the quality of the data.
Our study has limitations. The simulated data in our study is designed to give an indication of the CNR-dependency of clinical DSC data. However, a direct comparison between simulations and clinical data is challenging due to several factors. First, in the simulated data, the underlying vascular characteristics are kept constant and the variation in the measured values can be attributed directly to the added noise. In clinical data, however, the vascular characteristics may vary from voxel to voxel and a heterogenous distribution of tissue properties is unavoidable. Second, the measured CBV or VSI in clinical MRI can vary even though the underlying blood volume or vessel size is constant. This can be due to differences in the biophysical properties of the tissue (e.g. vessel size, vessel distribution, diffusion coefficient etc [20,32,33].), or instabilities of the MRI scanner system that does not change the observed CNR. The relation between the parameter and the underlying vascular properties is even more intricate regarding the vortex areas, as this parameter is not linked to one main tissue property, as CBV or VSI, but has been shown to be affected by several vascular properties, including blood volume fractions, MTT and the branching structure of the vessels [6][7][8]. Furthermore, in simulations, uncorrelated Gaussian noise was added to the complex MRI signal. However, in clinical data, noise from different sources can be spatially and/or temporally correlated and is not necessarily reflected in the CNR as defined here. In addition, leakage effects were not considered in the simulations. In clinical data, leakage corrections were performed, but may influence both the measured CNR and the parameters. Finally, a range of tumor related scenarios could produce low CBV and thereby low CNR, including radio-damaged tissue (pseudoprogression), necrotic regions, oedematous tissue, or tumors regions with low or no angiogenic activity. Careful assessment of CNR in these regions are therefore warranted. Nevertheless, while the results from our simulations should be interpreted with care, the findings can serve as an indication as to how CNR affects the measurement in clinical DSC imaging.
In conclusion, our simulations show that noise has different effects on the vascular parameters from combined GE-SE DSC imaging. We found that CBV measurements become overestimated when CNR is low, while VSI and vortex area measurements become underestimated. Moreover, the vortex area generally displayed a larger sensitivity to noise than the other parameters. Our clinical data suggest that CNR < 4 gives highly unreliable measurements, and the distribution of relative values in tumour regions can vary considerably depending on the cut-off value used for CNR.