A Fully Automated Method for Quantitative Cerebral Hemodynamic Analysis Using DSC–MRI

Abstract

Dynamic susceptibility contrast (DSC)-based perfusion analysis from MR images has become an established method for analysis of cerebral blood volume (CBV) in glioma patients. To date, little emphasis has, however, been placed on quantitative perfusion analysis of these patients, mainly due to the associated increased technical complexity and lack of sufficient stability in a clinical setting. The aim of our study was to develop a fully automated analysis framework for quantitative DSC-based perfusion analysis. The method presented here generates quantitative hemodynamic maps without user interaction, combined with automatic segmentation of normal-appearing cerebral tissue. Validation of 101 patients with confirmed glioma after surgery gave mean values for CBF, CBV, and MTT, extracted automatically from normal-appearing whole-brain white and gray matter, in good agreement with literature values. The measured age- and gender-related variations in the same parameters were also in agreement with those in the literature. Several established analysis methods were compared and the resulting perfusion metrics depended significantly on method and parameter choice. In conclusion, we present an accurate, fast, and automatic quantitative perfusion analysis method where all analysis steps are based on raw DSC data only.

Introduction

Dynamic susceptibility contrast (DSC) perfusion magnetic resonance (MR) imaging (MRI) is a relatively new technique to characterize brain tumors and in particular gliomas. Analysis of relative cerebral blood volume (rCBV) is used to differentiate low-grade from high-grade gliomas (Covarrubias et al, 2004) and to assess patient prognosis and tumor progression (Law et al, 2008). Although some investigators have reported quantitative analysis of tumor blood volume and perfusion values using DSC-based MRI (Law et al, 2006), quantitative perfusion analysis of this patient group has not received as much attention to date. The limited number of quantitative studies may reflect a combination of the increased complexity and user dependence associated with quantitative analysis methods, combined with an unclear need for such methods due to the already high diagnostic accuracy of qualitative methods.

An increasing need for quantitative perfusion analysis techniques may, however, become apparent with the move toward using MR-based perfusion analysis in longitudinal analysis of treatment response in primary brain tumor patients (Gerstner et al, 2008). Absolute CBV and cerebral blood flow (CBF) derived from DSC imaging have been monitored in tumoral, peri-tumoral, and normal-appearing brain tissue before, during, and after radiation treatment, and some preliminary studies already indicate that the tracer mean transit time (MTT) may provide important information on treatment response (Nagesh et al, 2007). Other studies have shown that early CBV changes during radiotherapy in gliomas predicted different physiological responses to treatment, and this information may help identify tumor sub-volumes that are radioresistant and might benefit from intensified radiation (Fuss et al, 2000). Given the abnormal properties of the neovasculature associated with malignant tumors, showing a combination of high vascular density and disorganized vessel architecture, it would be expected that the MTT would be altered in malignant tissue. Further, quantitative perfusion analysis may enable the monitoring of baseline perfusion in unaffected brain tissue, allowing inter- and intra-patient comparisons across MR machines and institutions.

It is therefore believed that DSC imaging-based quantitative perfusion analysis of patients with primary brain tumors may become increasingly important. Such analysis does, however, pose several methodological challenges related to, for example, identification of arterial input function (AIF), AIF partial volume (PV) effects (PVEs), non-linear dose response, and choice of deconvolution method. Although the technical literature addressing these issues is extensive, little has been done to test the clinical utility of quantitative DSC-based perfusion imaging in large patient materials. Clinical implementation of quantitative methods adds further challenges related to user friendliness, user independence, robustness, and processing speed. The aim of this work was thus to develop and validate an end-to-end automatic method for quantitative perfusion analysis using DSC images for routine clinical use. A wide range of deconvolution techniques and methodological variations have been proposed in the literature, and hence several deconvolution methods and models were tested. As a first validation step, the focus of the analysis in our study was perfusion measurements in normal-appearing white and gray matter in patients with confirmed glioma. Hence, the methodology presented includes automated algorithms for extraction of normal-appearing brain tissue using only raw DSC data, thereby eliminating the need for image co-registration and processing of multiple data sets.

Materials and methods

Patient Selection

Study approval was obtained from the Regional Medical Ethics Committee and patients were included only if informed consent was provided. A total of 101 patients (51 males, mean age 51 years, range 8 to 79 years) with confirmed glioma, selected from an ongoing prospective tumor perfusion study, were included in the analysis. All patients received a glioma diagnosis after MR perfusion imaging and subsequent surgery (resection or biopsy). All imaging was performed before surgery and any subsequent treatment.

MR Imaging

Imaging was performed at 1.5T (Siemens Sonata, Symphony, or Avanto; Siemens AG, Erlangen, Germany), using an eight-channel (Symphony/Sonata) or a 12-channel (Avanto) head-coil. DSC–MRI was performed using a gradient-echo, echo-planar imaging sequence acquired during contrast agent (CA) administration. The imaging parameters were the following: TR/TE 1,430 ms/46 ms, bandwidth 1,345 Hz/pixel (12 axial slices) or 1,720/48, bandwidth 1,500 Hz/pixel (14 axial slices), field of view 230 × 230 mm, voxel size 1.8 × 1.8 × 5 mm³, and inter-slice gap 1.5 mm. For each slice, 70 images were recorded at intervals equal to the repetition time. After approximately eight time points, 0.2 mmol/kg of the CA gadobutrol (Bayer Schering Pharma AG, Berlin, Germany) was injected at a rate of 5 mL/sec, immediately followed by a 20-mL bolus of saline (B Braun Melsungen AG, Melsungen, Germany) also at 5 mL/sec.

