Glioblastomas are aggressive brain tumors with dismal prognosis. One of the main bottlenecks for developing more effective therapies for glioblastoma stems from their histologic and molecular heterogeneity, leading to distinct tumor microenvironments and disease phenotypes. Effectively characterizing these features would improve the clinical management of glioblastoma. Glucose flux rates through glycolysis and mitochondrial oxidation have been recently shown to quantitatively depict glioblastoma proliferation in mouse models (GL261 and CT2A tumors) using dynamic glucose-enhanced (DGE) deuterium spectroscopy. However, the spatial features of tumor microenvironment phenotypes remain hitherto unresolved. Here, we develop a DGE Deuterium Metabolic Imaging (DMI) approach for profiling tumor microenvironments through glucose conversion kinetics. Using a multimodal combination of tumor mouse models, novel strategies for spectroscopic imaging and noise attenuation, and histopathological correlations, we show that tumor lactate turnover mirrors phenotype differences between GL261 and CT2A mouse glioblastoma, whereas recycling of the peritumoral glutamate-glutamine pool is a potential marker of invasion capacity in pooled cohorts, linked to secondary brain lesions. These findings were validated by histopathological characterization of each tumor, including cell density and proliferation, peritumoral invasion and distant migration, and immune cell infiltration. Our study bodes well for precision neuro-oncology, highlighting the importance of mapping glucose flux rates to better understand the metabolic heterogeneity of glioblastoma and its links to disease phenotypes.

Animals and cell lines

This study was performed in strict accordance with European Directive 2010/63 and the Portuguese law (Decreto-Lei 113/2013), following the FELASA (Federation of European Laboratory Animal Science Associations) guidelines and recommendations concerning laboratory animal welfare, and aligned with the ARRIVE (Animal Research: Reporting of In Vivo Experiments) guidelines. All animal experiments were performed at the Champalimaud Foundation Vivarium under project #05318, pre-approved by the competent institutional and national authorities: ORBEA (Champalimaud Foundation Animal Welfare Body) and DGAV (Direcção Geral de Alimentação e Veterinária), respectively. All the surgeries were performed under isoflurane anesthesia, and every effort was made to minimize suffering. A total of n=10 C57BL/6j male mice were used in this study, bred at the Champalimaud Foundation Vivarium, and housed with ad libitum access to food and water and 12h light cycles. GL261 mouse glioma cells were obtained from the Tumor Bank Repository at the National Cancer Institute (Frederick MD, USA): “GLIOMA 261”, sample number 0507815. CT2A mouse glioma cells were kindly provided by Prof. Thomas Seyfried at Boston College (Boston MA, USA). Both cell lines were grown in RPMI-1640 culture medium supplemented with 2.0 g/l Sodium Bicarbonate, 0.285 g/l L-glutamine, 10% Fetal Bovine Serum (Gibco) and 1% Penicillin-Streptomycin solution. The cell lines tested negative for mycoplasma contamination using the IMPACT Mouse FELASA 1 test (Idexx-BioResearch, Ludwigsburg, Germany).

Glioma models

Tumors were induced in previously described [66]. Briefly, intracranial stereotactic injection of 1 x10⁵ GL261 or CT2A cells was performed in the caudate nucleus (n=5 and n=5 mice, respectively); analgesia (Meloxicam 1.0 mg/Kg s.c.) was administered 30 min before the procedure. Mice were anesthetized with isoflurane (1.5-2.0% in air) and immobilized on a stereotactic holder (Kopf Instruments, Tujunga/CA, USA) where they were warmed on a heating pad at 37 ºC, while body temperature was monitored with a rectal probe (WPI ATC-2000, Hitchin, UK). The head was shaved with a small trimmer, cleaned with iodopovidone, and the skull exposed through an anterior-posterior incision in the midline with a scalpel. A 1 mm hole was drilled in the skull using a micro-driller, 0.1 mm posterior to the bregma and 2.32 mm lateral to the midline. The tumor cells (1x10⁵ in 4 μL PBS) were inoculated 2.35 mm below the cortical surface using a 10 µL Hamilton syringe (Hamilton, Reno NV, USA) connected to an automatic push-pull microinjector (WPI Smartouch^TM, Sarasota FL, USA), by advancing the 26G needle 3.85 mm from the surface of the skull (~1mm skull-to-brain surface distance), pulling it back 0.5 mm, and injecting at 2 μL/min rate. The syringe was gently removed 2 min after the injection had finished, the skin sutured with surgical thread (5/0 braided silk, Ethicon, San Lorenzo Puerto Rico) and wiped with iodopovidone. During recovery from anesthesia, animals were kept warm on a heating pad and given an opioid analgesic (Buprenorphine 0.05 mg/Kg s.c.) before returning to their cage. Meloxicam analgesia was repeatedly administered at 24- and 48-hours post-surgery.

