Automatic Detection and Quantitative DCE-MRI Scoring of Prostate Cancer Aggressiveness

Purpose To develop a robust and clinically applicable automated method for analyzing Dynamic Contrast Enhanced (DCE-) MRI of the prostate as a guide for targeted biopsies and treatments. Materials and methods An unsupervised pattern recognition (PR) method was used to analyze prostate DCE-MRI from 71 sequential radiotherapy patients. Identified regions of interest (ROIs) with increased perfusion were assigned either to the peripheral (PZ) or transition zone (TZ). Six quantitative features, associated with the washin and washout part of the weighted average DCE curve from the ROI, were calculated. The associations between the assigned DCE-scores and Gleason Score (GS) were investigated. A heatmap of tumor aggressiveness covering the entire prostate was generated and validated with histopathology from MRI-ultrasound fused (MRI-US) targeted biopsies. Results The volumes of the PR-identified ROI’s were significantly correlated with the highest GS from the biopsy session for each patient. Following normalization (and only after normalization) with gluteus maximus muscle’s DCE signal, the quantitative features in PZ were significantly correlated with GS. These correlations straightened in subset of patients with available MRI-US biopsies when GS from the individual biopsies were used. Area under the receiver operating characteristics curve for discrimination between indolent vs aggressive cancer for the significant quantitative features reached 0.88–0.95. When DCE-scores were calculated in normal appearing tissues, the features were highly discriminative for cancer vs no cancer both in PZ and TZ. The generated heatmap of tumor aggressiveness coincided with the location and GS of the MRI-US biopsies. Conclusion A quantitative approach for DCE-MRI analysis was developed. The resultant map of aggressiveness correlated well with tumor location and GS and is applicable for integration in radiotherapy/radiology imaging software for clinical translation.

inTrODUcTiOn Prostate cancer is the second most common cancer in American men (1). Clinical decisions, related to the need for prostate biopsy, define target areas for biopsy and regions that require attention in focally directed therapy, are multifactorial and complex. Delivery of targeted radiation dose to high-risk tumor areas in lieu of escalating dose to the entire prostate will reduce overall complication risks (2), a strategy that requires robust and quantitative imaging of potential tumor lesions (3,4). The use of multiparametric MRI (mpMRI) for prostate cancer is rapidly evolving because of its growing availability and ability to combine functional [perfusion via dynamic contrast enhanced (DCE-MRI) and diffusion via diffusion-weighted imaging (DWI)], and anatomical information [T2-weighted (T2w) MRI]. DCE-MRI is an established component of prostate mpMRI, and DCE-MRI alone has reported sensitivity and specificity ranges of 46-96 and 74-96%, respectively (5)(6)(7)(8)(9)(10). The role of DCE-MRI for tumor diagnosis and assessment in the current version of Prostate Imaging Reporting and Data System (PI-RADS) (version 2), however, is minimized (11). DCE-MRI is only considered when DWI in the peripheral zone (PZ) is indeterminate. The PI-RADS Steering Committee justified the reduced role of DCE-MRI in prostate cancer assessment by the lack of expert consensus, reflective of the difficulty of interpreting DCE-MRI sequences by eye. An automated and quantitative assessment has the potential to improve consistency in identifying high risk prostate volumes.
There has been considerable efforts to standardize the analysis of DCE-MRI (12,13) with pharmacokinetic modeling (14); however, there is significant variability in the calculated rate constants. For instance, K trans , the volume transfer coefficient that measures capillary permeability (14), calculated on the same DCE-MRI data by participants from several academic medical centers yielded a within-subject coefficient of variation of 0.59 (15). This modeling approach is impeded by difficulties in estimating the arterial input function (16), insufficient temporal resolution (17,18) and low signal-to-noise ratio (SNR). Previously, it has been shown that the application of pattern recognition (PR) techniques to DCE-MRI overcomes the challenges of low temporal resolution and low SNR (19)(20)(21). Here, the PR approach is extended by introducing novel quantitative features that allow for comparisons between patients. The method automatically: (i) delineates the region in the prostate with increased perfusion; (ii) allocates the dominant lesion to either the PZ or transition zone (TZ); (iii) assigns a DCE scorea zone-specific quantitative measure of aggressiveness; and (iv) generates a spatial map of tumor aggressiveness based on DCE score. The result is a color-coded map of aggressiveness on a pixel level, which may serve as a guide for targeted biopsies and treatments.

