QuantitativeT2: interactive quantitative T2 MRI witnessed in mouse glioblastoma

Abstract. The aim of this study was to establish an advanced analytical platform for complex in vivo pathologies. We have developed a software program, QuantitativeT2, for voxel-based real-time quantitative T2 magnetic resonance imaging. We analyzed murine brain tumors to confirm feasibility of our method for neurological conditions. Anesthetized mice (with invasive gliomas, and controls) were imaged on a 9.4 Tesla scanner using a Carr–Purcell–Meiboom–Gill sequence. The multiecho T2 decays from axial brain slices were analyzed using QuantitativeT2. T2 distribution histograms demonstrated substantial characteristic differences between normal and pathological brain tissues. Voxel-based quantitative maps of tissue water fraction (WF) and geometric mean T2 (gmT2) revealed the heterogeneous alterations to water compartmentalization caused by pathology. The numeric distribution of WF and gmT2 indicated the extent of tumor infiltration. Relative evaluations between in vivo scans and ex vivo histology indicated that the T2s between 30 and 150 ms were related to cellular density and the integrity of the extracellular matrix. Overall, QuantitativeT2 has demonstrated significant advancements in qT2 analysis with real-time operation. It is interactive with an intuitive workflow; can analyze data from many MR manufacturers; and is released as open-source code to encourage examination, improvement, and expansion of this method.


Introduction
The soft tissue contrast in magnetic resonance imaging (MRI) results from the water-filled biological microcompartments. Few MRI relaxation parameters, such as spin-lattice (T1) and spinspin (T2) relaxation times, are designed to measure the distribution of water protons. T2 is specifically sensitive to the extent of water compartmentalization. 1 T2 times decrease as tissue water loses magnetization to surrounding molecules. For example, white matter (WM) has a short T2 time, compared to gray matter (GM) and cerebrospinal fluid (CSF), since WM water is highly compartmentalized by myelin bilayers that promote magnetization exchange. Consequently, it is possible to identify distinct pools of water molecules in biological tissues according to their unique T2 times. When pathological conditions induce changes in tissue structures, significant deviations from normal T2s may occur.
A two-dimensional T2-weighted MR image, which is commonly used in clinical practice, is created by taking a snapshot of the T2 decay for each voxel of the imaging plane.
Quantitative T2 (qT2) imaging is more complex. It requires a Carr-Purcell-Meiboom-Gill (CPMG) imaging sequence 2 that applies multiple refocusing pulses to generate a multiecho T2 image. With sufficient data points, pure T2-weighted decays are computed and further analyzed to determine the combination of T2 times that contributed to the measured decay. 3 This allows scientists to create T2 distributions and discern microcompartmental structures. 4 Multiexponential T2 decay analysis has been applied to study cartilage, 5 blood, 6 muscle, 7 fat, 7 cervix, 7 liver, 8 breast, 9 and prostate. 10 In addition to medical applications, this method has been applied by food scientists to determine if cheese was made using raw or heat-treated milk, 11 and as a quality control tool for fat distributions in fried food. 12 In the brain, qT2 has been used to identify normal WM, 13 GM, 1,14 and tumor. 15 This method has detected myelin content in multiple sclerosis 16,17 and uncovered previously undetected water compartments in human brain pathology in vivo. 18,19 The typical qT2 analysis proceeds as follows: (1) a multiecho T2 MR image is acquired; (2) a region of interest (ROI) is specified on the image; (3) decay data within the ROI are averaged; and (4) a T2 distribution is created by fitting a sum of weighted exponential decays to the averaged decay data. 20 To make the fit more resilient to noise, a regularization routine is often applied.
Although qT2 showed promise for characterizing pathological structures, the following limitations have discouraged widespread clinical adoption of qT2. First, the placement of the ROI requires prior knowledge of where pathology is located. Second, the average signal of an ROI may mask important physiological variability of complex biological tissues. Third, the regularization causes underestimation of the myelin water fraction (MWF), and this error increases as the signal-tonoise ratio decreases. 20 Fourth, the analysis often requires multiple software programs with different interfaces and input requirements. The main objective of this study was to overcome these limitations and establish an advanced analytical platform for complex pathologies, in vivo.
For this purpose, we developed QuantitativeT2, an integrated tool for qT2 analysis. This new tool results from four novel contributions: (1) voxel-based multiexponential qT2; (2) voxelbased T2 distribution histograms; (3) interactive computation of quantitative parametric maps; and (4) numerical assessment of qT2 parameter maps for user-specified ROIs. We hypothesize that our method will improve the existing qT2 analysis. It should identify pathological tissues at early stages of disease with quantitative information on its progression. We analyzed an in vivo mouse model of human malignant glioma (MG) to assess its potentials. Diagnosis of MG is particularly challenging as it is diffuse, highly invasive, and patients rarely survive in the long-term. 21 Animal MG models are commonly used to study the disease progression in a controlled environment. Although increased T2 times have been observed in mouse gliomas using conventional T2-weighted MRI, this method could not detect gliomas at early stages. 22 Multiexponential T2 has been reported in rat glioblastoma, 23 but the contributing T2 values or underlying factors for multiexponential behavior were not explored. Our method allowed probing at the subvoxel level to quantify pathological T2 with specific information on tumor infiltration. We generated quantitative parametric maps that identified water compartment alterations caused by glioma invasion and compared these results to histology.

