Abstract
BACKGROUND AND PURPOSE: Clinical measurements of cerebral perfusion have been increasingly performed with multiecho dynamic susceptibility contrast–MR imaging techniques due to their ability to remove confounding T1 effects of contrast agent extravasation from perfusion quantification. However, to this point, the extra information provided by multiecho techniques has not been used to improve the process of estimating the arterial input function, which is critical to accurate perfusion quantification. The purpose of this study is to investigate methods by which multiecho DSC-MRI data can be used to automatically avoid voxels whose signal decreases to the level of noise when calculating the arterial input function.
MATERIALS AND METHODS: Here we compare postprocessing strategies for clinical multiecho DSC–MR imaging data to test whether arterial input function measures could be improved by automatically identifying and removing voxels exhibiting signal attenuation (truncation) artifacts.
RESULTS: In a clinical pediatric population, we found that the Pearson correlation coefficient between ΔR2* time-series calculated from each TE individually was a valuable criterion for automated estimation of the arterial input function, resulting in higher peak arterial input function values while maintaining smooth and reliable arterial input function shapes.
CONCLUSIONS: This work is the first to demonstrate that multiecho information may be useful in clinically important automatic arterial input function estimation because it can be used to improve automatic selection of voxels from which the arterial input function should be measured.
ABBREVIATIONS:
- AIF
- arterial input function
- ΔR2*
- change in effective transverse relaxation rate
- QM
- quality of merit
- R
- Pearson correlation coefficient
- RMSerror
- root-mean-square error
Dynamic susceptibility contrast–MR imaging is an established clinical technique for imaging perfusion characteristics of the brain. Clinically, there has been interest in methods for accurate perfusion quantification across pathologies, which has motivated the transition from single-echo to multiecho techniques.1,2 These techniques desensitize perfusion measurements from bias resulting from contrast agent extravasation as is common in clinical contexts such as high-grade brain tumors.3 Most interesting, along with the measurement of multiecho signals comes the possibility of using this additional information beyond leakage correction.
The problems associated with standard single-echo DSC–MR imaging techniques arise from the assumption that changes in the measured signal intensity solely reflect changes in T2* relaxation times associated with contrast agent passage through the vasculature. This assumption is not valid in pathologies characterized by a disruption of the blood-brain barrier, in which contrast agent will accumulate in the surrounding tissue affecting both T2* and T1 relaxation times significantly and simultaneously. Fundamentally, contrast agent leakage breaks down the assumed relationship between measured signal intensity and underlying contrast agent concentration because signal changes cannot be assumed to be dominated by changes in T2*. Clinically, this broken assumption means that single-echo DSC–MR imaging techniques lacking leakage correction yield unreliable estimates of tissue perfusion.3 Therefore, there has been significant interest in implementation of multiecho techniques for perfusion measurement because they are insensitive to the T1 effects of contrast agent extravasation.4⇓⇓⇓–8
Multiecho DSC–MR imaging data also provide new information beyond leakage correction. For example, this new information could be used to improve the estimation of the arterial input function (AIF). AIF estimation is necessary for quantitative determination of cerebral blood volume, cerebral blood flow, and mean transit time from DSC–MR imaging data, and accurate estimation is highly dependent on the choice of locations for this estimation.9⇓⇓–12 Therefore, clinical applications place a priority on automated methods for AIF estimation due to their ease of use and consistency across users.13 However, to date, there is no consensus on how to reliably extract AIFs from multiecho data or to automate this process.
The goal of this study was to investigate methods by which multiecho DSC–MR imaging data can be used to automatically avoid voxels whose signal decreases to the level of noise when calculating the AIF. Such voxels are difficult to automatically identify with single-echo acquisitions. We expect that automatically avoiding such voxels would result in improved estimates of the AIF and would thereby provide more accurate measures of perfusion that are clinically feasible due to their automated nature.
Materials and Methods
Data Acquisition
In the present study, 117 pediatric clinical dual-echo DSC–MR imaging perfusion datasets were acquired on a 3T clinical MR imaging scanner with a body transmit coil used for excitation and an 8-channel array coil used for signal reception. Data were retrospectively anonymized and analyzed in accordance with the study regulations specified by the local internal review board. In each patient, standard clinical doses of Gd-DTPA contrast agent (0.1 mmol/kg) were administered through manual bolus injection approximately 1 minute after the beginning of perfusion imaging, followed immediately by a saline injection. Images were continually acquired for 4 additional minutes following the injection, yielding 5 minutes of data per patient. All data were acquired with the following parameters: single-shot gradient-echo EPI, FOV = 224 × 224 mm, voxel dimensions = 2.94 × 3.02 × 4.5 mm, number of sections = 21, TE1/TE2/TR = 15.16 ± 0.28/45.4/1500 ms, sensitivity encoding factor = 2, flip angle = 90°.
Image Processing
For each patient, images from each TE were separately converted to maps of change in effective transverse relaxation rate (ΔR2*) according to Equation 1.14 Dual-echo-based ΔR2* maps were also estimated according to Equation 2.14
where Spre refers to signal measured before contrast injection, STE1 refers to the signal at the earlier TE, STE2 correspondingly refers to signal measured at the later TE, and ΔR2*,DE refers to the dual-echo based estimate of ΔR2*.
Estimation of Bolus Arrival Time
To account for variable injection times among patients, we estimated the bolus arrival time for each patient. In general, our strategy was to identify the earliest time point whose signal demonstrated a significantly greater rate of change compared with the entire time course. Specifically, the temporal derivative of the mean whole-brain signal at TE1 was estimated for each patient. The extreme studentized deviate15 was calculated for each time point to identify outliers (P < .001), with the expectation that the points along the drop in signal corresponding to bolus arrival would have significantly higher rates of change than other time points. Isolated outliers (ie, outliers lacking a corresponding outlier either preceding or following it in time) were disregarded due to the expectation that bolus passage would have a duration of >1.5 seconds, desensitizing this analysis to random noise. The earliest remaining time point identified was considered the point of bolus arrival.
Automated AIF Estimation
To identify the most suitable voxels for AIF determination, we implemented a system of ranking candidate voxels. All voxels within the brain were ranked according to various qualities of merit (QMs), with low ranks indicating that a particular voxel was more suitable for AIF use and high ranks indicating poorer suitability. A combined score for each voxel was then calculated as the sum of all QM ranks being considered, with low combined scores indicating better suitability for AIF determination. The 10 voxels with the lowest combined score were selected for each patient, and the AIF was calculated as the average dual-echo-based ΔR2* time course across all 10 voxels. Thus, the different QM combinations could be assessed because any combination of QMs could be used in calculating the combined score.
Five QMs were measured in total and were available for use in identifying voxels for AIF estimation. The first 2 represented well-established AIF characteristics, and the final 3 represented potential new methods of comparing TE1 with TE2 to improve the AIF estimation. The first 2 QMs considered across all cases were the slope of the rise (ie, the rate of increase) in ΔR2*,DE across 4.5 seconds following bolus onset (ie, “slope”), and peak height of the ΔR2*,DE time courses in units of milliseconds−1(ie, “height”). These QM metrics are routinely used for automated AIF detection with single-echo DSC–MR imaging data.10 The 3 additional QMs considered alongside the first 2 were the following: 1) the linear Pearson correlation coefficient between ΔR2* calculated from TE1 and TE2 across 70 seconds spanning 10 seconds before bolus onset to 60 seconds following bolus onset (R); 2) the difference in the peak height of ΔR2* calculated from TE1 and TE2 (peak height difference); and 3) the root-mean-square difference (RMSerror) between the ΔR2* curves calculated from TE1 and TE2 according to Equation 3, where n represents the number of imaging volumes acquired.
We considered 4 possible QM combinations for identification of voxels for AIF calculation. These combinations were the following: 1) slope and height; 2) slope, height, and correlation; 3) slope, height, and peak height difference; and 4) slope, height, and RMSerror. These combinations represented standard parameters commonly used in the automated analysis of single-echo data (ie, combination 1) and 3 different combinations that leveraged the new multiecho information.
Perfusion Estimation
With each AIF selection method described above, perfusion maps were calculated for each patient on the basis of the dual-echo-based ΔR2* time courses. The AIF was converted to units of contrast agent concentration, assuming a quadratic relationship for contrast agent inside a blood vessel16 and a transverse relaxivity of Gd-DTPA (r2*) at 3T of 87 mmol−1 ms−1.17 Standard estimation by using circular singular-value decomposition was used to calculate CBV, CBF (milliliters/100 grams/minute), and MTT for the entire brain.18,19
Results
All data were manually inspected for evidence of artifacts associated with manual contrast agent injections. Features used for exclusion were evidence of double bolus peaks (indicating delay between contrast and saline injections) and insufficient baseline periods (inaccurate bolus timing and thus insufficient baseline periods). Of 117 datasets collected clinically during the study period (2 years and 162 days), 93 were considered acceptable quality for inclusion in the study. Qualitative inspection of the remaining voxel time courses in and around the branches of the middle cerebral artery in these data showed evidence of signal saturation at later TEs and potential T1 contamination at earlier TEs. An example can be seen in Fig 1.
Examples of potential pitfalls when selecting voxels from which to measure the arterial input function. Region 1: region including voxels within the major arterial branches. Note the decreased estimates at TE2, resulting in a distorted dual-echo estimate of ΔR2*. Region 2: voxels immediately outside the major arterial branch showing enhancement, resulting in artifactually low estimates of ΔR2* at TE1. Taken together, these regions illustrate the need for dual-echo estimates to avoid T1 effects and the need to avoid voxels that saturate into the noise floor to ensure accurate estimation of ΔR2* from the dual-echo data.
In all patients, voxels across the brain were ranked with respect to each QM. To evaluate how each of the 3 new metrics (R, peak height difference, RMSerror) for AIF voxel selection varied, along with more well-characterized metrics (ΔR2* onset slope and peak height), we calculated the mean values for R, peak height difference, and RMSerror across voxels and subjects and plotted them as a function of ΔR2* slope and peak height rank (Fig 2). In voxels with the shallowest bolus onset slope and shortest ΔR2* peak height, we observed low correlation between the ΔR2*,TE1 and ΔR2*,TE2. This may be consistent with voxels more distant from major arteries that have a lower contrast-to-noise ratio.
Illustration of potential markers of signal saturation by using multiecho data. Each potential marker (R, peak height difference, RMSerror) is shown across 3 sample sections in a representative subject (left column). Note markers identifying major arterial branches that are known to have signal saturation effects. In addition, each saturation metric is shown as a function of both ΔR2* peak height and bolus onset slope (right column). Contours are drawn on smoothed versions of the underlying images. Note that each marker varies nonmonotonically along the diagonal, indicating that simply maximizing the ΔR2* peak height and rise slope leads to suboptimal voxel identification for use in AIF determination.
Moving along the diagonal from shallow slopes and low peak heights toward steeper slopes and higher peak heights in ΔR2*, we observed that the correlation between ΔR2*,TE1 and ΔR2*,TE2 increased, reaching a maximum, and then decreased (Fig 2, upper right). This finding indicates that in the voxels with the highest ΔR2* peak height and steepest ΔR2* bolus onset slope, there is a discrepancy between the ΔR2* estimates measured from TE1 and TE2 independently. This is consistent with voxels experiencing signal saturation to the level of noise at later TEs. Likewise, the difference between ΔR2* peak heights measured from TE1 and TE2 was low for voxels with the shallowest onset slope and shortest peaks and then increased as the onset slope and peak heights increased (Fig 2, center right). The root-mean-square error between ΔR2*,TE1 and ΔR2*,TE2 followed a pattern similar to that of the peak height difference, though values of RMSerror varied little across voxels with shallow onset slope and short peak height (Fig 2, lower right).
Figure 3 shows an example of the effect of AIF selection strategy on the resulting global AIF shapes. When we considered only traditional AIF characteristics like bolus peak height and rate of ΔR2* increase, ΔR2* peaks were notably heightened initially, followed by broadened and truncated peaks, suggesting inclusion of voxels whose signal decreases to the level of noise during bolus passage. When using multiecho information to avoid voxels potentially with truncation artifacts, we saw the highest AIF peaks when using the correlation between ΔR2*,TE1 and ΔR2*,TE2 as the additional criterion for ranking voxels. This result corresponds to values for CBF that are slightly lower than reported values in adults.20,21 AIF peak heights were shorter when using peak height difference between ΔR2*,TE1 and ΔR2*,TE2 as the metric for avoiding voxels whose signal saturates to the level of noise and were shorter still when using the mean square error between the estimates of ΔR2*.
A demonstration of the difference between methods for selecting arterial input function voxels in a patient demonstrating abnormal left frontal perfusion. The upper row shows the location of the 10 voxels selected for AIF estimation in red, according to the indicated temporal properties. A projection map of major vessel locations, generated from maps of ΔR2* peak slope, is underlaid for spatial reference. The second row shows AIF time courses. Note that in all cases, comparing TE1 with TE2 results in avoidance of voxels whose signals are truncated by the noise floor, though there are differences in actual AIF according to which comparison method is chosen. The third and fourth rows show perfusion maps resulting from the above AIF estimates.
Comparing AIF estimates across patients, we saw similar results with reduced peak heights when using peak height difference or RMSerror as opposed to the Pearson R as the method for quantifying truncation artifacts. Looking at the distribution of mean whole-brain CBF estimates across subjects, we saw lower CBF when using peak height difference quantifications compared with RMSerror and Pearson R. Likewise, we saw the lowest CBV estimates under that condition as well, likely reflecting the larger area under the curve of the estimated AIF, due to its elevated tail following bolus passage (Fig 4). This finding is consistent with the sample perfusion maps shown in Fig 3.
The upper row shows AIF curves across subjects (mean ± SD). The middle row shows normal-appearing white matter mean CBF histograms across subjects The lower row shows normal-appearing white matter mean CBV histograms across subjects. White matter measurements are made in individually drawn, normal-appearing white matter ROIs. PHD indicates peak height difference.
Discussion
Here we present a new approach to automate AIF determination for multiecho DSC–MR imaging–based perfusion studies. The information contained in multiecho data enables improved estimation of the arterial input function though removal of T1 effects and identification and avoidance of voxels that are likely to be dominated by noise at later TEs. Data recorded from separate TEs can be independently quantified and compared against each other to identify voxels whose data are corrupted at either TE by artifacts. The multiecho-based AIFs exhibit features that are more consistent with the expected AIF shape and temporal characteristics and yield perfusion estimates in agreement with those in prior studies.
The advantages of the methods described here are not tied solely to the specifics of the automated AIF selection method. The “conventional” method for identifying AIF voxels, which is based on prior automated techniques,13 focuses on prioritizing the steepness of onset slope and ΔR2* peak height, though additional features could also have been included. For example, metrics like the bolus peak width and first moment of the concentration time curve may also be useful as AIF voxel-selection criteria, though all these criteria have important caveats that must be considered.22 In fact, the comparative benefits of incorporating multiecho information into the process of AIF voxel selection may decrease relative to conventional methods as more selection criteria are incorporated. Nevertheless, the specific choice of conventional AIF criteria is partially independent of the broader concept that multiecho information can be leveraged to provide information about physiologic signals that are not available when analyzing single-echo signals alone. For example, identifying slightly truncated ΔR2* peaks in single-echo data is difficult and, in some cases, may not be possible. The method of using multiecho information to avoid truncation artifacts can be implemented across a wide variety of AIF calculation strategies.
An interesting trend in multiecho DSC–MR imaging studies has been to estimate the AIF by using ΔR2*,TE1 alone, attempting to avoid the complication of identifying and avoiding voxels with truncation artifacts at the later TE. However, several disadvantages to this approach may be underappreciated. While it is less common for signals at early TEs to decrease to the level of noise, they may still be attenuated, resulting in noisier and more error-prone measures of ΔR2* during bolus passage. Furthermore, signals at early TEs are potentially confounded by T1 effects as illustrated in Fig 1. These effects theoretically underlie signals at later TEs as well, though they are more likely to be overlooked due to signals being dominated by T2* effects that result in a more characteristic signal decrease following bolus passage. These effects only become realized after the T1 effects of contrast agent passage are removed.
Advances in acquisition strategy such as multiband excitations23 and improved non-Cartesian k-space imaging have resulted in acquisitions with shorter TRs and TEs than have typically been available. Reductions in the TR may increase T1-weighting in the data similar to, though not as extreme as, the effect studied in the context of 3D PRESTO (Principles of Echo-Shifting with a Train of Observations) acquisitions.24 In the case in which section-wise acceleration results in reduced TR, multiecho techniques for AIF determination and tissue perfusion estimation gain increased importance. Likewise, as alternate acquisition strategies are used, in particular alternate k-space trajectories, very short TEs may become possible for the initial echo of multiecho acquisitions. This possibility too will alter the relative T1- and T2*-weighting of the resulting images, lending increased importance to multiecho approaches. In the extreme case in which the weightings of echoes in multiecho acquisitions become significantly different from each other, the approach of measuring the temporal similarity of one echo with the other may need to be reconsidered.
Because multiecho information is used to better select AIF voxels through avoiding truncation artifacts, the overall location of AIF voxels changes. This change may have implications for inclusion of partial-volume averaging artifacts. As Fig 3 illustrates, when one uses only bolus onset slope and peak height as the criteria for identifying AIF voxels, voxels are mainly identified along the M1 and M2 segments of the MCA. However, when one uses any of the methods for avoiding truncation artifacts, selected voxels increasingly lie along the posterior arteries, more distal points along the M2 MCA segments, and along smaller arterioles. While this outcome is the result of the advantageous avoidance of truncation artifacts, there is the potential of introducing partial volume artifacts. The chosen voxels may have reduced ΔR2* peaks due to voxels containing a mixture of signal from supplying arteries and local microcirculation. While other methods exist for addressing the issues of partial volume averaging,25 here our use of the peak height and bolus onset slope of ΔR2* as additional selection criteria rank arterial voxels more favorably. If one looked simply at ΔR2* peak height, Figs 3 and 4 suggest that using the Pearson correlation coefficient between ΔR2*,TE1 and ΔR2*,TE2 is the best choice as the metric for quantifying (and avoiding) truncation artifacts while minimizing partial volume artifacts.
Conclusions
The results shown herein indicate that multiecho information can be used to improve AIF estimation in DSC–MR imaging data through avoidance of truncation artifacts and removal of T1 effects. These represent another example of the larger concept that multiecho data can be used to provide a more robust measure of perfusion in clinical settings.
Footnotes
Disclosures: Ashley M. Stokes—RELATED: Grant: National Institutes of Health National Cancer Institute (1R01CA158079).* Christopher C. Quarles—RELATED: Grant: National Institutes of Health National Cancer Institute (1R01CA158079).* *Money paid to the institution.
This work was partially supported by the National Institutes of Health grant R01-CA158079.
Paper previously presented in part at: Annual Meeting of the International Society of Magnetic Resonance in Medicine, April 20–26, 2013; Salt Lake City, Utah.
Indicates open access to non-subscribers at www.ajnr.org
REFERENCES
- Received September 24, 2015.
- Accepted after revision December 14, 2015.
- © 2016 by American Journal of Neuroradiology