Patients and Mri acquisition
An Institutional Review Board (IRB) approved a protocol for retrospective review of mpMRI exams from prostate cancer patients. The IRB waived the need for informed consent. A total of 71 sequential patients, presenting for evaluation for radiation treatment (RT) between 2012 and 2015 and who underwent mpMRI exams on a 3T Discovery MR750 (GE, Waukesha, WI, USA) were identified from the departmental database. mpMRI exam consisted in part of: (i) axial T2w-MRI of the pelvis: resolution 1.25 mm × 1.25 mm × 2.5 mm; field of view: 320 mm × 320 mm; slice thickness = 2.5 mm (no gap); 72 slices, and (ii) DCE-MRI-12 series of T1 weighted (T1w) at 30-34 s temporal resolution, TR = 3.77-4.05 ms, TE = 1.69-1.78 ms acquired following intravenous bolus injection of a paramagnetic gadolinium chelate-0.1 mmol of gadobenate-glumine (Bracco Diagnostics Inc., Princeton, NJ, USA) per kilogram of body weight. The contrast is administered with a power injector (Spectris, Medrad Inc., Warrendale, PA, USA) at 2 mL/s and followed by a 20-mL saline flush.

Quantitative Dce-Mri analysis
The analysis pipeline is presented schematically in Figure 1. The DCE-MRI series were uploaded in MIM (MIM, Cleveland, OH, USA) and all subsequent analyses are performed using MATLAB R2014a (MathWorks, Natick, MA, USA) plugins. Prostate, PZ, and sample gluteus maximus (GM) volumes were manually contoured ( Figure 1A). The individual steps of the analysis are described below: • Motion correction ( Figure 1B). Each imaging set in the DCE-MRI series was aligned to the preceding one by finding the affine transformation (shifting, scaling, and rotation) that maximized their mutual information (22). Pixels within the prostate contour and extension by 18.75 mm (15 pixels) were considered. For each pair of before and after-correction images, the sum of squared pixel differences was computed and motion correction was carried out only when this sum was smaller post-correction. • Identification of well-perfused (suspicious) regions of interest (ROI) in the prostate. Non negative matrix factorization (NMF) was applied to signal-vs-time curves of all pixels within the prostate ( Figure 1C) (19). Briefly, if D is the data matrix, containing the individual pixel's signal-vs-time curves in its rows (baseline corrected by the average of the pre-contrast points), then D can be represented as a product of k basic temporal contrast signatures S(t) and their weights W(X) in each pixel, i.e., D ~ W × S under the constraint that all elements of W and S are non-negative. k = 3 is estimated by Principal Component Analysis of D as the number of significant Principal Components (23). The signal in the i th pixel in data matrix D can be represented as: The well-perfused pattern, Sj, is automatically selected from the k NMF patterns as: where AUCt1-t2 (Si) is the area under the curve (AUC) of the pattern Sm between times t1 and t2. Let Wj be the well-perfused (wp) contribution map (Wwp) associated with the well-perfused Early AUFC Early AUFC (SROI) pattern Sj (Swp). Wwp is an intensity map of the "well-perfused" pixels in the data ( Figure 1D).