Image Analysis

AIF detection:

The AIF was automatically identified in each slice using K-means cluster analysis (Mouridsen et al, 2006a). The K-means cluster implementation was based on a method developed for automatic segmentation of brain tissue sub-classes in diffusion-weighted images (Hadjiprocopis et al, 2005), but for our purpose applied to DSC concentration–time curves (CTCs) rather than pixel intensities (Emblem et al, 2009). In short, CTCs were assumed to belong to one of five tissue sub-classes (white matter, gray matter, arterial blood, venous blood, or ‘other’), and multiple features were used to assign each CTC into the most probable tissue class. The following CTC features were used: first moment (FM) of the CTC time integral, peak height (P_max), time to peak height (TTP), and bolus arrival time (T_d).

Since the final cluster allocation for each CTC is based on an iterative procedure with random starting points, there is no guarantee that the final solution is unique and optimal. The cluster analysis procedure was therefore repeated 50 times and the final AIF for each image slice was selected from the lowest FM cluster (determined for each run) with the highest average P_max. The average of five CTCs in this cluster with the lowest FM was then selected as the final AIF

AIF PV correction:

Since arterial diameters are small compared with the voxel size used in DSC–MRI, the measured AIFs are attenuated and distorted by PVEs (Kjolby et al, 2009; van Osch et al, 2005). Given that the dose-response is the same in venous and arterial voxels, and that venous voxels without PVEs can be identified, the PV correction factor can be estimated from the area under the curve of the respective arterial and venous CTCs (Ziegelitz et al, 2009). However, in DSC acquisitions the first-pass T2* effect in voxels containing only blood is generally so large that complete signal loss occurs, with resulting venous CTC distortion (Knutsson et al, 2007). Therefore, the first-pass phase was not used to correct for PVEs. Instead, the pseudo-steady-state phase after the first pass was used under the assumption that the venous CA concentration is sufficiently low to avoid complete signal loss. Assuming that the dose-response is similar in arteries and veins in the steady-state phase, and that the measured voxel-wise change in T2* relaxation rate is proportional to CA concentration, the PV correction factor, f_pv, can be determined as:

$\begin{array}{c}{C}_{\mathrm{a}\mathrm{\_}\mathrm{s}\mathrm{s}}(t)=f_{\mathrm{p}\mathrm{v}}\cdot C_{\mathrm{v}\mathrm{\_}\mathrm{s}\mathrm{s}}(t)\end{array}$

(1)

where C_{a_ss}(t) and C_{v_ss}(t) are the steady-state concentrations in the identified arterial and venous voxels, respectively. It was assumed that the concentration levels in both the arteries and veins followed similar decay curves (or remained approximately constant) so that the mean concentration levels were used and the PV correction factor was then determined as the ratio of the mean steady-state concentrations:

$\begin{array}{c}{f}_{\mathrm{p}\mathrm{v}}=\frac{\overline{C_{\mathrm{a}\mathrm{\_}\mathrm{s}\mathrm{s}}}}{C_{\mathrm{v}\mathrm{\_}\mathrm{s}\mathrm{s}}}\mathrm{a}\mathrm{n}\mathrm{d}{C}_{\mathrm{a}}(t)=\frac{1}{f_{\mathrm{p}\mathrm{v}}}{C}_{\mathrm{a}}^{\mathrm{P}\mathrm{V}}(t)\end{array}$

(2)

where C_a^PV is the CA concentration in the AIF in the presence of PVEs.

Venous output functions (VOFs) were identified using the same cluster analysis approach as described above, but using a large CTC FM as the selection criterion. The steady-state levels in the AIFs and VOFs were estimated by first fitting a γ-variate function to the respective CTCs. The CTC tails were then identified as the portion of the curves with amplitude above the corresponding value of the γ-variate-fitted curve (starting at the last time-point), as shown in Figure 1B. A Fourier-based, low-pass filter was applied to the identified AIF and VOF before identification of the PV factor to improve the stability of the method. A single PV correction factor was estimated for all slices based on the average steady-state concentration levels obtained from the identified arterial and venous CTCs. For each subject, the validity of the PV correction factor was assessed (from the criterion that 0 < f_PV ≤ 1) and any failure to detect at least one pixel representing a valid AIF in each slice was reported during the analysis.

Figure 1

Estimation of PV correction factor. (A) The AIFs automatically detected in each slice and (B) the corresponding average AIF across all slices together with the average VOF, which was also based on VOFs determined individually in each slice (not shown). The corresponding γ-variate-fitted curves are shown as dotted lines and the intersection points between the AIF/VOF and the γ-variate-fitted curves were used to determine the steady-state portion of the respective curves as indicated. Finally, the PV correction factor was determined from the difference in the steady-state levels of the AIF and VOF, respectively, as indicated by the arrow.

Quantitative perfusion analysis:

Quantitative analysis of tissue perfusion and blood volume was performed using established tracer kinetic models (Ostergaard et al, 1996; Wu et al, 2003) expressing the relationship between the tissue response and the AIF in terms of the convolution integral:

$\begin{array}{c}{C}_t(t-t_{\mathrm{d}})=\frac{K_{\mathrm{h}}}{\rho }F\cdot C_{\mathrm{a}}(t)\otimes R(t-T_{\mathrm{d}})\end{array}$

(3)