Animal Model
The mouse MG model was established by implanting immunocompromised mice with patient-derived brain tumor initiating cells (BTICs), a subpopulation of brain tumor cells that maintain the ability to self-renew, proliferate, and give rise to differentiated daughter cells that repopulate the tumor. 24 BTICs, when implanted into the brains of severe combined immunodeficiency (SCID) mice, result in highly invasive tumors that closely resemble MG in humans. 21,24,25 In this study, six-weekold female SCID mice (Charles River Laboratories, Ontario, Canada) were housed in groups of three and maintained on a 12-h light/dark schedule at a temperature of 22 AE 1°C and a relative humidity of 50 AE 5% for a period of 14 weeks. Food and water were available ad libitum. For the experimental group (n ¼ 3, malignant glioma or MG group), the MG tumor model was established on day 3 by inoculating each mouse intracerebrally, in the right putamen with 1 × 10 4 BTICs collected from human surgical specimens. The control group (n ¼ 3, control or C group) contained three, age-matched, normal SCID mice. All animal and human tissue protocols were reviewed and approved by University of Calgary Animal Care Committee.

MRI Protocol
All mouse brains (C and MG groups) were scanned in vivo on day 93 following the BTIC inoculation. The MR imaging was carried out using a 9.4 Tesla (T) Bruker Avance system (Bruker, Billerica, Massachusetts). Four axial slices were imaged from each brain using a CPMG sequence. With a repetition time of 3000 ms, 128 echoes were collected with 5.5 ms echo spacing and 4 averages. The slices were 0.75 mm thick, covering a 1.92 cm 2 field of view. The image matrix contained 128 × 128 pixels. Each pixel represented signals from a voxel of 0.15 × 0.15 × 0.75 mm volume. MR data were saved in raw data format (8 bit unsigned). For each animal, the MR slice containing the largest cross-section of the mouse brain was selected for further analysis.

T2 Decay Analysis
Multiexponential decay analysis was carried out following welldescribed techniques 26,27 applied most recently by our group 20 for ROI-based analysis. Briefly, if one assumes that the T2 decay within a particular voxel is multiexponential, and can be decomposed into a summation of monoexponential functions, then the signal, y, measured at each of N CPMG echoes can be modeled as (1) where t i are the CPMG echo times, M is the number of T2 bins used to model the T2 decay (described in a subsequent paragraph), and s j are the relative weightings for each monoexponential function. Contrary to the conventional practice, we solved this equation to determine individual s j weightings-once for each voxel in the MR slice. This intensive computation was executed in C language by representing Eq. (1) as an N × M matrix [with N rows corresponding to the echo measurement times (t i ) and M columns corresponding to the T2 bins] as shown in Eq. (2). With preset T2 bins, the unknowns, s j , were solved using iterative non-negative least-square (NNLS) algorithm by minimizing error between the model and the measured echoes. : (2) In this study, the echo magnitude reached the noise floor at the 97th echo. Therefore, N was set to 96 and echoes 97 to 128 were not used to avoid analysis artifacts that can be caused by the Rician noise floor. 28 In addition, M ¼ 120 logarithmically spaced T2 bins ranging from 8.25 ms (1.5× shortest echo time) to 1056 ms (2× longest echo time) were used to model the T2 decay in each voxel. The T2 − s j combinations for every voxel were stored independently as T2 distributions, and also summed together to create a T2 distribution histogram for the entire MR slice.