• Segmentation of well-perfused ROI. Five methods to segment
Wwp were tested for identification of the ROIwp, the well-perfused ROI. The first was based on Otsu thresholding (24); the others were based on a parameter β = 40, 50, 60, and 70, describing how "pure" is the pattern in a pixel (20). For instance, if β = 60, the ROIwp is defined by pixels for which the weight of the well-perfused pattern Wwp is >60% of the total sum of the weights (Eq. 1). Simultaneously, the assignment of ROIwp to PZ vs TZ is made by using variable fraction (10, 15, and 20%) of ROIwp in PZ. ROIwp in PZ and TZ were denoted pzROIwp and tzROIwp, respectively ( Figure 1E). The volumes in PZ outside of pzROIwp were considered normal appearing tissues PZ (NATPZ) and TZ (NATTZ). • Quantitative feature extraction. The signal-vs-time curve SROI is reconstructed using Eq. 1 for pixels within pzROIwp or tzROIwp using k = 3; SROI is an approximation of the average of the curves within the region. By using the NMF reconstruction this signal is effectively de-noised (19). SROI and the average curve from the muscle volume SGM were used for quantitative analysis (Figure 1F). Six quantitative features (25,26), summarized in Table 1 and illustrated in Figure S1 in Supplementary Material, were computed based on these two signals ( Figure 1G). The features are divided in three groups: (i) Features extracted directly from the DCE-MRI intensity curve; (ii) Features based on the fitting of SROI and SGM to a bi-exponential model; and (iii) Features calculated by fitting these curves to a late enhancement linear model. The parameters from the bi-exponential model are: initial static intensity s0, plateau intensity sm, start of enhancement t0, time-to-peak τ, and washout slope wobiexp ( Figure S1 in Supplementary Material). The wash-in slope wi was computed as (sm − s0)/τ. The late enhancement linear model computed the slope and intercept of a linear fit of the intensities between DCE270 and DCE330. This slope is referred to as wolinear. • Map of aggressiveness. Let F be the value of a quantitative feature for the curve SROI in the identified volume ROIwp. Using the distribution of F in associations with GS, a spatial map of tumor aggressiveness was computed ( Figure 1H) For convenience, the feature map is represented using a 1 to 10 scale, grouping the Wwp i F/ω in 10 bins. Feature range was bound by 5th and 95th percentile of the feature values in the entire population of subjects with positive biopsies. DCE-feature scores are allocated by uniformly dividing the bounded interval in ten regions. Pixels with feature values smaller than the lower bound or larger than the upper bound were assigned a DCE-featurescore of 1 and 10, respectively.

statistical analysis
Spearman correlation coefficient (ρ) of the DCE-scores with the highest GS from the biopsy session was computed (patient level analysis). The GS were grouped in a clinically relevant categories, separating GS7 (3 + 4) and GS7 (4 + 3) (27). To study the ability of the features to discriminate between indolent (GS = 6) and aggressive (GS ≥ 7) lesions, the association between the quantitative features and these two classes was investigated. The significance of the median difference between groups was evaluated using Kruskal-Wallis test, followed by Dunn's posttest. In addition, the AUC of the receiver operating characteristics (ROC) analysis was computed. A subset of the patients underwent MRI-ultrasound fused targeted (MRI-US) biopsies. The statistical test above were carried out for the DCE-scores of each biopsy location and biopsy GS (biopsy level analysis). Statistical analysis was performed using MATLAB. All tests were two-sided. Significance was set at a p-value <0.05.

Patients
The flow diagram of inclusion and exclusion criteria of the analyzed patients is presented in Figure S2 in Supplementary Material. Nine patients (9/71, 12%) were excluded because: four patients lacked complete DCE-MRI series and five had suboptimal contrast data due to acquisition artifacts or incomplete contrast administration. The clinical characteristics of the analyzed cohort are summarized in Table S1 in Supplementary Material. The highest GS from the biopsy session was used.

Motion correction
For each pixel in the dataset, DCEN was aligned to DCEN − 30, (N = 90, 120, …, 330) and motion was quantified as the mean absolute pixel gray scale-intensity difference between DCEN and DCEN − 30, before and after correction. For 27 patients (44%), the movement correction procedure resulted in essentially no change (change of pixels squared differences before and after: median = 0.02%; range = −2.38 to 1.51%). For 35 patients (56%), the procedure yielded significant improvement: median = 9.93%, range = 2.61-46.84%. The procedure is illustrated for the patient with the largest correction ( Figure S3 in Supplementary Material).

identification and Delineation of Well-Perfused (suspicious) rOi in the Prostate
All five segmentation methods resulted in significant correlations between the volumes of ROIwp and GS. The Spearman correlation ρ of the volume of ROIwp and GS (6, 3 + 4, 4 + 3, ≥8), its p-value, and the number of patients assigned to PZ or TZ, or rejected due to volumes less than 0.5 cc are presented in Table 2. While in general, with increasing β, the correlations with GS increased, the resulting smaller ROIwp volumes caused a large number of patients to be rejected due to ROIwp < 0.5 cc. β = 60 and 10% PZ vs TZ were selected. The median size of the well-perfused region ROIwp across the patients was 1.95 cc (range 0.23-12.21 cc). Following 0.05 cc cleanup for disconnected small volumes, the resulting volumes ROIwp were on average 18.19% smaller and ranged in size between 0.16 and 11.52 cc. Eight patients with ROIwp <0.5 cc were eliminated from further processing. For the remaining 54 patients, the mean volume of ROIwp was 2.57 cc (median = 1.84 cc, range 0.51-11.52 cc). Forty-one (76%) of the tumors were in the PZ (pzROIwp). The mean volume of pzROIwp was 1.08 cc (median = 0.80 cc, range 0.15-5.01 cc) and the mean volume of tzROIwp was 3.09 cc (median = 1.28 cc, range 0.80-11.52 cc). The distributions of pzROIwp volumes in association with GS divided in three (6, 3 + 4, >4 + 3) and four groups   Figure 2A. In both cases, Spearman correlation coefficients were significant.