In vivo Studies

Longitudinal MRI

GBM-bearing mice were imaged every 5-7 days on a 1 Tesla Icon MRI scanner (Bruker BioSpin, Ettlingen, Germany; running ParaVision 6.0.1 software), to measure tumor volumes. For this, each mouse was placed in the animal holder under anesthesia (1-2 % isoflurane in 31% O₂), heated with a recirculating water blanket, and monitored for rectal temperature (36-37 ºC) and breathing (60-90 BPM). Tumor volume was measured with T2-weighted ¹H-MRI (RARE sequence, ×8 acceleration factor, repetition time TR = 2500 ms, echo time TE = 84 ms, 8 averages, 1 mm slice thickness, and 160×160 µm² in-plane resolution), acquired in two orientations (coronal and axial). Each session lasted up to 30 min/animal. 

End-point MRI and DMI

GBM-bearing mice with tumors ≥35 mm³ (longitudinal MRI assessment) were scanned on a 9.4T BioSpec MRI scanner (Bruker BioSpin, Ettlingen, Germany; running under ParaVision 6.0.1), using a ²H/¹H transmit-receive surface coilset customized for the mouse brain (NeosBiotec, Pamplona, Spain), as described before [31]. Before each experiment, GBM-bearing mice fasted 4-6h, were weighed, and cannulated in the tail vein with a catheter connected to a home-built 3-way injection system filled with: 6,6′-²H₂-glucose (1.6M in saline); Gd-DOTA (25 mM in saline); and with heparinized saline (10 U/mL). Mice were placed on the animal holder under anesthesia (as in 2.3.1). Coilset quality factors (Q) for ¹H and ²H channels were estimated in the scanner for each sample based on the ratio of the resonance frequency (400.34 and 61.45 MHz, for protons and deuterium, respectively) to its bandwidth (full width at half-minimum of the wobbling curve during the initial tuning adjustments): 175±8 and 200±12, respectively. Mice were imaged first with T2-weighted ¹H-MRI (RARE sequence, x8 acceleration factor, 3000 ms TR, 40 ms TE; 2 averages, 1 mm slice thickness, 70 µm in-plane resolution) in two orientations (coronal and axial). Then, the magnetic field homogeneity was optimized over the tumor region based on the water peak with ¹H-MRS (STEAM localization: 6x6x3 mm volume of interest, i.e. 108 µL) using localized 1^st and 2^nd order shimming with the MapShim Bruker macro, leading to full widths at half-maximum (FWHM) of 28±5 Hz.

DMI was performed using a slice-FID chemical-shift imaging pulse sequence, with 175 ms TR, 256 spectral points sampled over a 1749 Hz window, and Shinnar-Le Roux RF pulse [67, 68] (0.42ms, 10kHz) with 55º flip angle, to excite a brain slice including the tumor: 18×18 mm field-of-view, and 2.27 mm slice thickness. After RF pulse calibration (using the natural abundance semi-heavy water peak, DHO), DGE-DMI data were acquired for 2h23min (768 repetitions), with i.v. bolus of 6,6′-²H₂-glucose (2 mg/g, 4 µL/g injected over 30 s; Euroisotop, St Aubin Cedex, France). Data were sampled with an 8×8 matrix and 4-fold Fourier interpolated [69], rendering a 560 µm in-plane resolution. A reference T2-weighted image was additionally acquired with matching field-of-view and slice thickness, and 70 µm in-plane resolution.

