An Efficient Framework for Accurate Arterial Input Selection in DSC-MRI of Glioma Brain Tumors.

INTRODUCTION: Automatic and accurate arterial input function (AIF) selection has an essential role for quantification of cerebral perfusion hemodynamic parameters using dynamic susceptibility contrast magnetic resonance imaging (DSC-MRI). The purpose of this study is to develop an optimal automatic method for arterial input function determination in DSC-MRI of glioma brain tumors by using a new preprocessing method.
MATERIAL AND METHODS: For this study, DSC-MR images of 43 patients with glioma brain tumors were retrieved retrospectively. Our proposed AIF selection framework consisted an effcient pre-processing step, through which non-arterial curves such as tumorous, tissue, noisy and partial-volume affected curves were excluded, followed by AIF selection through agglomerative hierarchical (AH) clustering method. The performance of automatic AIF clustering was compared with manual AIF selection performed by an experienced radiologist, based on curve shape parameters, i.e. maximum peak (MP), full-width-at-half-maximum (FWHM), M (=MP/ (TTP × FWHM)) and root mean square error (RMSE).
RESULTS: Mean values of AIFs shape parameters were compared with those derived from manually selected AIFs by two-tailed paired t-test. The results showed statistically insignificant differences in MP, FWHM, and M parameters and lower RMSE, approving the resemblance of the selected AIF with the gold standard. The intraclass correlation coefficient and coefficients of variation percent showed a better agreement between manual AIF and our proposed AIF selection than previously proposed methods.
CONCLUSION: The results of current work suggest that by using efficient preprocessing steps, the accuracy of automatic AIF selection could be improved and this method appears promising for efficient and accurate clinical applications.

Introduction

Dynamic susceptibility contrast-enhanced magnetic resonance imaging (DSC-MRI) is a common technique for quantifying perfusion parameters of the brain, i.e. cerebral blood flow (CBF), cerebral blood volume (CBV), and mean transit time (MTT), which have widespread applications in brain tumor grading and treatment monitoring [1-4]. Quantification of DSC-MRI can be performed using relative or absolute quantification. The former is performed by directly assessing the curve shape descriptors, while the latter requires accurate estimation of the parameters based on the kinetic theory and therefore, allows for patient-specific evaluation, regardless of the protocol and scanner being used.

Nonetheless, one of major problems in absolute quantification of DSC-MR images is accurate determination of the arterial input function (AIF), which describes the concentration of contrast agent during its passage from the arteries feeding the brain tissue. In clinical applications, AIF is typically obtained manually by drawing a region of interest around the middle cerebral artery (MCA) or internal carotid artery (ICA) [5-9]. Nonetheless, manual AIF selection procedure is time-consuming, subjective and it depends on the experience of the radiologist. Thus, developing automatic and accurate AIF selection method is important to achieve the reliable diagnosis.

In this context, Murase et al [8] and Mouridsen et al [10] proposed fuzzy c-means (FCM) and k-means clustering approaches, respectively. Due to poor reproducibility of FCM and k-means clustering techniques, Yin et al proposed agglomerative hierarchy (AH) technique [11], which showed high reproducibility compared with FCM and k-means clustering. However, evaluation of the AH method was carried out on healthy subjects, which may lack generalization for clinical applications, like the assessment of brain tumors. Two of such problems can be outlined: (1) the similarity of the curve shapes within tumorous and arterial regions, which causes incorrect selection of AIF curves within tumorous regions instead of arteries, (2) selection of truncated, noisy or non-ideal curves mistakenly as AIF.

To address these problems, in the current study, we have implemented an efficient AIF selection framework, consisting of a new preprocessing step to exclude tumorous, truncated, noisy and non-ideal curves followed by AH clustering to increase accuracy and precision of automatic AIF selection.

Material and Methods

In this retrospective study, DSC MR images of 43 patients who have been histopathologically-confirmed with glioma brain tumors (20 low-grade glioma, and 23 Glioblastoma Multiform (GBM)) were used for evaluation of methods.30 of data were downloaded from the Cancer Imaging Archive [12] and 11 of them were acquired from Imam Khomeini hospital. The study approval was obtained from the institutional ethics committee and the subjects were included only if they provided their informed consent.