associations of Quantitative Features with Prostate cancer aggressiveness (Patient level analysis)
Six quantitative features ( Table 1) were computed based on the signal SROI for the patients with suspicious area in PZ (pzROIwp). The correlations between these features and the highest GS (6, 3 + 4 and >3 + 4) are presented before and after normalization with the signal from the muscle reference SGM (Figures 2B,C).
Note that the correlations with GS improve substantially following normalization with the muscle signal, reaching significance (p < 0.05) in all cases but one. When GS is divided in four groups: 6, 3 + 4, 4 + 3, 8 − 10 the same five features remained significant ( Figure 2D). It should be noted that none of the features were significant for tzROIwp (data not shown). In addition, comparisons including the distribution of the quantitative features in NATPZ and NATTZ we carried out ( Figure S4 in Supplementary Material). Their corresponding AUC for discrimination between NATPZ and pzROIwp ranged between 0.62 and 0.99. AUCs of DCE-scores in tzROIwp were also high. In view of the small number of lesions in TZ, however, these results should be treated with caution.

correlation of Quantitative Features with Mri-Us Biopsies (Biopsy level analysis)
Sixteen patients underwent MRI-US biopsies and had total of 35 positive biopsies ( Figure S5 in Supplementary Material). After  confirmation of the location of the actual biopsy, the region with needle tracks was back-projected onto T2w. The overlap between the biopsy region and the well-perfused region ROIwp was denoted as usROIwp. Quantitative features were computed for each biopsy in PZ (29/35) using the average signal in usROIwp. The associations of the sixfeatures for PZ tumors and the targeted biopsy GS (divided in three and four groups) are presented in Figures 3A,B. The same four features were significant using three, four GS groups and when tumors were split into indolent (GS = 6) and aggressive (GS ≥ 7) [wash-in ratio reached marginal (p = 0.051) significance] (Figure 3C). The features resulted in AUC for classifying indolent vs aggressive PZ lesions ranging between 0.85 and 0.95 ( Figure 3D).

Map of aggressiveness
Using Wwp and the DCE-score F from ROIwp, each pixel in the prostate was scored for aggressiveness using 1 to 10 scale. In Figures 4A,B, two axial slices of the prostate are shown with indicated regions representing the intersection of the identified well-perfused region ROIwp and the biopsy locations. A schematic representation of the prostate with the biopsy needle tracks (yellow) in the target volumes (1 in the green, 2 in the red, and 1 in the blue target) is depicted in Figure 4C. On pathology review, the biopsies were assigned GS = 6 (green target) and GS = 7 (for both biopsies in the red target). Map of aggressiveness was generated for early AUFC ratio feature. In Figures 4D,E, the map is overlaid on the patient axial MRIs slices. The areas of the positive biopsies are clearly identified. The color scheme of the maps is given on the right. There are no pixels with DCE-scores 9 and 10, indicating lack of areas with high GS. The highest DCE-score on the maps (left to right) was 6 and 8, consistent with the biopsy findings. The MRI-US dataset was utilized to confirm the accuracy of the global map of aggressiveness. On one hand, the map of aggressiveness featureWwp was generated using only Wwp and F. On the other  hand, the feature value for each biopsy ROI was estimated in the previous section using usROIwp. In Figure S6 in Supplementary Material, the x-axis corresponds to the feature computed using usROIwp. The y-axis is the mean value of the map of aggressiveness over the region of the biopsy, mean [featureWwp(mrusROI)].
There was an excellent concordance between these two estimates for each biopsy ROI (n = 29) for each feature: the correlation coefficients ranged between 0.85 and 0.99; the slope for the linear fit ranged between 0.87 and 0.98. This graph also confirms the co-localization of the "hot" spots of the map of aggressiveness with the biopsy location.