Quantitative Parametric Mapping
A T2 distribution is evaluated by gmT2 and WF measures. The gmT2 is mean T2 time on a log scale. 29 The graphic user interface (GUI) of QuantitativeT2 is equipped with two slider bars that allowed selection of a T2 range in the T2 distribution histogram. For a range between T2 low and T2 high , gmT2 of each voxel for the cross-sectional image was calculated using the following formula: (3) The gmT2 values were arranged in a parametric map, superimposed on a gray-scale T2 scan and displayed on upper left panel (Figs. 1 and 2). WF value indicates the relative signal strength of a T2 range with respect to the entire T2 distribution. 1 QuantitativeT2 computed WF values individually for each voxel using the following formula: where T2 low and T2 high defined the selected T2 range, and T2 min and T2 max were the minimum and maximum T2 times in the T2 distribution histogram. Areaðs j Þ referred to the T2 distribution histogram area, made up of discrete s j weightings, between the defined T2 times. The WF map was then superimposed on a gray-scale T2 scan and displayed in upper right panel ( Figs. 1 and 2). The gray-scale scans are used for anatomical reference. We implemented flexible color mapping and window/leveling algorithms for improved visualization. A sketch function was implemented for both parametric maps for the user to draw an arbitrary ROI for further analysis. A local T2 distribution histogram was computed from the T2 data only within the user-defined ROI. These distributions were scaled for display such that the area under the distribution curve was 1 and represented 100% of the signals.

Histological Analysis
The animals were immediately sacrificed following MRI. The brain from each animal was collected, fixed with 4% paraformaldehyde, paraffin-embedded, and cut into 5-μm-thick slices.
For each mouse, a whole brain tissue section corresponding to the MR scan segment was stained for nestin, using an antibody against this human neural cell progenitor protein 25 (Sigma-Aldrich, Oakville, Ontario, Canada). Hematoxylin was used as a nuclear counter stain. A second tissue section within the MR scan slice was stained with hematoxylin and eosin (H&E) (EMD Chemicals, Gibbstown, USA; Sigma-Aldrich, Oakville, Ontario, Canada). The quantitative analysis of the nestin stained histology sections was performed using the ACIS III automated cellular imaging system. It consisted of five steps: (1) the histology slide was digitized using ACIS III software; (2) the digital section was compared against the segmented gmT2 map of the same mouse, and histological regions corresponding to the T2 bands were identified; (3) within the histological region for each T2 band, four circular ROIs, each with a diameter of 500 μm, were specified using ACIS III; (4) the percentage fraction of brown area (nestin stained, tumor cells) was calculated for each ROI; and (5) the average percentage was computed for the four ROIs corresponding to each T2 band. Data were collected from three control mice and three experimental mice by repeating the same process.

QuantitativeT2 Architecture and Workflow
QuantitativeT2 was developed and tested on an Apple Mac Pro running Mac OS X 10.6. The Mac Pro was equipped with 3 GB of 1066 MHz DDR3 ECC SDRAM and powered by two 2.93 GHz Quad Core 45-nm Xeon W3540 processors. QuantitativeT2 used the C programming language for the computationally intensive voxel-based qT2 and the Objective C programming language for data management. The GUI was developed using Objective C and Interface Builder. The software components were organized according to the model view controller (MVC) architectural pattern and were compiled by Xcode (version 3.2). MVC isolated the application logic from data input and presentation. This allowed independent development and maintenance of the data processing and user interface units. Our matrix-fitting algorithm was written in the C programming language; it incorporated the open-source NNLS code developed and validated by Deguet. 30 Analysis in QuantitativeT2 followed a simple three-step approach.
Step 1: the system performed voxel-based qT2 decay analysis for all voxels within the brain region. This produced a T2 distribution histogram, or a signature, that illustrated the general T2 properties of the tissue cross-section. The pathological T2 values may be visible in the distribution.
Step 2: the user selected a T2 range from the histogram using slider bars. The system computed and displayed gmT2 and WF maps for the selected range in real time. These maps helped localize regions with pathological T2.
Step 3: the user specified an ROI guided by the gmT2 map values and the extent of heterogeneity on the WF map. The mean and standard deviation of WF and gmT2 data within the ROI were then computed and saved.

