Glioblastoma (GBM) is the most common and aggressive intracranial tumor in adults. Despite decades of research, survival remains dismal, with a mean survival time of 15 months and a mean 5-year survival rate of <10% (1–3). The standard therapy consists of maximal safe surgical resection, followed by radiation and then chemotherapy with temozolomide, and this protocol has changed very little since its inception (2). Although many patients initially respond, GBM inevitably recurs in most patients (1).
Given its poor prognosis, GBM is a perfect candidate for targeted therapy using precision medicine; thus, understanding the genetic profile of GBM is vital. However, known spatial and temporal tumor heterogeneity is recognized as a significant obstacle for accurate diagnosis impeding the use of targeted therapies. GBM is among the first tumor types for which in-depth genomic analyses were performed (4), and integrated genomic analyses performed in a large cohort of patients uncovered complex clonal heterogeneity (4), intertumoral heterogeneity (5–10), and intratumoral heterogeneity (11–14). Recent studies have revealed that customized treatments targeting specific mutations may not be efficacious in altering patient outcomes (15). Although in this case the cause of treatment failure is unclear, previous studies have shown that between 67% and 91% of posttreatment failures no longer harbored the target mutation (16, 17). This is probably because of clonal expansion of cells outside the targeted mutation. GBM is also known to change its molecular signature with treatment (18). These findings reflect the highly dynamic genetic evolution present in GBM, confounding the elucidation of its complete genetic landscape.
Multiregional tumor sampling is ideal in GBM and other tumor types, but it is often impractical and potentially risky for routine clinical application (19, 20). Furthermore, tumor sampling during different treatment phases is often clinically unindicated. However, precise mechanisms underlying tumor heterogeneity and its functional significance can be evaluated in preclinical models for discovery of driver mutations, discovery of mechanisms of drug resistance, identification of genetic and imaging biomarkers, and development of clone-specific cancer therapies.
Here, we devised and successfully showed a magnetic resonance imaging (MRI)-guided stereotactic biopsy method for isolation and functional profiling of tumoral heterogeneity to assess temporal and geographical genomic variations during disease progression in an intact murine GBM model. Longitudinal tumor biopsy samples were obtained at different stages of tumor growth and in spatially distinct locations. RNA sequencing (RNA-Seq) of biopsy tissue showed unique temporal and spatial gene expression patterns.
The p53 mutant GBM mouse model was obtained from Dr. Yuan Zhu Children's National Medical Center (Washington DC). The genotype and maintenance of the strain have been described previously (21). Animal experiments were approved by the University Committee on the “Use and Care of Animals” (UCUCA) at the University of Michigan.
Mice were anesthetized using a 1%–2% isofluorane–air mixture throughout the duration of the MRI experiment. Then, 5 minutes before anesthesia administration, 50 μL of 0.5M gadolinium (Gd)–diethylenetriaminepentaacetic acid (DTPA) (Magnevist; Bayer Healthcare Pharmaceuticals) was injected intraperitoneally. Once anesthetized, mice were securely fastened in the animal handling MRI cradle with a bite bar and stererotactic bilateral ear bars. An external fiducial marker consisting of Gd contrast agent was secured to the coil ∼1 cm left of the mouse skull. MRI scanning was performed using a 7 T Direct Drive MR system (Agilent, CA) using a volume radiofrequency transmit coil and a 2-channel fixed-tuned actively decoupled surface receive coil with the mouse handling cradle including bilateral ear bars and a bite bar for securing the mouse head (Rapid MR International). During MRI, heated air was introduced into the magnet bore (World Precision Instruments) to maintain animal body temperature. T1-weighted spin echo images were acquired for each animal using a spin echo sequence with the following parameters: repetition time/echo time = 510/15 milliseconds, field of view = 20 × 20 mm2, matrix size = 128 × 128, slice thickness = 500 μm, number of slices = 25 interleaved sections, and 2 signal averages per phase encode step.
MRI Biopsy Validation Study
Phantoms were constructed using agarose gel blocks. Bubbles with varying volumes of air (5, 10, and 20 μL) were injected into the gel before it solidified using a Hamilton syringe. The MRI radiofrequency coil bed with an associated fiducial marker was modified to allow for secure attachment to a Kopf stereotactic device to provide for a consistent and reproducible coordinate system. MRI-guided biopsy was performed on the gel blocks as described below to validate targeting accuracy. Instead of extracting “tissue” (eg, gel), we injected a 1:10 dilution of 0.5M Gd–DTPA (Magnevist) into the agarose gel at the calculated location of the intended target volume (eg, air bubble). If the air bubble was successfully filled with the contrast agent, the air would be displaced and a bright signal in the T1 magnetic resonance (MR) image would be used to verify success versus failure. Gel-based target biopsies were performed on 15, 12, and 12 bubbles with 5-, 10-, and 20-μL volumes, respectively. The accuracy of the method was determined as follows: [(success)/(total number of target volumes)] × 100%. Success was determined as the ability of the biopsy needle to completely fill the target volume as determined by MRI.
Intracranial Biopsy Procedure
After acquisition of prebiopsy MR images, the mouse was anesthetized with a ketamine/xylazine (0.1/0.02 mg/kg) mixture. The tumor was localized with 3-dimensional (3D) Cartesian coordinates on the basis of the distance from the external fiducial marker. The coil was secured to a coil-holder with an attached Kopf stereotactic coordinate system. The mouse was prepped for surgical biopsy and the eyes were protected using eye lubrication. The scalp was sterilized with a topical antiseptic solution and allowed to dry. A 1.5-cm longitudinal incision was made across the midline beginning just posterior to the nasion. The scalp and pericranium were gently retracted and the skull was exposed with cotton-tip applicators. A drill was attached to the stereotactic system (Foredom K-1070 with MH-170 drill with stereotaxic attachment-Model 1474; with Kyocera #70 105-0280.400 drill bits). The attached drill was used to localize the skull burr hole site based on the MRI-determined Cartesian coordinates. After the burr hole was made, a fine aspiration needle (Inrad Model 54722) was secured to the stereotactic instrument and placed into the burr hole. The needle was then advanced to the MRI-determined depth at the proximal/superior edge of the tumor. The inner cannula was removed and the needle was advanced toward the center of the tumor. The needle was attached to a mini-bore extension set (Abbott 32072) connected using a locking aspiration syringe. A syringe was used to withdraw with 2.5 mL of pressure for 3 seconds before the pressure was re-equalized, and the needle was withdrawn from the skull. The biopsy sample was collected in a sterile 2-mL centrifuge tube. Samples were either snap frozen or manually dissociated and cultured at 37°C.
After the biopsy, the burr hole was covered with bone wax (#W810, Ethicon), and the incisions were approximated and closed with Vetbond (#1469SB, Vetbond). To minimize pain, 100 μL (1 mg/mL) of Rimadyl (Pfizer Animal Health) was subcutaneously administered following completion of the surgical procedure. The mouse was kept secured in the coil holder and the MRI was repeated. T1-weighted images were acquired as described above to confirm the accuracy of the biopsy and to rule out significant complications.
A second late-stage biopsy was performed when the tumor volumes reached ∼100 μL. Two distinct biopsy locations were identified and similar biopsy procedures were carried out as described above. The notation Mx[y] was used to denote mouse number [x] and tumor biopsy sample number [y]. At the conclusion of the final biopsy and MR scans to confirm biopsy locations, the animals were sacrificed.
Tumor biopsy tissues were used for RNA extraction. Total RNA was extracted using TRIzolTM (#10296010; Life Technologies) and ethanol precipitation as per the manufacturer's instructions. RNA-Seq library preparation and sequencing was completed by DNA sequencing core at the University of Michigan. Specifically, the RNA library was prepared with NEBNext® rRNA Depletion Kit (#E6310L, NEB). Further, 30–45 million 50nt-single-end reads were obtained from each sample with an estimated 20–30× sequencing depth.
RNA-Seq Data Analysis
Single-end reads of RNA-Seq data were mapped to the mouse transcriptome (GRCm38.p4 /mm10) using BWA, duplicate reads were marked by Picard, and base quality scores were recalibrated using GATK, by referring to the PRADA pipeline (22), followed by proceeding with the mapped reads to calculate the read counts and RPKM (reads per kilobase of transcript per million mapped reads) values for 39,017 expressed genes (Ensembl Gene ID) using UCSC/mm10/ensGene transcript database, with associated gene names annotated based on Ensembl genome assembly (Biomart). Differential expression of genes was calculated by applying quantile-adjusted conditional maximum likelihood method that fits read counts with negative binomial models and obtains dispersion estimates, followed by an exact test to determine differential expressions based on false discovery rate (FDR ∼ adjusted P value) (23). Genes were considered to be differentially expressed when the change was greater than 1.5-fold, and the adjusted P value (FDR) was <0.05%.
Pathway and network analyses of differentially expressed genes in early- versus late-stage biopsies were performed using ingenuity pathway analysis (IPA) software (Qiagen). The molecule activity predictor in IPA was used to predict the upstream or downstream activation or inhibition of a given pathway. The P value of the enrichment score was used to evaluate the significance of the overlap between observed and predicted gene sets. All statistical analyses were performed and graphs were obtained under R (version 3.2.3) environment.
Tumor purity was evaluated on RPKM (reads per kilobase of transcript per million mapped reads) values of gene expression profiling using ESTIMATE, based on the enrichment of gene signatures in stromal and immune cells (24). The molecular classification of samples was performed using ssGSEA.
To enable temporal and geographical multisampling of tumor tissue in a murine GBM model, we devised a biopsy apparatus (Figure 1A). For tumor localization, a 3D Cartesian coordinate diagram with external fiducial markers was used to guide biopsy (Figure 1B). To validate the stereotactic biopsy efficiency and accuracy, an MR phantom consisting of agar gels was developed. Three different volumes of air bubbles (void spaces) were used as representative tumor sizes. Once the locations of the voids were determined by MRI within the 3D coordinate system on the basis of the distance from the external fiducial marker, Gd was stereotactically injected into spatially distinct voids via biopsy needles. MRI T1 images were taken before and after Gd injection to determine the accuracy of the biopsy method (Figure 1C). Further, 12–15 mock biopsies were performed for each size, and the accuracy was shown to be 100, 92, and 73% at volumes of 20, 10, and 5 μL, respectively (Figure 1D). These results show that sampling of tumors as small as 5 μL can be routinely achieved with a high degree of accuracy using this MRI-guided biopsy method.
Next, we demonstrated this method for stereotactic MRI-guided multiregional sampling of spontaneous murine p53 mutant GBM tumors with tumor volumes of 10–20 μL. Animals were secured in a stereotactic biopsy station (Figure 1A), and the tumor was localized on the basis of the distance from the external fiducial marker (Figure 1B) using MRI. MR images were easily acquired both before and after stereotactic biopsy sampling (as depicted in Figure 2A) by transferring the MRI animal cradle to a holder with an attached stereotactic coordinate system. After surgical preparation, the Cartesian coordinate system was used to localize the tumor, and a fine aspiration needle attached to an articulating apparatus was used to withdraw the biopsy tissue specimen. The tissue was either snap frozen for RNA extraction or manually dissociated for cell culture (data not shown). T1-weighted MRI was repeated after biopsy to confirm the accuracy of the biopsy and the sampling location. Three p53 mutant GBM animals were used to perform an early-stage biopsy. Further, 2 geographically distinct locations were sampled during the late-stage biopsies for each animal at a tumor volume of ∼100 μL (Figure 2, B and C) to investigate intertumoral genomic heterogeneity.
The analysis of spatiotemporal genomic heterogeneity was accomplished by RNA-Seq-based expression profiling of early- and late-stage biopsies, which revealed distinct temporal and spatial gene expression patterns. We compared both early- and late-stage biopsies of 3 individual mice, and we depicted the top 10 genes upregulated (Figure 3A) and downregulated (Figure 3B) in late-stage biopsies. Pathway-enrichment analysis of differentially expressed genes during tumor progression further identified unique molecular patterns. Upregulated and downregulated pathways in late- versus early-stage biopsies are depicted (Figure 3C). Overall, the pathway-enrichment analysis of differentially expressed genes in early- versus late-stage biopsy samples showed significant differences in molecular patterns during tumor progression. Moreover, spatially varying differences in expressed genes were also observed in late-stage tumors, revealing that intertumoral heterogeneity in this GBM mouse model could be detected using this method.
Transcriptional profiling of GBM has been previously performed in large-scale clinical studies, and heterogeneous genotypes consisting of subclones with varying tumorigenic, proliferative, differentiation, and drug responsiveness have been shown (4, 5, 18). To identify whether such spatial (geographical) gene expression patterns were present in a preclinical p53 mutant GBM model, 2 spatially distinct locations of late-stage tumors using the aforementioned criteria and IPA pathway analysis were compared (Figure 3D). Interestingly, transcription factor OTX2 was differentially regulated between the 2 spatially distinct regions. OTX2 has been implicated as an oncogene in medulloblastoma where it is highly expressed and amplified in a subset of these tumors (25). OTX2 regulates cell cycle and differentiation, and its regulated genes may serve as potential therapeutic targets (25). Differentially regulated pathways between regions are depicted (Figure 3E). Up- and/or downregulated genes were found, and although detailed studies need to be conducted to investigate the roles of each gene in tumor growth and survival, some genes identified have been indicated to play a role in overall patient survival and tumorigenesis. For example, Ighg1 has been shown to affect overall survival in GBM patients, indicating that it may be associated with more aggressive tumor characteristics in the late stages of tumor growth.
Molecular characterization of tumors to guide treatment selection is hampered by intratumor heterogeneity (19, 26). Longitudinal and multifocal monitoring of clonal dynamics may provide for improved patient care (19). In this study, we devised an experimental approach for image-guided tumor sample acquisition from isolated GBM biopsies for sequencing to enable comparative analyses of gene expression patterns and tumor growth. This method allowed for longitudinal tissue sampling over time and in geographically specific regions to gain further insight into the mechanisms of heterogeneity, thus allowing for linking phenotypic and genomic characteristics.
Although phenotypic analyses are not shown here, neurosphere cultures may also be established from GBM biopsies (data not shown). Single-cell expansion and functional assays may be performed from these biopsies and are a topic for future investigations. Although cell to cell differences were not investigated in this study, this approach would also allow for single-cell isolation and investigation of clonal heterogeneity. These studies can easily be expanded to include DNA extraction and the analysis of genetic variations by exome or deep sequencing approaches. We believe that these types of sequencing analyses may aid in the identification of genetic and imaging biomarkers.
Currently, MRI-based stereotactic intracranial procedures are routine in clinical practice. However, an accurate intracranial biopsy technique has not been shown in mouse models. Here, we established an innovative and accurate stereotactic image-guided murine biopsy method that is capable of accurate and repeated biopsy collection from an individual animal. Tissues collected can be used to investigate genomic dynamics during disease progression and modulation during treatment intervention. Spatial and temporal variations in the molecular landscape of preclinical models may provide unique insights into the development and interrogation of novel therapeutic strategies. This information provides opportunities for evaluating genomic responses during treatment to gain a further biological understanding of drug/target modulation and mechanisms of drug resistance and clonal evolution. Elucidation of the precise mechanisms of existing or emerging tumor heterogeneity and its functional significance in preclinical models provides opportunities to discover driver mutations, emergence of drug resistance, identification of genetic and or imaging biomarkers, and development of clone-specific cancer therapies. Here, we devised a method to isolate and functionally profile intratumoral heterogeneity by using MRI-guided biopsies to assess temporal and geographical differences during disease progression in a murine GBM model. We successfully collected tumor biopsy samples at different stages of tumor growth and individual locations. RNA-Seq results showed unique temporal and spatial gene expression patterns in GBM. This approach may allow for new opportunities for assessing tumor tissue samples from discrete spatial locations over time to improve our understanding of the role of spatial heterogeneity and genomic instability in adaptation mechanisms and therapy resistance. In the future, this approach may be adapted to investigate treatment-induced changes and adapted to study other primary and metastatic tumor types. In summary, the combination of MRI-guided tumor sampling with advanced genome sequencing can facilitate identification of mechanisms leading to heterogeneity and development of clinically useful measures of heterogeneity, which may guide identification of new therapeutic strategies.