DiscUssiOn
DCE-MRI plays an ancillary role in tumor assessment to T2-weighted (T2w) and DWI in the contemporary recommendations for prostate cancer review and scoring in PI-RADSv2. The diminished role of DCE-MRI is confounded by: (i) difficulties for radiologists to quantitatively assess multiple DCE-MRI series via image viewing; and (ii) lack of reproducibility of the quantitative parameters when modeling approaches are used. In addition, other factors, such as chronic inflammation and differential characteristics of the prostate zones further complicate the DCE-MRI interpretation. The thesis of this work is that DCE-MRI can reliably identify and consistently score prostate cancer lesions if the region of well-perfused tissue is properly delineated and sensible quantitative features, using internal signal normalization are used. The presented approach considers the signal-vs-time curves from each pixel in the dataset and an unsupervised PR technique is used to identify the pixels with the characteristic for tumor perfusion pattern. Quantitative features are calculated from a representative signal-vs-time curve from these pixels. Several of these features, especially after normalization with signal from the muscle, were significantly correlated with biopsy GS, and differentiated indolent (GS = 6) from aggressive tumors (GS ≥ 7). For instance, Litjens et al. (28) found that the τ-induced volume and wash-out had a large relative importance in their ROC analysis for classification of cancer (normal/benign). Vos et al. (10) found that wash-in and wash-out were statistically different between low-grade (GS = 6) and high-grade [GS > 7(4 + 3)] in the PZ. The Area under the Receiver Operating Characteristics Curve (AUROC) for indolent vs. aggressive in this work were 0.72 for wash-in and 0.65 for wash-out. In comparison, after normalization by muscle, the AUROC for wash-in and wash-out of 0.85 and 0.88 are reported here. Chen et al. found a significant correlation between wash-out and GS, with an AUROC for indolent/aggressive classification of 0.88 (26). Significant correlation and coincidentally the same AUROC of 0.88 is reported here. It should be noted that unlike these previous reports, where few significant features are found, the majority of the evaluated features here were significantly correlated with GS, both on patient and biopsy level analysis. This is attributed to the robustness of well-perfused ROI selection and fitting a single DCE curve of high SNR, rather than pixel-by-pixel in the ROI (10). This reduces the errors associated with noise, data acquisition and motion artifacts.
Here, a novel pipeline for automatic generation of a wellperfused ROI is presented. The method performs well in datasets acquired at low temporal resolution. Most current research and PI-RADS guidelines suggest that high-temporal resolution be used. While this is a requirement for compartment model analysis, the straight-forward features extracted here are significantly correlated with GS even at low temporal resolution (Δt = 30 s).
This study has several limitations. First, this is a retrospective study with all inherent limitations of the design. Second, the application of the thresholds require contours of PZ and TZ. The use of prostate atlas will decrease significantly the need of manual contouring (29). And finally, the generalizability of the determined thresholds to other MRI sequences, vendors, magnetic field strengths, and coils (endorectal vs body) should be investigated. Until then, it should be assumed that the thresholds are valid for DCE data, acquired under identical conditions.
In conclusion, a quantitative approach for DCE-MRI analysis is developed and the resultant map of aggressiveness is integrated in radiotherapy/radiology imaging software as a guide for targeted biopsies and treatments. In RT of prostate cancer, dose escalation for prostate cancer has been shown to reduce biochemical failure (30). While dose escalation also has been shown to reduce the need for androgen deprivation in intermediate to high risk patients (4,31,32), when the entire prostate is dose escalated, there is an increased risk of secondary adverse events. Dose escalation only to determinate prostate habitats has the potential reFerences to improve tumor control with less toxicity than when the entire prostate is dose escalated.

eThics sTaTeMenT
This study was carried out in accordance with the recommendations of the Institutional Review Board (IRB) of the University of Miami. The IRB approved a protocol for retrospective review of mpMRI exams from prostate cancer patients. The IRB waived the need for informed consent.
aUThOr cOnTriBUTiOns RS, AP, and NP conceived the study; AP, MA, and SP provided patients; NP, RS, FC, and OK collected and assembled data; NP, RS, BM, FM, RC, and AP conducted data analysis and interpretation; NP, RS, and AP wrote the manuscript; All authors reviewed the manuscript for intellectual content and approved the final version.

FUnDing
This publication was supported by Grants R01CA189295 and R01CA190105 from the National Cancer Institute.