where ⊗ denotes convolution, C_t(t) is the CTC in tissue, C_a(t) is the PV corrected arterial CTC, F is a scaling factor proportional to CBF, K_h is a factor correcting for the difference in hematocrit in large vessels and capillaries (K_h = 0.71), and ρ is tissue density (ρ = 1.04 g/mL) (Rempp et al, 1994). R(t) is the tissue residue function describing the fraction of the instantaneous bolus of tracer present in the tissue at time t. T_d is the delay in the tracer arrival so that:

$\begin{array}{c}R(t-T_{\mathrm{d}})=\left\{\begin{array}{ll}0& t<T_{\mathrm{d}}\\ 1& t=T_{\mathrm{d}}\\ 0& t\to \mathrm{\infty }\end{array}\right.\end{array}$

(4)

Equation (4) assumes that the AIF is not significantly dispersed (Calamante, 2005) and T_d can then be estimated from the time-shift of the peak amplitude of R(t), relative to t = 0. CBF is then given by:

$\begin{array}{c}\mathrm{C}\mathrm{B}\mathrm{F}=\frac{K_{\mathrm{h}}}{\rho }F\end{array}$

(5)

Two different methods to estimate CBV an MTT (henceforth referred to as method-I and method-II) were compared. In method-I, CBV was estimated from the assumption that fractional blood volume in a voxel scales with the area under the CTC, normalized to the area under the curve of the AIF (Ostergaard et al, 1996):

$\begin{array}{c}{\mathrm{C}\mathrm{B}\mathrm{V}}^{\mathrm{I}}=\frac{K_{\mathrm{h}}}{\rho }\frac{\int C_{\mathrm{t}}(t)\mathrm{d}t}{\int C_{\mathrm{a}}(t)\mathrm{d}t}\end{array}$

(6)

From the central volume principle (Lassen and Perl, 1979) the corresponding MTT is then given as:

$\begin{array}{c}{\mathrm{M}\mathrm{T}\mathrm{T}}^{\mathrm{I}}=\frac{{\mathrm{C}\mathrm{B}\mathrm{V}}^{\mathrm{I}}}{\mathrm{C}\mathrm{B}\mathrm{F}}\end{array}$

(7)

In method-II, MTT was determined independently from the area under the residue function (Vonken et al, 1999) as follows:

$\begin{array}{c}{\mathrm{M}\mathrm{T}\mathrm{T}}^{\mathrm{I}\mathrm{I}}=\int R(t)\mathrm{d}t\end{array}$

(8)

and CBV is then given by the central volume principle as follows:

$\begin{array}{c}{\mathrm{C}\mathrm{B}\mathrm{V}}^{\mathrm{I}\mathrm{I}}=\frac{K_{\mathrm{h}}}{\rho }F\int R(t)\mathrm{d}t\end{array}$

(9)

Equation (3) was solved to obtain CBF, MTT^II, CBV^II, and T_d using four different singular value decomposition (SVD)-based methods:

• Method (a): A standard truncated SVD (sSVD), as originally proposed by Ostergaard et al (1996).

• Method (b): A block-circulant SVD (cSVD) extension of sSVD, which has been shown to be less sensitive to tracer delays (Wu et al, 2003).

• Method (c): Oscillation-limited cSVD (oSVD), which is an iterative method repeating the cSVD deconvolution process until the oscillation in the residue function is below a threshold defined by the oscillation index (OI) (Wu et al, 2003).

• Method (d): Iterative Tikhonov regularization SVD method (iTrSVD) using the L-curve criterion to determine the optimal regularization parameter (Hansen and Oleary, 1993).

For methods (a) and (b) the singular value threshold (P_SVD) was varied between 0.025 and 0.4 (fraction of maximum singular value) and for method (c) the OI was varied between 0.1 and 0.9. The range of both P_SVD and OI was chosen to cover the entire range of values reported to have been used in previous studies. The deconvolution analysis was performed on the raw CTC data without any prior spatial or temporal smoothing. The signal–noise ratio (SNR) was calculated from the baseline portion of each DSC data set as SNR = SI₀/σ, where SI₀ is the average signal intensity in baseline and σ is the SI₀ standard deviation (s.d.).

Tissue clustering:

The gray and white matter were automatically segmented using cluster analysis based on the temporal characteristics of the CTCs. Before this classification, areas of edema secondary to brain tumor were identified by a simple one-feature cluster analysis eliminating the class with the highest signal intensity in the baseline (pre-contrast) images. This approach assumed that the baseline DSC images were sufficiently T2-weighted to render pathological tissue hyper-intense in the baseline images. Once areas of abnormal appearance were removed, the remaining CTCs were classified using K-means cluster analysis with three classes and a single feature (area under the CTC). The cluster with the highest area under the curve was assumed to represent blood, the middle cluster gray matter, and the lowest cluster white matter.

Non-brain tissue voxels were removed before tumor extraction and tissue classification using a two-step procedure. First, the noise level in the baseline DSC images was determined automatically using the method of Otsu (1979). A seed point was then automatically selected centrally in the brain and a seed-growing algorithm was used to select all connecting voxels above the noise level. Due to the relatively poor resolution of the DSC images, the orbital structures as well as the skull appeared disconnected from the brain and were thus consistently removed with this approach.

All image-processing methods were written in C+ + (CodeGear C+ + Builder; Embarcadero Technologies, San Francisco, CA, USA) and integrated in a modified version of the software package nordicICE (NordicImagingLab AS, Bergen, Norway).

Statistical Analysis

Possible correlation between the perfusion metrics and patient age was tested using stepwise linear regression analysis and analysis of variance, with tumor grade and gender included as possible confounding factors. Regression analysis was performed for the whole patient population as well as separately for males and females and the regression coefficient, β, was determined. Gender differences in the perfusion metrics were tested using two-sample t-tests with assumed equal variance, and differences between deconvolution methods were tested using paired t-tests. Both age and gender effects were tested for all deconvolution methods. Significance level was set to α = 0.05 for all tests.

Since identification of both AIF and VOF (and hence PV correction factor), as well as the gray/white-matter tissue classes, was based on iterative cluster methods with random initial seed points, stability analysis was performed to assess the uncertainty in the resulting perfusion values introduced by the cluster analysis. Ten subjects were selected randomly and the complete perfusion analysis was repeated 50 times in each of these subjects. Mean and s.d. values for CBF and CBV^II for the combined gray- and white-matter tissue classes were obtained using cSVD deconvolution, with a fixed value for P_SVD. The mean relative uncertainty in CBF and CBV^II was then calculated as the average of the individual s.d.s normalized to the mean value for each perfusion parameter for each of the 10 subjects. This test was performed both with and without PV correction to assess any additional uncertainty introduced by the automatic PV correction method.

Statistical analysis was performed using Minitab15 (Minitab, State College, PA, USA).

Results

Of the 101 cases investigated, there were six WHO grade-I pilocytic astrocytomas, 19 grade-II diffuse astrocytomas, 5 grade-II oligodendrogliomas, 6 grade-II oligoastrocytomas, 3 grade-II ganglioglioma, 6 grade-III anaplastic astrocytomas, 2 grade-III anaplastic oligodendrogliomas, 5 grade-III anaplastic oligoastrocytomas, and 49 grade-IV glioblastomas. There was no difference in mean age between the gender groups (mean age 51.8 years in males and 50.6 years in females; P = 0.74).

AIF/VOF Detection and Tissue Cluster Analysis

For all patients, an AIF and a valid PV correction factor could be identified in all image slices using the automated method. Figure 1A shows the resulting AIFs for each slice in a sample subject with glioblastoma (male, age 40 years). Figure 1B shows the corresponding mean AIF and the mean VOF averaged across all slices. The average (± s.d.) correction factor was f_PV = 0.30 ± 0.06 (range 0.13 to 0.49). Figure 2 summarizes the tissue clustering results from the same patient as in Figure 1. The average (± s.d.) number of voxels across all patients classified as gray and white matter was 7,800 ± 1,877 and 13,167 ± 1,986, respectively, and the corresponding gray/white-matter volume ratio was 0.59 ± 0.12. The total analysis time for one subject, including automatic AIF/VOF determination, tissue clustering, and generation of all perfusion maps on a standard desktop PC, was about 15, 20, 80, and 90 secs when applying sSVD, cSVD, oSVD, and iTrSVD, respectively.

Figure 2

A sample case showing segmentation of tumor and classification of non-affected white/gray matter in a patient with glioblastoma (40-year-old male). Identification of gray- and white-matter tissue classes and blood was based on cluster analysis of the relative CBV (rCBV) values, resulting in the three tissue segments shown in panel D as overlays on the rCBV map shown separately in panel A. The three cluster classes are shown in red (high CBV), green (medium CBV), and blue (low CBV). The tissue cluster with the highest intensity in the baseline DSC scans (B) were excluded from the identified tissue classes and, as shown, the excluded region corresponds well with the area affected by the tumor according to the hyper-intense region in the T2-weighted image (C). The medium-CBV class corresponded well with the expected gray-matter regions and the low-CBV class with white-matter regions. Panel E shows the corresponding histogram distributions of the rCBV values (using the same color code as in panel D) for all slices for the respective tissue classes.

Dependence of Perfusion Metrics on Deconvolution Method

CBF in both gray and white matter decreased monotonously in a bi-linear manner with increasing P_SVD and cSVD deconvolution gave lower CBF values than sSVD deconvolution at equal P_SVD for all values of P_SVD, as shown in Figure 3A. Using oSVD deconvolution, CBF increased monotonously with increasing OI threshold as shown in Figure 3B. The mean CBF obtained with iTrSVD deconvolution is shown for comparison in both figures. The mean combined white- and gray-matter CBF estimated with the four different SVD methods converged at around 50 mL/100 g/min and this perfusion value was obtained using P_SVD = 0.1, P_SVD = 0.2, and OI = 0.06 for cSVD, sSVD, and oSVD deconvolution, respectively. Table 1 summarizes the mean values of CBF, CBV, MTT, and T_d obtained with the four different deconvolution methods and two kinetic models using the regularization parameters specified above.

Figure 3

Dependence of CBF estimates on SVD regularization parameter. The figures show the change in estimated CBF for the combined gray- and white-matter tissue classes as function of singular value cut-off (P_SVD, as fraction of maximum singular value) for (A) sSVD and cSVD, and (B) as a function of the OI used for oSVD. In both figures, the mean CBF obtained with iterative iTrSVD is shown for comparison (dotted line). Note that all SVD methods converge to the same CBF value at about 50 mL/100 g/min.

Table 1

Comparison of perfusion metrics obtained for 101 subjects with different deconvolution methods and kinetic models

Method	Tissue	CBF	CBVII	MTTII	T d
sSVD (0.2)	WM	30.41 ± 9.22	2.04 ± 0.62	4.32 ± 0.80	0.89 ± 0.91
	GM	65.13 ± 17.32	3.88 ± 0.91	3.85 ± 0.79	0.29 ± 0.30
	GM/WM	2.16 ± 0.17	1.92 ± 0.19	0.89 ± 0.07	0.31 ± 0.14
cSVD (0.1)	WM	30.15 ± 8.18	1.98 ± 0.51	4.07 ± 0.88	0.86 ± 0.99
	GM	58.13 ± 15.75	3.49 ± 0.84	3.78 ± 0.93	0.32 ± 0.31
	GM/WM	1.93 ± 0.10	1.79 ± 0.19	0.93 ± 0.09	0.93 ± 0.09
iTrSVD	WM	26.36 ± 11.33	1.91 ± 0.54	5.00 ± 1.09	1.37 ± 0.81
	GM	64.31 ± 25.70	3.68 ± 0.92	4.17 ± 1.10	0.77 ± 0.47
	GM/WM	2.47 ± 0.25	1.96 ± 0.22	0.83 ± 0.08	0.54 ± 0.16
oSVD (0.06)	WM	27.10 ± 11.34	1.89 ± 0.65	4.60 ± 1.21	0.96 ± 1.38
	GM	56.64 ± 19.77	3.42 ± 1.04	3.90 ± 1.06	0.40 ± 0.65
	GM/WM	2.13 ± 0.26	1.84 ± 0.22	0.85 ± 0.08	0.39 ± 0.19
		CBVI	CBVII	MTTI	MTTII
cSVD (0.1)	WM	2.30 ± 0.58 a	1.98 ± 0.51	4.81 ± 0.99 a	4.07 ± 0.88
	GM	4.57 ± 1.05 a	3.49 ± 0.84	4.99 ± 1.03 a	3.78 ± 0.93
	GM/WM	2.01 ± 0.17	1.79 ± 0.19	1.04 ± 0.07	0.93 ± 0.09

Abbreviations: CBF, cerebral blood flow; CBV, cerebral blood volume; cSVD, block-circulant SVD; GM, gray matter; iTrSVD, iterative Tikhonov regularization SVD method; MTT, mean transit time; oSVD, oscillation-limited cSVD; sSVD, standard truncated SVD; SVD, singular value decomposition; T_d, bolus arrival time; WM, white matter.

CBF (± s.d., mL/100 g/min), CBV (mL/100 g), and MTT (secs) obtained with the different deconvolution methods: sSVD (0.2) with singular value threshold, P_SVD = 0.2; cSVD (0.1) with P_SVD = 0.1; iTrSVD; oSVD (0.06), with oscillation index, OI = 0.06. The superscripts (I) and (II) used for MTT and CBV denote which kinetic method (method-I or method-II) was used to estimate the parameter. T_d is tracer arrival delay in seconds.

a

Significant difference (paired t-test; P < 0.05) in the corresponding perfusion parameter obtained using method-I versus method-II.

Figure 4 shows the changes in CBF, CBV^II, and MTT^II and T_d as a function of P_SVD using cSVD deconvolution. Similar trends were observed using sSVD. CBV^II decreased slightly with increasing P_SVD, but not to the same extent as for CBF, whereas MTT^II increased almost linearly with increasing P_SVD. T_d decreased with increasing P_SVD, especially in white matter.

Figure 4

Correlation between perfusion metrics and applied singular-value threshold (P_SVD). Change in (A) CBF, (B) CBV, (C) MTT, and (D) tracer delay time (T_d) as a function of P_SVD for gray matter (GM) and white matter (WM) using cSVD.

CBV^I was significantly higher in both gray matter and white matter as compared with the corresponding values of CBV^II, irrespective of the deconvolution method used to estimate CBV^II. The corresponding values of MTT^I were also significantly higher both in gray and white matter as compared with the same values obtained for MTT^II.

Correlation of Perfusion Metrics with Age, Gender, and Tumor Grade

The correlation between CBF and patient age was generally higher using cSVD as compared with sSVD and the highest correlation was obtained using cSVD (0.1). A significant correlation between CBF and age was only observed using cSVD or sSVD deconvolution and not using iTrSVD or oSVD (for any OI value). The age and gender correlations obtained with cSVD (0.1) deconvolution are summarized in Table 2 and the linear regressions of age versus CBF, CBV^II, and MTT^II are shown in Figure 5. CBF decreased with age and MTT^II increased with age in both gray and white matter. No correlation was found between CBV^II or CBV^I and age or gender in any tissue type. CBF was significantly higher, and MTT^II significantly shorter, in females as compared with that in males both in both gray and white matter. The correlation between age and CBF/MTT^II was generally more pronounced in females than in males. No correlation was found between tumor grade and any of the perfusion metrics for normal-appearing gray or white matter. CBF was higher—and MTT was shorter—in females than in males for all deconvolution methods except for oSVD with OI < 0.06.

Figure 5

Dependence of perfusion metrics on age. Linear regression showing correlation between age and CBF (top row), CBV (middle low), and MTT (bottom row) for white matter (left column) and gray matter (right column). Significant age dependence was found for CBF and MTT in both WM and GM, but not for CBV, in any tissue type. The highest age correlation (as shown in this figure) was observed using cSVD with a singular value threshold (P_SVD) of 0.1.

Table 2

Age and gender dependence of CBF, CBV, and MTT

Tissue	Gender	CBF	CBVII	MTTII
Mean (± s.d.)
WM	Female (n = 50)	33.52 a ± 8.10	2.06 ± 0.54	3.77 a ± 0.79
	Male (n = 51)	26.84 ± 6.88	1.89 ± 0.48	4.36 ± 0.87
	All (n = 101)	30.15 ± 8.18	1.98 ± 0.51	4.07 ± 0.88
GM	Female	64.66 a ± 16.08	3.62 ± 0.90	3.49 a ± 0.89
	Male	51.72 ± 12.58	3.37 ± 0.77	4.06 ± 0.89
	All	58.13 ± 15.75	3.49 ± 0.84	3.78 ± 0.93
β (± s.e.)
WM	Female	−0.20 ± 0.06	0.00 ± 0.00	0.001 ± 0.006
	Male	−0.05 ± 0.06	0.00 ± 0.00	0.011 ± 0.008
	All	−0.14 ± 0.04	0.00 ± 0.00	0.010 ± 0.005
GM	Female	−0.43 ± 0.11	0.00 ± 0.00	0.023 ± 0.007
	Male	−0.17 ± 0.11	0.00 ± 0.00	0.017 ± 0.006
	All	−0.33 ± 0.06	0.00 ± 0.00	0.018 ± 0.005
P(H:β = 0)
WM	Female	0.001	0.10	0.143
	Male	0.394	0.89	0.180
	All	0.002	0.31	0.046
GM	Female	< 0.001	0.57	0.010
	Male	0.14	0.46	0.018
	All	< 0.001	0.76	0.001

Abbreviations: CBF, cerebral blood flow; CBV, cerebral blood volume; cSVD, block-circulant SVD; GM, gray matter; MTT, mean transit time; SVD, singular value decomposition; WM, white matter.

Estimated CBF (± s.d., mL/100 g/min), CBV (mL/100 g), and mean transit time, TTP (secs) in WM and GM as a function of age and gender. MTT and CBV were estimated using method-II. The values shown were obtained using cSVD with singular value threshold, P_SVD = 0.1; cSVD (0.1) gave the highest age correlation of all the regularization parameters and SVD methods investigated. The regression coefficient, β (± s.e.), was obtained from linear regression.

a

Significant difference (P < 0.05) compared with the same parameter in males (two-sample t-test; P < 0.05). The numbers in bold indicate regression coefficients, which are significantly different from null (ANOVA; P < 0.05). P(H:β = 0) is the probability that null hypothesis (β = 0) is true.

SNR and Stability Analysis

The average whole-brain SNR measured from the pre-contrast dynamic signal was 45.6 ± 14.5. SNR was higher in white matter (47.4 ± 15) than in gray matter (45.5 ± 15) (paired t-test; P < 0.001). The uncertainty in the perfusion estimates as measured by the stability analysis was 7.3% and 8.0%, respectively, for CBF and CBV^II when PV correction was applied. The uncertainty dropped to 3.6% for both CBF and CBV^II without PV correction.

Discussion

We present a fully automated, user independent tool for quantitative analysis of hemodynamic parameters using DSC images. The proposed methods provide global gray- and white-matter perfusion analysis, including slice-wise identification of AIFs, PV corrections, and tissue segmentation, computed in about 10 to 15 secs on a standard desktop computer without any user interaction. Contrary to other reported DSC-based methods (Shin et al, 2007), all analysis is based solely on the DSC images, eliminating the need for time-consuming MR protocols and complex image co-registration routines.

The average CBV and CBF values obtained from normal-appearing white and gray matter in the 101 glioma patients investigated are in good agreement with previously reported values for patients and normal subjects obtained using different MR-based techniques (Larsson et al, 2008; Petersen et al, 2009; Shin et al, 2007). The observed decrease in CBF with age and the gender differences in CBF and MTT are also in agreement with previous MRI and [¹⁵O]H₂O positron emission tomography studies (Borghammer et al, 2008; Leenders et al, 1990; Shin et al, 2007). Contrary to some previous studies, no age-related changes in CBV were observed. The absence of any change in CBV with age is consistent with the observed age-related increase in MTT, suggesting that decrease in CBF with age is caused by a true reduction in tissue blood flow and is not just an indirect effect of reduced CBV, which (from the central volume principle) would have left MTT unchanged.

In spite of the good agreement between our results and those of many previous studies, large discrepancies exist in the reported normal cerebral CBF and CBV values, as well as age and gender variations in these values, obtained using different MR techniques and also using different imaging modalities, and a gold standard seems hard to establish. The large deviation in the reported normal perfusion values probably reflects both true inter- and intra-patient physiological variations (Leenders et al, 1990), as well as methodological differences between the different imaging techniques used to measure cerebral perfusion. The methodological challenges are particularly apparent in MR-based perfusion analysis and an extensive literature exists in addressing the various problems relating to AIF selection (Bleeker et al, 2009; Calamante et al, 2004; Mouridsen et al, 2006a), PVEs in the AIF (Kjolby et al, 2009; van Osch et al, 2005), non-linear dose response (Calamante et al, 2007; Kiselev, 2005), and choice of optimal deconvolution technique and model parameters (Calamante, 2005; Calamante et al, 2003; Knutsson et al, 2004; Murase et al, 2001; Wu et al, 2003). Since the aim of this study was to establish a user-independent tool for robust automated quantitative perfusion analysis, several of the most commonly used deconvolution techniques were tested and compared.

Our results raise some important methodological issues related to quantitative DSC perfusion imaging. The significant dependence of CBF and MTT on the choice of SVD singular value threshold (P_SVD) as observed in our study (especially at low P_SVD values) is in good agreement with results from numerical simulations (Murase et al, 2001), and correct choice of P_SVD is clearly important. Although CBF varied significantly with the choice of regularization parameter values, we found that the average CBF values obtained with the different deconvolution methods converged at about 50 mL/100 g/min for the combined gray- and white-matter tissue classes. Interestingly, this CBF value was obtained using P_SVD = 0.1 (for cSVD), which also gave the highest correlation between CBF and patient age. At lower P_SVD values CBF varied significantly for small changes in P_SVD and the s.d. in the perfusion estimates also increased. The optimal choice of P_SVD is clearly dependent on the SNR in the raw data (Wu et al, 2003) and the average SNR in this study (acquired at three different 1.5-T systems) was about 45, with an inter-scan SNR variation of about ± 30%. A method has previously been proposed, which aims at finding the optimal P_SVD for each voxel based on the local noise characteristics (Liu et al, 1999). We did not test this approach but an iterative SVD method with no free parameters using Tikhonov regularization (iTrSVD) was investigated. This method was found to give perfusion values comparable to those obtained with standard truncated SVD methods (cSVD and sSVD), but with a tendency to give a larger gray matter/white-matter CBF ratio. This was a result of the iterative regularization algorithm tending to apply a higher degree of regularization in white matter as compared with gray matter (data not shown). Further, the iTrSVD method proved to be less sensitive in detecting age variations in CBF and

MTT, which may imply that the method generally applies too much SVD regularization in both tissue types. We used an iterative Tikhonov regularization algorithm that aims at minimizing the residual norm while keeping the residual function solution norm as small as possible (Hansen and Oleary, 1993). Others have expanded this method further, minimizing the first-order derivative of the solution norm, thereby minimizing oscillations in the resulting residue function estimate (Calamante et al, 2003; Larsson et al, 2008). This method was not tested in this work due to the anticipated resulting increase in the already long processing time of the iTrSVD method (90 secs) as compared with that of the standard non-iterative SVD techniques (15 to 20secs). Further, since oscillations in the residue function, R(t), will mainly be centered on the zero value (Mouridsen et al, 2006b), a minimization of the solution norm (in effect minimizing ∫R(t)²dt) will also effectively result in a solution with little oscillations. An alternative method to reduce oscillations in R(t) is to define an OI and repeat the deconvolution interactively until the oscillations in the resulting R(t) are below the defined limit (Wu et al, 2003). The resulting perfusion values obtained with this approach (oSVD) were also comparable to those obtained with standard SVD, and the resulting CBF depended significantly on the choice of OI, thereby introducing a new free parameter that needs to be determined and optimized for a given SNR profile. The expected age-dependent variations in CBF and MTT were not detected using oSVD with any of the OI thresholds investigated. It is therefore concluded that standard truncated SVD methods may be the most appropriate method for quantitative perfusion analysis of normal brain tissue, given that the optimal singular value threshold has been determined. This conclusion may not hold under conditions of highly abnormal flow patterns or in situations where the CA bolus is significantly delayed or dispersed, due to the inherent sensitivity of SVD-based deconvolution to differences in tracer arrival times between tissue voxels and voxels defining the AIF (Calamante et al, 2006). Different methods have been proposed to overcome this problem and two block-circulant methods (cSVD and oSVD) proposed by Wu et al (2003) were tested, although abnormal perfusion patterns were not expected in our data since care was taken to exclude abnormal-appearing tumor tissue from the tissue clusters included in the analysis. The cSVD method was found to give systematically lower CBF estimates at a given P_SVD as compared with standard SVD (sSVD), but the dependence of CBF on P_SVD was identical for the two methods. Further, the estimated tracer delays were not significantly different using block-circulant SVD (cSVD and oSVD) versus standard SVD (sSVD and iTrSVD). All deconvolution methods gave larger tracer delays in white matter as compared with that of gray matter, as previously reported (Christensen et al, 2008).

Two different kinetic models for deriving CBV and MTT values were compared. In the most commonly used approach (model-I), only CBF is determined by deconvolution and CBV is determined independently from the tissue/AIF CTC time integrals. MTT is then given by the central volume principle: MTT^I = CBV^I/CBF. In model-II both CBF and MTT are determined from the residue function after deconvolution and CBV is then the dependent parameter given by the product of CBF and MTT^II. For the two models to provide equal estimates of CBV and MTT (Equations (6) and (9)), the following criterion must be satisfied:

$\begin{array}{c}F\int R(t)\mathrm{d}t=\frac{\int C_{\mathrm{t}}(t)\mathrm{d}t}{\int C_{\mathrm{a}}(t)\mathrm{d}t}\end{array}$

(10)

Model-I was found to give significantly higher CBV values as compared with model-II for all deconvolution methods and P_SVD values, suggesting that Equation (10) does not hold true experimentally. This result is consistent with a previous study of stroke patients (Perkio et al, 2002) where it was concluded that method-I may lead to over-estimation of CBV due to recirculation effects in the CTCs, and further that this method is more sensitive to noise as compared with method-II since the SVD deconvolution process effectively filters out noise components by eliminating small singular values dominated by noise. Interestingly, the value of CBV^II was found to be almost independent of P_SVD. This can be explained by the fact that reduction in CBF with increasing P_SVD is offset by a corresponding increase in MTT^II so that CBV^II (equal to the product of CBF and MTT) is almost constant, as shown in Figure 4.

Since the aim of our work was to develop a quantitative perfusion method that can readily be implemented in a clinical setting, it was found highly desirable to base all analysis solely on raw DSC images to simplify the analysis and speed up processing. Automatic correction for PVEs in the AIF is a critical component in the processing pipeline, and the aim to only use raw DSC images prevented the use of more common PV correction schemes based on additional T1-weighted images (Wirestam et al, 2007) or quantitative T1 maps (Shin et al, 2007). The venous first-pass signal response is generally distorted in DSC images due to excessive T2* effects in intravascular voxels unaffected by PVEs, and complex, user-dependent correction regimes are needed to restore the true VOF (Knutsson et al, 2007). These saturation effects could in theory be reduced by reducing the echo time, but this would also negatively affect the dynamic signal change in parenchymal voxels, reducing the sensitivity of the method. The PV correction factor was therefore determined from the tail of the CTC in the identified AIF and a corresponding VOF in each slice. This was found to give a much more robust estimate of the PV correction factor as compared with using the first-pass venous response when implemented as an automated routine.

Different from most previously published methods, both AIF and VOF were identified in each slice to eliminate the need for automatic identification of a given vessel segment based on anatomical location, which would significantly complicate image processing. The slice-wise AIF determination may have advantages in that the AIF may better reflect the true input function for each slice, minimizing artificial delays and dispersion effects (Calamante et al, 2006). The disadvantage of this approach is that the AIF may not be well-defined in all slices, especially in the basal part of the brain. In spite of this potential limitation, the slice-wise AIF/VOF determination approach was found to be robust and reproducible, and providing physiologically plausible perfusion estimates with inter-patient variations comparable to the values obtained using different modalities and perfusion analysis methods.

We obtained an average PV correction factor of 0.3 in our study, with a variation across subjects of about 20% (in terms of the relative s.d.), suggesting that the method provides a reasonable estimate of PVEs across subjects. Assessing the overall inter-patient variation in the perfusion estimated, there is, however, clearly a significant variability, which cannot be explained by gender- or age-dependent effects. Some of this variability is likely due to imperfections in the PV estimation and AIF determination. One limitation with the proposed method is a lack of correction for possible shape distortion in the AIF due to, for example, non-linearity in the arterial dose-response. Such distortions, not observable from the steady-state signal, could lead to broadening of the residue function (and hence the MTT). Further, significant CA-induced T1-shortening effects could potentially occur in the steady-state phase, which may lead to under-estimation of the steady-state CA concentrations. This effect would, however, affect both arteries and veins, and the net effect on the AIF/VOF steady-state ratio is therefore expected to be small.

A limitation in the study design was that a healthy, age-matched control group was not included. Only patients with known primary brain tumor were investigated and although we implemented a method to automatically exclude pathological tissue from the analysis, there is no guarantee that a portion of the normal-appearing tissue masks may include tissue with abnormal hemodynamic properties. There was, however, no correlation between the CBV or CBF measured in the apparent unaffected tissue and tumor grade. Such a correlation would be expected if a significant portion of the tissue masks included tumor tissue, since high-grade gliomas are known to give significantly higher CBF and CBV values as compared with low-grade gliomas (Law et al, 2006). Other studies have also shown that normal-appearing regions of interest do not show altered CBV values compared with healthy controls in patients with primary brain tumors (Wenz et al, 1996). A further potential limitation in the methodology is that gray- and white-matter tissue clustering was based solely on the dynamic properties of the CTC, without any reference to structural images. Identification of the different tissue classes is commonly based on signal intensity differences in heavily T1-weighted images, but structure-based identification of gray matter would not account for long-range T2* effects occurring in large vessels in the DSC images, leading to significant over-estimation of CBF and CBV in gray matter. Large vessel domination could potentially have been reduced by using a spin echo-based EPI sequence rather than a gradient echo-based sequence, but this would have resulted in significant reduction in temporal resolution or brain coverage since longer echo times are needed in spin echo–EPI sequences to maintain sufficient CA sensitivity. Gray- and white-matter identification based on the dynamic CTC information was shown to reliably identify two tissue classes unaffected by large intravascular T2* effects, as shown in Figure 2. The fact that voxels that would structurally be identified as gray matter were excluded (due to intravascular T2* effects) from the gray matter class in our analysis, is reflected in the lower gray/white-matter volume ratio (0.6) compared with that in the literature (Wenz et al, 1996). This approach also eliminates the need for co-registration of structural images and DSC images, which would add complexity to the analysis procedure.

In summary, a fast and fully automated quantitative perfusion analysis method is presented where all analytical steps are based on raw DSC data only. The resulting hemodynamic parameters extracted automatically from normal-appearing gray- and white matter, are in good agreement with the literature values. The gender- and age-related changes in CBF were best observed using cSVD with P_SVD = 0.1, and no age correlation was observed when applying the more complex iterative SVD regularization schemes. This may suggest that simple truncated SVD deconvolution approaches may be best suited for analysis of cerebral perfusion in the absence of significant alternations in cerebral hemodynamic properties. The automated perfusion analysis package can readily be implemented in a clinical workflow, thereby allowing intra- and inter-patient monitoring of normal perfusion values of all subjects where a DSC perfusion examination is acquired. Further, given the results obtained using normal brain tissue, the method presented may be valuable for quantitative assessment of alterations in hemodynamic parameters in tumors in response to therapy, and studies to investigate this are currently being initiated.