TCGA DSC-MR imaging (T2*-based perfusion MRI) had been carried out on 3 Tesla MR Scanner (GE Medical Systems) with the following parameter specifications: TE/TR= 40/1900 ms, flip angle=90, FOV= 230 × 230 mm2, matrix=128 × 128, slice thickness= 5 mm and for each slice, 60 or 95 images were recorded at intervals equal to the repetition time.

Perfusion images of Imam Khomeini hospital were imaged with 3 tesla Siemens MR scanner with: TE/TR= 45/2340 ms, flip angle=90, FOV= 230 × 230 mm, matrix=128 × 128, slice thickness= 5, for each slice 50 images were recorded. Image acquisition had been performed before and immediately after injection of 0.2 mmol/kg the Gd-DTPA with a rate of 5 ml/sec.

Image Analysis

An Overview of DSC-MRI Quantification

Hemodynamic parameters, i.e. CBV, CBF and MTT were quantified based on the indicator dilution theory, from the acquired T2* images. First, the slices containing major arteries, such as MCA and ICA, were chosen. For each voxel, signal intensity changes over time were obtained [2] from:

[1]

where S(t) is signal intensities at time t, S is the mean baseline signal intensity,T* is the Gradient echo transverse relaxation time and k is the unknown proportionality factor with an assumed value of k=1. A linear relationship was assumed between contrast agent changes and transverse relaxation rates (1/T* ) [10] and the curves of concentration agent changes were obtained based on the following formulation:

[2]

Where C (t) is the contrast agent concentration changes over time t; S (t) is the signal intensity changes over time t; T is echo time; and k is a constant that relates to the contrast agent relaxivity and other unknown factor parameters, and is typically set to 1. CBV that it can be obtained by using concentration-time curve and the selected AIF [1,2]:

[3]

Where ρ=1.04 (g/ml) is defined as the brain density and KH=0.73 refers to hematocrit differences in large and small vessels. Regarding to indicator dilution theory, CBF is measured by deconvolution of concentration-time curve and AIF [2,13]:

CBF×R(t)=C [4]

R (t) is the residue function that describes the fraction of the contrast agent present within the tissue at time t. MTT is measured using central volume theorem [14]:

[5]

Gamma variate fitting can be used based on the following equation to isolate the first pass and eliminate recirculation [10,15]:

C(t)=A(t-t exp(-(t-t [6]

α is the measure of inflow velocity, β is the washout velocity, t is the bolus arrival time and A represents the maximum peak of concentration-time curve. Initial values should be defined beforehand for more reliable convergence. An alternative is the simplified gamma function curve fitting proposed by Chan et al [16] which shows lower errors and can be executed faster than the previous methods. Therefore, here, the simplified gamma variate function was used to fit on concentration-time curves for isolating the first pass [17]:

C(t' )=C [7]

t'=(t-t

where t is the bolus arrival time; C the maximum peak of concentration curve; α is the inflow velocity steepness measure; t is the time to reach the maximum peak of concentration-time curve.

Automatic Determination of AIF

The main discriminating characteristics of AIF concentration-time curve from tissue concentration-time curves are the maximum peak (MP), full-width at half-maximum (FWHM), area under the curve (AUC), and time to peak (TTP). Arterial curve has higher MP, larger AUC, lower TTP, and narrower FWHM in comparison with tissue concentration-time curves. In this work, the following steps were carried out for automatic AIF determination (Figure 1):

Figure1

Flowchart showing the automatic AIF determination processes carried out to DSC-MRI perfusion clinical data in the proposed study.

1- Motion Correction: Involuntary motions like breathing and heartbeat, which cause misalignment of consecutive images during DSC-MRI, were corrected. All of the images were aligned to the first pre-contrast primarily-acquired image based on a rigid transformation using SPM (available at ; version SPM99) and INRIAlign 1.01 (available at INRIAlign.html) [18,19]. As smoothing causes errors in DSC-MRI quantification, no smoothing was applied.

2- Exclusion of Tumorous Regions: As the curves belong to the tumorous regions mimic arterial curves, the possibility of choosing false curve increase. Hence, the curves of the tumorous regions should be excluded before AIF selection. Mean signal intensity in the baseline for the tumorous region part is higher than other parts of the brain tissue. By using k-means clustering, signal intensity curves were classified to five clusters. Thus, the cluster with the highest mean baseline value was removed as the tumorous cluster (Figure 2) [1,20].

Figure2

Tumorous region removal. Two sample of signal intensity curves: (a) tumorous region and (b) white matter. Signal intensity in tumorous region has higher signal intensity in baseline part in comparison with white matter. The subplots (c) and (d) shows two Sample of T* images before and after tumorous regions removal.

3- Noise Reduction: Signal intensity changes during the first pass of contrast agent could hit the noise floor level, which results in truncated concentration-time curves [21] and can cause significant systematic errors in CBV and CBF parameter quantifications [22]. As mentioned above, concentration-time curves consist of three parts, including baseline, first-pass and steady-state. Typical concentration-time curve ascends from baseline part to the maximum peak, then it descends from the maximum peak to the steady state. In this simple precondition, if in these ascending and descending paths, contrast agent concentration reaches zero, this point is considered to be contaminated by noise. Therefore, average of the concentration values of a point before and one after the noisy point is replaced with the value of the noisy point to reform the curve (Figure 3).

Figure3

Truncated concentration curve correction. A Sample of truncated concentration curves were shown in subplot (A), truncated time point highlighted by red circle. Mean value of after and before concentrations were replaced in truncated time point (red circle). Correction of this curve displayed in subplot (B).

4- Elimination of Distorted Curves: Due to various scanning artifacts, such as partial volume effects, physiologic pulsations, and shifts in voxels, some highly distorted concentration curves may be resulted [10,11,15,23,24]. These curves induce errors in the automatic AIF selection procedures. For exclusion of these irregular curves, we propose a new criterion. Typical concentration-time curves are bell-shaped and form an ascending-descending path between the baseline and steady-state portions. If a horizontal line is intersected with the bell-shaped portion of the concentration-time curve, two intersecting points will result. However, noise-distorted curves show fluctuations in an ascending-descending path. Compared with a typical concentration-time curve, if a horizontal line is intersected with a distorted concentration-time curve, more than two intersection points result. So, for discarding these distorted curves, concentration-time curves were intersected separately by 5 horizontal lines with values in the range of 0.5 to 0.9 of the maximum peak (0.1 unit intervals). If the horizontal lines and concentration-time curves return more than 2 intersecting points, this curve will be regarded as a distorted and irregular curve and is discarded from further processing.

5- Selection of Candidate AIF Curves: Through this step, by using arterial curve characteristics, such as higher MP, early TTP, and lower FWHM, non-arterial curves were discarded and candidate AIF curves were maintained. This was performed by taking the average MP of the remaining curves from the previous steps. The MP of each curve was compared with the average MP and if it was lower than the average, the curve was excluded. This procedure was repeated based on TTP and FWHM and if TTP and FWHM of the curves were larger than their average, they were discarded. Ultimately, several candidate AIF curves would be maintained.

6- Imposing Partial-Volume Criterion: Bleeker et al [9] proposed a criterion to detect of the curves contaminated with the partial volume effect (PVE). A simplified form of this criterion was employed by Yin et al, [11,23,24], which will be described briefly here. First, the area under the gamma-variate function fitted on the first-pass of a concentration-time curve was calculated (AUC). Second, initiation of the steady-state was assumed as the first time-point which was lower than three-tenths of a MP, the ten succeeding time points after the initial point of the steady-state portion were integrated (AUC). Then, the mean ratio of AUC to AUC was calculated, and if the result fell outside the accepted range (mean ratio±20 %), the curve was considered to be affected by PVE and was discarded.

7- Clustering: Finally, AH clustering was applied to classify the remaining curves into five clusters, including white matter, gray matter, arterial, venous and other curves, based on their similarity. AIF cluster and curve were selected by the measure of M=MP/ (TTP × FWHM).

The flow chart of our proposed automatic AIF selection procedure (Proposed-AIF) is illustrated in Figure 1.

Comparison with Automatic AIF Selection Method by Yin et al

One of the most recent and efficient methods for automatic AIF selection has been proposed by Yin et al. on normal subjects [11,15,23,24]. In our study, we evaluated their method by clinical brain tumor data.

As it mentioned similarity of tumorous and AIF curves caused Yin et al method which was performed using normal data, get trapped by tumorous curves and AIF wrongly was obtained in the tumorous region. To solve this problem, at first, the tumorous region was excluded as it explained in our proposed method. Afterward, the following preprocessing steps were used in Yin et al proposed method exactly applied:

1- In the first step, the area under concentration curves was calculated and 90 percent of them, with the lowest value were discarded as tissue curves.

2- Then, Roughness measurement used to discard noisy curves:

[8]

25% of the curves with largest integral were excluded as noise curves.

3- PVE-affected curves were excluded using the condition proposed by Bleeker et al [9] which is described in step 6 of the previous section. However, they used gamma fitting of Equation 6 to separate first pass which had longer time expectation time and errors in fitting.

AH clustering methods were applied with cluster number of five, and automatic AIF acquired by Yin et al method (to be called Yin-AIF hereafter) was applied and AIF was selected by AH clustering.

There were some difference between our proposed and Yin et al method. At first, using tumorous data, which its tumorous curves were challenging issue. Furthermore, in our proposed method noisy, distorted and non-arterial curves excluded with new preprocessing steps which completely different with Yin et al preprocessing steps. Moreover, In PVE condition we used simplified gamma fitting which was more accurate and faster than Yin et al used method.

In both methods, final AIF obtained by using the same AH clustering method to determine the preprocessing steps effect on AIF selection.

Statistical Analysis

For evaluation of the performance of each clustering method, an experienced radiologist was asked to insert regions of interest in (ROIs) around middle carotid artery (MCA) and internal carotid artery (ICA). Among the selected AIFs, the representative AIF for that patient (True-AIF) was selected based on curve shape characteristics like higher MP and AUC, and lower FWHM and TTP. For automatic and True AIFs, shape parameters i.e. FWHM, TTP, MP, AUC and M (M=MP/ (TTP × FWHM)) were calculated.

The difference between automatic AIF and True-AIF was computed based on the root mean square error:

[9]

Where n is the number of dynamic scanning volumes.

For comparing mean values of AIF shape parameters obtained by manual and automatic methods, Student’s paired t-test analysis was performed and p-values less than a level of 0.05 were considered to be statistically significant. In ideal situations, the mean values of AIFs were derived from the manual method, as the gold standard of comparison [8,11,15,23,25], and the automatic method is exactly the same. Therefore, to accept an automatically selected AIF, the mean values of True-AIF and automatic AIF should show insignificant differences (p-value>0.05).

Furthermore, the intra-class correlation (ICC) coefficient was performed to quantify the degree of agreement between automatic- and True- AIF. In ICC analysis, if automatic-AIF is similar to True-AIF, the ICC should be close to 1. Another measure for concordance of our automatically selected AIF and True-AIF in comparison with Yin-AIF is the coefficient of variation (CV), which can be calculated by the ratio of standard deviation to the mean (STD/Mean) of each shape characteristic parameter. For each parameter, decreased CV is desirable.

Results

Gamma-Variate Fitting

For evaluating the accuracy of gamma-variate fitting methods, complicated (Equation 6) and simplified gamma-variate fitting (Equation 7) were respectively used in our proposed method and the one proposed by Yin et al. The mean value of MP, TTP, and FWHM parameters were extracted from 4311 concentration curves. These curves were generated in the sixth step (Bleeker condition) of automatic AIF procedure for 43 data. By using the TTP point, a 14-second time window was defined to include the first-pass part of the concentration curve and decrease steady state and baseline effects on calculation RMSE between gamma-variate fittings and original concentration curves. Also, execution time for each fitting was computed.

According to Table 1, MP and TTP in the simplified form of gamma fitting were exactly equal to mean values of original concentration curves and the FWHM mean value was close to the mean FWHM of original concentration-time curves. RMSE in simplified gamma fitting decreased as expected. Figure 4 shows a good agreement between the original concentration-time curve and simplified gamma fitting. Moreover, this method had required execution time compared with the initial gamma fitting (Table 1). Statistical analysis showed significant difference for FWHM, RMSE and execution time between complicated gamma-variate fitting and simplified gamma fitting. Hence, using simplified gamma fitting was more efficient in terms of computation burden and accuracy than the complicated method.

Table 1

Comparison between mean values of gamma variate fittings shapes parameter obtained using two kind of gamma fitting and original curve.

	MP	TTP	FWHM	RMSE	Execution Time
Complicated gamma fitting	0.0433	21.96	5.599	0.0052	0.209
Simplified gamma fitting	0.0381	22.56	4.695	0.0044	0.031
Original Curve	0.0381	22.56	4.381	-	-

Figure4

Two gamma fitting method accuracy. A sample of gamma variate function fitted (blue curves) on concentration curve (black curves). In subplot (A) complicated gamma variate fitting was carried out in last studies and in subplot (B) Chan et al simplified gamma variate were used in our proposed method were illustrated. Its intuitively obvious gamma fitting used in proposed method has better fitting.

Automatic AIF Selection

According to the procedures explained, AIFs were obtained separately based on our new proposed method and Yin’s approach. As described before, the True-AIF shows discriminating characteristics like earlier TTP, narrower FWHM, higher MP and M. The shape parameters for the two automatic method and True-AIF were compared.

The results of calculating CV for True-AIF, our proposed and Yin’s methods are indicated in Table 2. CV was lower for the four parameters, i.e. MP, TTP, FWHM, and M, based on our proposed method in comparison with Yin’s method. However, the CV for AUC parameter based on Yin’s method was less than our proposed method. Overall, the decrease in CV for four of the descriptive parameters indicates that our proposed method resulted in less standard deviation in estimation of AIF parameters.

Table 2

Coefficient of variation percent illustrated for our proposed, Yin-based and True AIFs.

	MP	AUC	TTP	FWHM	M
True-AIF	32	38	23	37	50
Our proposed method	31	86	30	41	68
Yin-based AIF	38	37	31	60	170

Table 3 shows the mean values for the calculated parameters based on our proposed and Yin’s methods. It is apparent that the differences between the mean values of parameters generated by the two automatic methods are significantly different. ICC analysis value for MP increased from 0.41 to 0.7, which indicates good agreement between MP of our proposed method and True-AIF. TTP values obtained from both automatic methods are in good agreement with each other and there is no significant difference (P=0.66). However, in Table 4, ICC showed better performance of Yin’s method with respect to our proposed method, for TTP, as ICC for our proposed method and Yin’s method were 0.66 and 0.86, respectively. The mean FWHM for AIF selected by new method was closer to True-AIF and the difference between automatic methods is significant (P=0.001). ICC values for FWHM increased from 0.35 to 0.44. AUC for both automatic methods are in agreement with True-AIF and there was no significant difference (P=0.759). Moreover, ICC value for both methods were approximately the same. The M value for our proposed method is in line with True-AIF compared with Yin’s method; in addition, this difference was statistically significant (P=0.037). The increase in ICC values for M was remarkable and the values changed for Yin’s and our proposed methods from 0.08 and 0.58, respectively. RMSE for automatic and True-AIF differences, notably decreased in our proposed method and statistical analysis of pair t-test showed the significant difference (P=0.002). Overall, the mean values of the AIF shape parameters for True-AIF and our proposed method are more in agreement compared to Yin’s method and True-AIF.

Table 3

Comparison of the AIFs shape parameters obtained from different automatic methods and the true AIF.

	MP	AUC	TTP	FWHM	M	RMSE
True-AIF	0.111	0.434	24.23	3.492	0.00155	-
Our proposed method	0.109	0.418	24.95	3.64	0.0016	0.0176
Yin-based AIF	0.095	0.403	25.70	2.59	0.0036	0.0261
P-value	0.002	0.759	0.367	0.001	0.037	0.002

Table 4

Intraclass correlation coefficient agreement measure between automatic method and True AIF.

	MP	AUC	TTP	FWHM	M
Our proposed method	0.70	0.48	0.86	0.44	0.58
Yin-based AIF	0.41	0.41	0.66	0.35	0.08

In Figure 5, comparison of the AIF derived from manual and automatic methods is indicated for a patient. As it can be observed from this figure, Yin’s method mistakenly selected noisy curves instead of AIF which is different from the AIF selected by the manual method. But our proposed method was capable to identify the same AIF as the one selected by the manual method.

Figure5

AIF selection in automatic methods. The plots of (A) and (B) are two sample of AIF were selected by Yin method and the plots of (C) and (D) are the AIFs which selected by new proposed method for the same data. AIF obtained by proposed method are more ideal and it was the same of True AIF.

Discussion

Determination of AIF has an essential role in absolute quantification of DSC-MRI perfusion maps, where errors in accurate AIF selection can lead to significant errors in estimation of the perfusion maps [23]. Manual selection of AIF depends on the expertise of the radiologist to identify large arteries, such as MCA and ICA and to choose the correct AIF curve by assessing arterial curve characteristics [7,9,10,21]. Therefore, manual selection of AIF is tedious, subjective, and lacks reproducibility. To overcome these limitations, several automatic methods have been proposed for AIF selection.

In the automatic method proposed by Mouridsen et al [10], two new preprocessing steps were used to exclude non-arterial curves including tissue, noisy and truncated curves. Arterial curve had higher CBV in comparison with tissue curves, including gray matter and white matter, so they calculated area under the concentration curves (AUC), which is proportional with CBV, and 90% of smallest AUC values were considered as tissue curves and discarded. Furthermore, 25% of the curves with largest AUC were eliminated as noisy curves. Then, they used K-means clustering to classify the remaining curves to reach the final AIF. However, K-means clustering has low reproducibility. Yin et al [11] used the two preprocessing steps suggested by Mouridsen for discarding tissue and noisy curves. They also calculated the Bleeker condition for excluding PVE-affected curves. Instead of K-means, they examined AH clustering which has higher reproducibility, to classify the remaining concentration curves accurately. The result showed that AH clustering could determine AIF much better than K-means. One of the limitations of the work carried out by Yin et al was using healthy volunteer data to evaluate their method.

Our evaluation of Yin’s method on DSC-MR images of brain tumor patients revealed several issues. One of the major problems was similarity of the curves of arterial and the tumorous region that caused incorrect selection of tumorous curves instead of arterial curves. Yet after removing tumorous regions, AIF selection by Yin’s method did not improve and the AIF obtained by Yin’s method showed significant differences with the True-AIF, which was selected manually by an expert radiologist. Evaluating the discarded and AIF curves in Yin’s method indicated that many arterial curves were incorrectly excluded as noisy curves and some truncated and noisy curves were chosen as AIF. As a result, using preprocessing steps to exclude tumorous, noisy and truncated curves were necessary to determine AIF with high accuracy.

As a solution to this problem, in the current study, truncated curves were corrected and a new condition was used to exclude noisy curves. Furthermore, AIF curve has a higher maximum peak, area under curves and lower FWHM and time to peak. Yin and Mouridsen just used area under curve feature of AIF and excluded 90% of the curves. Instead, in this work, we assessed more descriptive features of AIF, including maximum peak, time to peak and FWHM. As a result, in our new condition, by using more AIF curve features, the higher number of non-arterial curves were excluded.

The condition proposed by Bleeker et al for excluding PVE-affected curves depends on correct isolation of the first-pass of the perfusion curve by gamma-variate fitting. Investigation of the gamma-variate fitting used in Yin’s method showed noticeable fitting errors more computational burden. Instead, in this study, we used a simplified form of the gamma-variate fitting. This simplified gamma-variate fitting was faster and statistical analysis showed higher accuracy. As a result, using this fitting can increase the accuracy of excluding PVE-affected curves. Finally, AH clustering was applied and AIF was automatically determined.

Our new proposed method with improved preprocessing steps showed more accurate results for AIF selection in comparison with Yin’s method and the errors remarkably decreased. This improvement in the results can be due to adopt new preprocessing steps for correcting truncated curves, eliminating noisy curves with new conditions, using more arterial concentration curve characteristics, such as maximum peak, FWHM, and time to peak for discarding more non-arterial curves, and employing a new simplified gamma-variate function for correction of recirculation peak.

However, our work had some limitations. First, images of some patients had been acquired with fewer slices, which complicated identification of the slice containing large arteries, i.e. MCA and ICA. Second, our study population was limited to 43 brain tumor patients. Generalizing this method to clinical applications, requires assessments of the technique on larger patient populations including more diverse diseases, like acute stroke, arterial stenosis and other brain disorders.

Conclusion

In the present study, we developed a new automatic AIF selection framework for quantification of DSC-MR images of brain tumor patients by using new preprocessing steps and AH clustering method. It was shown that the AIF obtained by our method is in agreement with manually-selected AIF, and has the potential to produce robust and accurate quantification of perfusion.