Results
QuantitativeT2 has introduced a comprehensive voxel qT2 analysis method, by combining novel and existing features, which is intuitive and effective for disease diagnosis. Figure 3 shows proton density (PD) weighted and T2-weighted scans of C1 and MG1 as obtained from MR scanner. These scans are commonly used in clinical practice. Although informative, it is difficult to interpret these scans and isolate pathological tissues. Figures 1 and 2 show the analysis by QuantitativeT2 from raw MR data to parametric maps and numerical results for C1 and MG1. Figure 4(a) displays the T2 distribution histogram computed from C1. Most of the signal was concentrated in a single peak centered at 47.6 ms with tightly packed values in the 30.4 to 60.8 ms range. A second peak corresponding to CSF was observed at 1056 ms. Figure 5(

Tumor Segmentation by T2-Specific Mapping
Figures 4(b) and 4(c) display the gmT2 and the WF maps computed for the low T2s for mouse C1. The gmT2 map selectively displayed tissues corresponding to the T2 range specified in Fig. 4(a). It displayed the average T2 value, within the T2    Other regions of the brain corresponding to higher T2 values (outlined in purple) displayed a heterogeneous WF distribution in Fig. 5(e).

Comparison of qT2 with Histology
For validation purpose, the results obtained by QuantitativeT2 were compared against histology. The histology stains of the control group C1 to C3 revealed no tumor cells (data not shown). Figure 8(a) shows that the four color-coded regions in the gmT2 map of the glioma-bearing mouse brain, MG1, correspond with varying tumor cell densities, when visually compared with Fig. 8(b). Here density of brown stain represents tumor cell population. At 40x magnification [Figs. 8(c)-8(f)], the increase in tumor cell population with increased T2 times is striking. The histological analysis also confirmed that the BTICs were initially implanted at the central mass of the pathological tissue (highlighted in orange in the gmT2 map). In Figs. 8(g)-8(j), the nuclei are stained blue and the extracellular matrix stained pink. Increasing tumor cell density was observed with the increase in T2 for the regions corresponding to T2 bands 1, 2, and 3. Degeneration of the extracellular matrix was observed in the regions corresponding to the highest T2 time, T2 band 4 [ Fig. 8(j)]. Similar results were observed in brain sections stained for nestin or H&E from all three mice with glioma.
The histology sections stained for nestin were later quantified and analyzed using the ACIS III automated cellular imaging system for both control and tumor-bearing mice. For each whole brain section, histological regions corresponding to the T2 bands listed above were identified and the average  percentages of brown area (tumor cell cytoplasm) were calculated for these regions. Figure 9 demonstrates the nestin stain percentages in normal brains for the control groups (C1 to C3) and BTIC tumor-bearing brains for the experimental group (MG1 to MG3). The T2 time increased, from 30.4 to 84.3 ms as the nestin stain% increased, and therefore, the tumor cell density increased. However, nestin% decreased in brain regions with the longest T2 times (87.8 to 149.2 ms). The data trends of nestin percentages were very strongly correlated (p < 0.001) according to Spearman's rank correlation coefficient 31 among the MG tumor-bearing brains, MG1, MG2, and MG3.

Discussion and Conclusion
Our voxel-based analysis method has enhanced the spatial resolution of qT2 analysis compared to ROI-based approaches. It now permits visual interpretation of the general T2 properties of all brain tissues in the scanned slice. Long T2 times are readily identified in the histogram, allowing easy localization of the corresponding anatomical regions. It has enabled the identification of T2 ranges of interest exhibiting subtle differences-these ranges would likely remain unidentified using the traditional ROI-based approaches. QuantitativeT2 has combined interactive computation and presentation of gmT2 and WF parametric maps. Two different properties-water microcompartment structure by gmT2 and partial voluming among these compartments by WF-are now revealed simultaneously by this arrangement. We have also improved the contrast and appearance of parametric maps by flexible color mapping and window leveling algorithms. The local ROI analysis tool has been effective for the isolation and characterization of specific tissues by gmT2 and WF measurements. Overall, we have improved the interpretation of qT2 results with a faster, easier, and intuitive workflow.
QuantitativeT2 is free and released as an open-source code. It can process the raw MR data from many manufacturers. It integrates mathematical processing, visual displays, and interactive analysis functionalities in one platform, such that no other software is required for data processing. It is also time-efficient; it completes qT2 analysis of 128 × 128 scan matrix in <2 s on a desktop/laptop computer running Mac OS X 10.6 or newer. Yoo et al. recently described a method for qT2 analysis that leveraged multicore CPU and graphical processing units to reduce computation time. 32 Their work is highly complementary, and different from ours, in several important ways. Their work is geared toward rapid, noninteractive, batch processing of multiple exams and is focused on production of static MWF maps. Our tool is designed for a dynamic workflow and interactive mapping of gmT2 and WF.
In the literature, qT2 MRI has been shown effective for analyzing multiple neurological conditions, such as tumor, 15 multiple sclerosis, 16 and phenylketonuria. 19 In the current study, we have analyzed MG, one of the most common primary central nervous system tumors in adults. In conventional MRI, a long T2 is observed in regions of increased tissue water and blood volume, as well as in areas where brain tissue texture has been lost from gliosis. 33 In a rat glioblastoma model, both short T2 (20.7 AE 5.4 ms) and long T2 (76.4 AE 9.3 ms) were observed in tumor ROIs, whereas a single T2 component (48.8 AE 2.3 ms) was observed for normal GM. 23 In the latter study, the ROI-based approaches used resulted in coarser spatial resolution and reduced detail in T2 histograms. In contrast, our voxel-based analysis method has shown that the T2 times between 30 and 150 ms are related to cellular density and the integrity of the extracellular matrix. The genotypic and phenotypic changes in MG are expected to result in high cellular density and matrix degeneration, 34 which is in accordance with our results. qT2 mapping has been recommended for improved monitoring of glioblastoma patients under therapy. 35 The analytical advancement from the conventional MRI to our analysis approach is clearly noticeable by comparing Fig. 3 with Fig. 5. The PD and T2 MRI in Fig. 3, which are applied in clinical settings, identify the pathological region but fail to isolate the tumor from normal tissues. Gray-scale intensity is the only information available for further diagnosis. For the same MR slice, voxel-based multiexponential MRI clearly delineates the tumor and provides quantitative information on tumor infiltration. The gmT2 and WF maps in Fig. 5 provide anatomical information specific to the progression of MG, which may assist in cancer staging. Likewise, QuantitativeT2 may also be useful for analyzing other neurological conditions and qT2 MRI data in general.
Due to technical constraints on data availability for human brain pathology, our method has been validated in a mouse model of MG up until the present date. The commonly used qT2 sequence scans one slice at a time, which is inadequate for patient screening. Currently available human qT2 data are acquired in research facilities and not clinical centers, which lack infrastructure to obtain the time-consuming multiecho sequences from all patients. Recently, promising progress has been made using new pulse sequences that accelerate multiecho Table 2 The average gmT2 and WF computed for four ROIs corresponding to the four T2 bands in segmented T2 distribution histograms of the malignant glioma (MG) induced mouse brain scans MG1 to MG3 in Figs. 7(a)-7(c).
data acquisition, while providing information from contiguous slices. 36,37 Combining the scans of multiple slices, these new sequences also produce significantly large data. The new pulse sequences may encourage the clinical application of qT2 only if there is an efficient analysis tool-QuantitativeT2-to process the data and visually integrate the qT2 analysis results. In fact, QuantitativeT2 has been validated against existing qT2 analysis methods and normal human data. 38 One limitation of MRI is its inherent low sensitivity level for identifying individual tumor cells. It can be argued that a very low concentration of tumor cells can remain unnoticed by voxelbased qT2 or T2-based segmentation as observed in the corpus callosum of the contralateral cerebral hemisphere in Fig. 8(b). This was due to the overwhelming effect of the normal tissue for a voxel volume (0.15 × 0.15 × 0.75 mm). However, the heterogeneous WF map of this region in Fig. 5  of WF for the prominent tissue and highlights the importance of further analysis. Perhaps the quantitative ROI characterization [shown in Fig. 5(f)] would provide more information on the tumor cells in that region. The strong correlation among the MG data nestin percentages proves it beyond doubt that the MG tumor infiltration follows a unique behavior in all three mice presented here as shown in Fig. 9. The tumor infiltration is related to T2 alterations, which are detectable by voxel-based qT2 analysis, when the tumor cells are detectable within MRI hardware capacity.
In conclusion, we developed QuantitativeT2-a software package for simple and real-time voxel-based qT2 analysis. QuantitativeT2 brings significant advancement to existing qT2 MRI technology. It provides a novel and interactive workflow for identification and characterization of T2 distribution histograms along with gmT2 and WF parameter maps. With it, we identified structural changes in glioma progression in relation to T2 variations. Being open-source, QuantitativeT2 should promote the use of qT2 analysis in research as well as in diagnostics. 30