Finally, animals underwent DCE T1-weighted ¹H-MRI (FLASH sequence, 8º flip-angle, 16ms TR, 4 averages, 150 repetitions, 1 slice with 140 µm in-plane resolution and 2.27 mm thickness, FOV size and position matching the DGE-DMI experiment), with i.v. bolus injection of Gd-DOTA (0.1 mmol/Kg, injected over 30 s; Guerbet, Villepinte, France). Animals were then sacrificed, brains were removed, washed in PBS, and immersed in 4% PFA.

MRI/DMI Processing

T2-weighted ¹H-MRI

T2-weighted MRI data were processed in ImageJ 1.53a (Rasband, W.S., ImageJ, U. S. National Institutes of Health, Bethesda, Maryland, USA, https://imagej.nih.gov/ij/, 1997-2018). For each animal, the tumor region was manually delineated on each slice, and the sum of the areas multiplied by the slice thickness to estimate the volume, which was averaged across the two orientations acquired (coronal and axial).

DGE-DMI

DGE-DMI data were processed in MATLAB^® R2018b (Natick, Massachusetts: The MathWorks Inc.) and jMRUI 6.0b [70]. Each dataset was averaged to 12 min temporal resolution and noise regions outside the brain, as well as the olfactory bulb and cerebellum, were discarded, rendering a 4D spectral-spatial-temporal matrix of 256×32×32×12 points. After automated phase-correction of each spectrum, the 4D matrix was denoised with  a tensor PCA denoising approach [32]. For this, a [8 8 8] window and tensor structure [1 2:3 4] were used for patch processing the spectral, spatial, and temporal dimensions with, whereas the a priori average standard deviation of the noise in each spectrum (calculated σ²) was used to avoid deleterious effects of spatially-correlated noise [34]. Then, these denoised spectra were analyzed voxel-wise by individual peak fitting with AMARES (similarly to the single-spectrum analysis reported previously in [52]), using a basis set for DHO (4.76 ppm: short- and long-T2 fractions [18]) and deuterium-labelled: glucose (Glc, 3.81 ppm), glutamate-glutamine (Glx, 2.36 ppm), and lactate (Lac, 1.31 ppm); relative linewidths referenced to the estimated short-T2 fraction of DHO, according to the respective T2 relaxation times reported by de Feyter et al [18]. The natural abundance DHO peak (DHOi) was further used to select and quantify both original and denoised spectra: SNR_DHOi >3.5 and 13.88 mM reference (assuming 80 % water content in the brain and 0.03 % natural abundance of DHO), respectively. Metabolite concentrations (CRLB<50%; otherwise discarded) were corrected for T1 and labeling-loss effects, according to the values reported by de Feyter et al (T1, ms: DHO, 320; Glc, 64; Glx, 146; Lac, 297) [18] and de Graaf et al (number of magnetically equivalent deuterons: DHO, 1; Glc, 2; Glx, 1.2; Lac, 1.7) [71], respectively. Thus, the concentration of each metabolite (m) at each time point was estimated as (Eq 1):

(Eq 1)

Area = peak area; Area0 = average peak area before injection; d = number of magnetically equivalent deuterons corrected for labelling-loss effects; C = T1 correction factor (1-exp(-TR/T1)); and Conc_ref = reference DHO concentration.

The time-course changes of ²H-labelled metabolite (Glc, Glx and Lac) concentrations were fitted using a modified version of the kinetic model reported by Kreis et al [19], to estimate the maximum rate of Glc consumption (total, V_max) for Glx synthesis (mitochondrial oxidation, V_glx) and Lac synthesis (glycolysis, V_lac), and the confidence intervals for all estimated parameters:

(Eq 2)

The coupled differential equations describing the concentration kinetics of each metabolite were:

(Eq 3)

(Eq 4)

                                                   (Eq 5)

where: k_g, apparent rate constant of glucose transfer between blood and tumor (min⁻¹); k_glx, apparent rate constant of Glx elimination (min⁻¹); k_lac, apparent rate constant of lactate elimination (min⁻¹);  Glc concentration in plasma (mM); a₁, the Glc concentration after the bolus injection (mM); and k_p, the effective rate constant of labeled glucose transfer to tissue (min⁻¹). As reported previously [31], the following parameters were fixed: fraction of deuterium enrichment (f), at 0.6 [19]; constant for glucose uptake (k_m), at 10 mM [72, 73]; and the extravascular-extracellular volume fraction (v), at 0.22 – average estimation from DCE-T1-weighted MRI analysis (Table 1 - Supplementary file 1a). All the other parameters were fitted without any restrictions to their range. Metabolic rate maps were displayed and analyzed pixel-wise using cut-off points defined by 5 times their respective confidence intervals.

DCE T1-weighted MRI

DCE T1-weighted MRI data were processed with DCE@urLab [74], as before [31]. First, ROIs were manually delineated for each tumor and the time-course data was fitted with the Extended Tofts 2-compartment model [75], to derive the volume transfer constant between plasma and tumor extravascular-extracellular space (k^trans), the washout rate between extravascular-extracellular space and plasma (k_ep), and the extravascular-extracellular volume fraction (v_e). Then, each dataset was reprocessed by down-sampling the original in-plane resolution to match the DGE-DMI experiment (0.56×0.56×2.27 mm³), and fitting the time-course data pixel-wise with the Extended Tofts 2-compartment model to derive k^trans, k_ep, and v_e maps (pixels with root-mean square error >0.005 discarded).

Histopathology and Immunohistochemistry

Whole brains fixed in 4% PFA were embedded in paraffin and sectioned at 30 different levels on the horizontal plane, spanning the whole tumor area. 4 µm sections were stained with H&E (Sigma-Aldrich, St. Louis MO, USA), digitized (Nanozoomer, Hamamatsu, Japan), and analyzed by an experimental pathologist blinded to experimental groups, according to previously established criteria [31]. Then, QuPath v0.4.3 built-in tools [76] were used to highlight different tumor regions: Tumor ROIs, corresponding to the bulk tumor, were delineated first with “create threshold” and then manually corrected; P-Margin ROIs, including areas of peritumoral infiltration, were delineated with “expand annotations” by expanding 100 µm the tumor margin toward the adjacent brain parenchyma; Infiltrative ROIs, corresponding to specific infiltrative regions, were manually annotated. Between 3 to 6 sections of each tumor were also immunostained for Ki67 (mouse anti-ki67, BD, San Jose CA, USA; blocking reagent, M.O.M ImmPRESS kit, Vector Laboratories, Burlingame CA, USA; liquid DAB⁺, Dako North America Inc, Carpinteria CA, USA), digitized (Nanozoomer, Hamamatsu, Japan), and analyzed with QuPath built-in tools [76] for Tumor and P-Margin ROIs, defined as detailed above. Thus, Ki67^+/- cells were counted semi-automatically to determine the total number of cells, the cell density, and the proliferation index (% Ki67⁺ cells) as the average across slices for each ROI, and respective Tumor/P-Margin ratios. This procedure was repeated for each animal. In addition, one histologic section corresponding to each DGE-DMI slice was immunostained for Iba-1 (rabbit anti-Iba-1, Fujifilm Wako PCC, Osaka, Japan; NovolinkTM Polymer, Leica Biosystems, UK; liquid DAB+, Dako North America Inc, Carpinteria CA, USA), digitalized (Philips UFS v1.8.6614 slide scanner) and analyzed in QuPath. Tumor region and peritumoral margin regions were automatically annotated as outlined above, and Iba-1 positive staining was quantified across all annotations using the threshold tools, adjusted for each slide to account for variations in staining intensity, to calculate the percentage of Iba-1 positive area: (Iba-1+ area / total annotation area) * 100.

Statistical analyses

Data were analyzed in MATLAB^® R2018b (Natick, Massachusetts: The MathWorks Inc.) using the two-tailed Student’s t-test, either unpaired (comparing different animal cohorts) or paired (comparing the same animal cohort in different conditions). Differences at the 95% confidence level (p=0.05) were considered statistically significant. Correlation analyses were carried out with the Pearson R coefficient. Error bars indicate standard deviation unless indicated otherwise.