Abstract
Purpose: To identify mechanisms and mediators of resistance to antiangiogenic therapy in human glioblastoma.
Experimental Design: We carried out microarray gene expression analysis and immunohistochemistry comparing 21 recurrent glioblastomas progressing during antiangiogenic treatment with VEGF neutralizing antibody bevacizumab to paired pretreatment tumors from the same patients.
Results: Microarray analysis revealed that bevacizumab-resistant glioblastomas (BRG) had two clustering patterns defining subtypes that reflect radiographic growth patterns. Enhancing BRGs (EBRG) exhibited MRI enhancement, a long-established criterion for glioblastoma progression, and expressed mitogen-activated protein kinases, neural cell adhesion molecule-1 (NCAM-1), and aquaporin 4. Compared with their paired pretreatment tumors, EBRGs had unchanged vascularity and hypoxia, with increased proliferation. Nonenhancing BRGs (NBRG) exhibited minimal MRI enhancement but had FLAIR-bright expansion, a newer criterion for glioblastoma recurrence since the advent of antiangiogenic therapy, and expressed integrin α5, laminin, fibronectin1, and PDGFRβ. NBRGs had less vascularity, more hypoxia, and unchanged proliferation than their paired pretreatment tumors. Primary NBRG cells exhibited more stellate morphology with a 3-fold increased shape factor and were nearly 4-fold more invasive in Matrigel chambers than primary cells from EBRGs or bevacizumab-naive glioblastomas (P < 0.05).
Conclusion: Using microarray analysis, we found two resistance patterns during antiangiogenic therapy with distinct molecular profiles and radiographic growth patterns. These studies provide valuable biologic insight into the resistance that has limited antiangiogenic therapy to date. Clin Cancer Res; 18(10); 2930–42. ©2012 AACR.
Initial successes of angiogenesis inhibitors were tempered by their failure to produce enduring anticancer responses. Here, we describe 2 phenotypes of glioblastomas developing resistance to antiangiogenic therapy with distinct molecular and radiographic profiles. The first exhibited enhancing nodular MRI growth, expressed mitogen-activated protein kinases and aquaporin 4, and exhibited unchanged vascularity and hypoxia, and increased proliferation compared with paired pretreatment tumors. The second had nonenhancing poorly circumscribed MRI growth, expressed α5β1 integrin, laminin, and fibronectin, and exhibited less vascularity, more hypoxia, and unchanged proliferation compared with paired pretreatment tumors. The translational significance of these findings is 2-fold. First, our data reveal heterogeneous tumor responses to hypoxia induced by antiangiogenic therapy–-some respond with neovascularization and proliferation, others with perivascular invasion. Second, targeting mediators of these responses that we identified could disrupt resistance to antiangiogenic therapy, allowing these treatments to fulfill their therapeutic promise.
Introduction
Recognition of the role of VEGF in developing the vascularity of glioblastomas, which contributes to their growth and treatment resistance, has led to clinical trials of humanized monoclonal VEGF antibody bevacizumab as monotherapy or combined with DNA damaging agents such as irinotecan in glioblastoma patients (1–4). Two trials showing efficacy of bevacizumab monotherapy (3, 4) led to the 2009 U.S. food and drug administration (FDA) approval of bevacizumab for recurrent glioblastoma, making bevacizumab just the third FDA-approved glioblastoma treatment in nearly 4 decades. Randomized trials stemming from these results are studying bevacizumab in newly diagnosed glioblastomas, potentially allowing bevacizumab to join standard treatment regimens for newly diagnosed and recurrent glioblastomas.
Unfortunately, as with other cancers (5), the response to antiangiogenic therapy in glioblastoma is often transient, with 40% to 60% radiographic progression rates after initially successful bevacizumab treatment in phase II clinical trials (1, 2). Glioblastomas progressing during bevacizumab therapy can exhibit nonenhancing FLAIR-bright growth (6) or restricted diffusion (7) on MRI. Although these growth patterns were initially considered common after antiangiogenic therapy, subsequent analyses have shown them to occur in less than half of glioblastomas progressing during antiangiogenic therapy (6, 8, 9). Thus, imaging after resistance to antiangiogenic therapy suggests heterogeneous resistance mechanisms and illustrates the biology of antiangiogenic therapy resistance, as FLAIR bright nonenhancing growth is thought to represent tumor infiltration, whereas restricted diffusion is believed to represent hypoxia. This pattern of increased hypoxia and invasiveness has also been described in preclinical models of VEGF blockade (10–12). Uncircumscribed growth after antiangiogenic therapy often limits the benefit of surgery (13), and tumor hypoxia after antiangiogenic therapy reduces response to available chemotherapies.
Preclinical studies suggest that tumors become resistant to antiangiogenic therapy by transcriptional reprogramming allowing tumor cells to grow while the antiangiogenic target remains inhibited (14). This paradigm represents a departure from resistance to traditional DNA damaging chemotherapy, which typically involves gene mutations. Because antiangiogenic therapy resistance reflects transcriptional changes more readily generated than mutations characterizing traditional chemotherapy resistance (14), these responses may occur to some extent in all tumors treated with antiangiogenic therapy, with tumors with the greatest transcriptional changes exhibiting antiangiogenic therapy resistance.
To identify mediators of glioblastoma resistance to antiangiogenic therapy, we carried out comprehensive microarray transcriptional analysis, immunohistochemistry (IHC), and Matrigel invasion assays comparing bevacizumab-resistant glioblastomas (BRG) to their paired primary tumors, allowing us to define changes occurring in tumor cells and the microenvironment as patients' tumors progressed from bevacizumab responsive to bevacizumab resistant.
Materials and Methods
Case selection
Review of the UCSF Brain Tumor Research Center (BTRC) database identified glioblastomas meeting 2 inclusion criteria: (1) after initial responsiveness, tumor radiographic progression during bevacizumab therapy required surgery, with response and progression defined per Response Assessment in Neuro-Oncology (RANO) criteria (Supplementary Methods; ref. 15); and (2) paired pretreatment and bevacizumab-resistant tissue was available for comparison.
MRIs
Every 4 to 6 weeks during treatment, patients underwent MRIs with T1 post-gadolinium images and T2-weighted fluid attenuated inversion recovery (FLAIR) sequences (Supplementary Methods). FLAIR and T1 gadolinium-enhanced MRI scans revealing radiographic progression during bevacizumab treatment were loaded into aidScans software (AnyIntelli; Ukraine) for calculating volumes by an observer blinded to other analyses.
Immunohistochemistry
Immunostaining is described in Supplementary Methods. Vessel densities were counted from 10 representative 40× fields of vWF immunostainings by 2 observers blinded to treatment group. Stainings were quantified by converting images to binary using ImageJ software (NIH, Bethesda, MD; ref. 16).
Matrigel invasion assay
Matrigel invasion assay is described in the Supplementary Methods.
Microarray analysis
RNA was extracted from paraffin blocks with the RecoverAll Total Nucleic Acid Isolation Kit (AM1975; Ambion, Inc.), using 60 to 120 μm of paraffin-embedded slices. RNA of sufficient quality and concentration was hybridized to multiple chips (7 distinct lots) to reduce batch or chip effect. Microarray analysis used the whole-Genome DASL Assay with HumanRef-8 BeadChips (Illumina, Inc.), a platform for formalin-fixed, paraffin-embedded extracted samples. To validate primary dataset results, 3 additional paired specimens were analyzed on Whole Human Genome Oligo Microarray 4 × 44K (Agilent) and Human HT-12 v4 Expression BeadChip Kits (Illumina). Two chips whose samples clustered differently from others on cluster dendrogram heatmaps were excluded because of chip effect. After eliminating samples with insufficient gene-profile signal detection, 9 pairs underwent further analysis. Additional microarray details, microarray statistics, and gene set enrichment analysis (GSEA) are described in Supplementary Methods. Microarray data was deposited in the ArrayExpress database (accession no. E-MEXP-3296).
Real-time reverse transcriptase PCR
Real-time reverse transcriptase PCR (RT-PCR) is described in Supplementary Methods.
Statistics
Microarray statistics are described above. IHC was analyzed by paired t tests. Kaplan–Meier analysis was used to compare survival between patient groups. Wilcoxon rank sum test was used for nonparametric comparisons of paired PCR and IHC data. Student t test was used for parametric comparisons. Fisher exact test was used to compare proportions. Comparisons of more than 2 groups used ANOVA, with post hoc testing for further paired analysis. P < 0.05 was statistically significant.
Results
Patient characteristics
Twenty-one glioblastoma patients met our inclusion criteria, having a period of initial bevacizumab responsiveness followed by radiographic progression during bevacizumab treatment, along with archived tissue from surgeries before and after bevacizumab treatment (Supplementary Table S1). These 21 patients underwent BRG resection an average of 31 days after their last bevacizumab dose (range = 14–49 days). Nine patients received bevacizumab monotherapy, whereas 12 received bevacizumab combined with carboplatin, temozolomide, or CPT-11 (Supplementary Table S1).
Radiographic volumetric analysis suggests two BRG subtypes
We quantified the portion of each BRG that was nonenhancing by measuring the percent of the FLAIR-bright volume that was enhancing on the MRI revealing tumor progression during bevacizumab treatment and found a significant pattern of clustering into 2 groups, a tendency not existing in the pretreatment specimens (Fig. 1A). We called the BRGs (n = 9) with a low percentage of FLAIR-bright volume exhibiting gadolinium enhancement (average = 5%; range = 1%–10%; Fig. 1A) nonenhancing BRGs (NBRG). We called BRGs (n = 12) with a higher percentage of FLAIR-bright volume exhibiting gadolinium enhancement (average = 70%; range = 50%–97%) enhancing BRGs (EBRG; Fig. 1A and B; Supplementary Figs. S1 and S2).
To investigate the hypothesis that the nonenhancing FLAIR-bright growth seen in NBRGs reflected infiltrative growth, we analyzed cases in which tumor–brain interface was part of the specimen. NBRGs (n = 5) had more infiltrative borders than EBRGs (n = 6; Fig. 1B–D) with nearly 5-fold greater distance of furthest discontinuous tumor invasion site from the continuous tumor edge (P < 0.001; Fig. 1D). Two NBRGs also contained perivascular tumor cell invasion (Fig. 1C), a phenomenon previously reported in VEGF knockout murine gliomas (10). EBRGs were treated with bevacizumab for the same duration before progressing as NBRGs (mean 23 weeks each; P = 1.0).
Patients with EBRGs versus NBRGs were equally likely to receive bevacizumab as monotherapy or with another agent (P = 0.4). EBRG patients had equal median survival as NBRG patients after starting bevacizumab (P = 0.2; Supplementary Fig. S3A) and from time of progression during bevacizumab treatment (P = 0.1; Supplementary Fig. S3B). EBRG patients had comparable overall survival (OS) from time of diagnosis as NBRG patients (median 92 weeks NBRGs; 85 weeks EBRGs; P = 0.3).
Because bevacizumab-treated glioblastomas sometimes exhibit restricted diffusion on MRI (7), we investigated BRG diffusion weighted imaging. Nine of 21 BRGs exhibited restricted diffusion. EBRGs and NBRGs were equally likely to exhibit restricted diffusion (P = 0.4). For patients with and without restricted diffusion, OS measured from diagnosis (median 79 vs. 93 weeks; P = 0.9) and from progression on bevacizumab (median 20 vs. 28 weeks; P = 0.4) were comparable.
Characterizing the microenvironment of bevacizumab-resistant glioblastomas
We used IHC of tissue from BRGs and their paired pretreatment specimens to characterize their hypoxia and vessel density, 2 tumor microenvironment characteristics potentially affected by antiangiogenic therapy. Compared with their paired specimens before bevacizumab treatment, NBRGs exhibited reduced vessel density, from 26 vessels per high powered field (hpf) down to 10 vessels/hpf (P < 0.05) and increased hypoxia, evidenced by more than 60% increased tissue expressing hypoxia marker CA9 (ref. 17; P < 0.01) and increased HIF-1α staining (P < 0.05) compared with before bevacizumab treatment (Fig. 2; Supplementary Figs. S4–S6). Conversely, compared with their paired specimens from before bevacizumab treatment, EBRGs exhibited unchanged vessel density (P = 0.2) and unchanged hypoxia, as evidenced by CA9 (P = 0.9) and HIF-1α staining (P = 0.4; Fig. 2; Supplementary Figs. S7–S9). Like EBRGs, control paired primary and recurrent bevacizumab-naive glioblastomas exhibited unchanged vessel density and hypoxia (P = 0.5–0.9; Supplementary Fig. S10).
Molecular profiling of common glioblastoma alterations in bevacizumab-resistant glioblastomas
Profiling BRGs for common glioblastoma genetic alterations, including epidermal growth factor receptor amplification, p53 mutations, isocitrate dehydrogenase 1 (IDH1) mutations, and loss of PTEN revealed no tendency for any alteration in pretreatment specimens to be associated with a resistance pattern and no changes in any alteration after bevacizumab resistance (P > 0.05; Supplementary Table S2).
Characterizing tumor cells in bevacizumab-resistant glioblastomas
We analyzed proliferation, density, morphology, and invasiveness of tumor cells before and after bevacizumab resistance. Using paired analysis from before and after bevacizumab resistance, NBRGs exhibited unchanged tumor cell staining with proliferation marker Ki-67 (P = 0.7), whereas EBRGs exhibited increased Ki-67 staining (P < 0.05; Fig. 3A and B). EBRGs and NBRGs both exhibited nonsignificant increases in cell density after bevacizumab resistance (P = 0.2–0.4; Fig. 3C). Like BRGs, cell density did not change in recurrent versus primary bevacizumab-naive glioblastomas and, like NBRGs, Ki-67 staining did not change in recurrent versus primary bevacizumab-naive glioblastomas (P = 0.8–0.9; Supplementary Fig. S10).
Furthermore, cultured NBRG cells exhibited round enlarged cytoplasm with abundant actin-positive extensions, whereas bevacizumab-naive glioblastoma and EBRG cells exhibited more polarized stellate morphology with fewer actin-positive extensions (Fig. 3D). These differences were quantified by measuring cell dendricity, a unitless value also called inverse shape factor that describes 2-dimensional shape (ref. 18; Supplementary Methods). Dendricity was elevated in NBRG cells compared with bevacizumab-naive or EBRG cells (P < 0.05; Fig. 3E).
These morphologic differences were further investigated by immunostaining glioblastoma multiforme (GBM) tissue for cadherins, adhesion molecules whose altered expression can impact cell morphology (19). Epithelial cadherin (E-cadherin) was not expressed by 10 bevacizumab-naive glioblastomas, 8 NBRGs, or 10 EBRGs (data not shown), consistent with rarity of glioblastoma E-cadherin expression reported elsewhere (20). Neuronal cadherin (N-cadherin), whose expression can alter cell morphology, increase invasiveness and occurs during epithelial to mesenchymal transitions (EMT) in cancer (21), was expressed by 20% of bevacizumab-naive glioblastomas, 20% of EBRGs, and 100% of NBRGs (P < 0.001; Supplementary Fig. S11).
Given the potential impact of these cell morphology differences on tumor cell migration, which involves reorganizing the cytoskeleton to contract through tight spaces, or invasion, which involves extracellular matrix degradation, and our finding of greater tumor cell infiltration at the brain–tumor interface in NBRGs, we used Boyden chambers coated with control or Matrigel inserts to measure the migration or invasion, respectively, of primary cells from NBRGs, EBRGs, and cells taken by image-guided biopsies from nonenhancing or enhancing regions of bevacizumab-naive glioblastomas and found nearly 4-fold more invasiveness in NBRG cells compared with EBRG cells or cells from nonenhancing versus enhancing regions of bevacizumab-naive glioblastomas (P < 0.05; Fig. 4), with unchanged migration between the 4 groups (P = 0.6; Fig. 4).
Microarray analysis of bevacizumab-resistant glioblastomas and their paired pretreatment specimens verifies two distinct phenotypes
Microarray analysis was used to analyze BRGs and their paired pre-bevacizumab treatment tumors to identify gene expression changes associated with resistance. The 9 paired specimens (4 NBRGs and 5 EBRGs) from which RNA could be extracted were analyzed by microarrays. Unsupervised clustering of all 18 samples revealed no clustering tendency (Fig. 5A). An unbiased gene selection approach in which differential gene expression data for each of the 9 pairs was clustered on high variance probes (top 98th percentile and above, 491 probes encompassing 478 unique genes) revealed 2 clustering patterns (Fig. 5B and C). Analysis of the 2 clustering patterns revealed one to contain NBRGs and the other to contain EBRGs (Fig. 5B and C), suggesting that these clustering patterns reflected differences in BRG gadolinium enhancement. After adjusting for multiple testing, 43 genes were differentially expressed in NBRGs versus EBRGs and 146 in EBRGs versus NBRGs with adjusted P < 0.05, including NBRG expression of integrin α5, laminin, fibronectin1, and PDGFRβ, and EBRGs expressing mitogen-activated protein kinases (MAPK) 4 and 10, neural cell adhesion molecule-1 (NCAM-1), and aquaporin 4 (oncologically pertinent genes and gene subsets in Table 1; complete list of differentially expressed genes in Supplementary Tables S3 and S4). GSEA revealed that, relative to EBRGs, NBRGs exhibited significant (Q value adjusted for false discovery < 0.01) increased expression of 245 of 13,375 gene sets, including those describing wound healing, collagen production, and chemokine- and cytokine-mediated inflammation, whereas, relative to NBRGs, EBRGs exhibited significant upregulation of 61 gene sets, including those describing central nervous system (CNS) development (Table 1).
Upregulated genes by subtype . | Log2 (fold change) . | Adjusted P . | Ref. . |
---|---|---|---|
EBRG genes (n = 5; relative to 4 NBRGs) | |||
MAPK4 | 3.3 | 0.01 | |
MAPK10 | 2.1 | 0.01 | |
NCAM-1 | 1.8 | 0.03 | |
Aquaporin 4 | 3.3 | 0.02 | |
NBRG genes (n = 4; relative to 5 EBRGs) | |||
Integrin α5 | 1.4 | 0.03 | |
Fibronectin1 | 2.2 | 0.01 | |
Laminin, alpha 4 chain | 1.2 | 0.04 | |
Platelet-derived growth factor receptor β (PDGFRβ) | 1.6 | 0.03 | |
EBRG genes (n = 5; relative to paired pre-EBRGs) | |||
Angiopoietin 2 | –1.3 | 0.0001 | (45) |
Integrin β1 | –1.0 | 0.0015 | (46) |
Integrin α8 | –1.0 | 0.0010 | (47) |
NBRG genes (n = 4; relative to paired pre-NBRGs) | |||
Fibronectin1 | 1.4 | 0.0006 | (40) |
Neutrophin 3 | 1.5 | 0.001 | (43) |
Integrin α5 | 1.0 | 0.002 | (44) |
TWIST1 | 1.0 | 0.002 | (40) |
PDGFRβ | 1.1 | 0.002 | (28) |
CXCL12 | 1.1 | 0.002 | (42) |
Interleukin 6 | 1.1 | 0.002 | (48) |
GSEA | Q | ||
EBRGs (n = 5; upregulated relative to 5 NBRGs) | |||
CNS development | 1 × 10−5 | ||
NBRGs (n = 4; upregulated relative to 4 EBRGs) | |||
Regulation of wound healing | 3 × 10−3 | ||
Collagen biosynthetic process | 9 × 10−3 | ||
Collagen metabolic process | 3 × 10−3 | ||
Collagen fibril organization | 2 × 10−3 | ||
Inflammation mediated by chemokine and cytokine signaling pathways | 2 × 10−3 | ||
EBRGs (n = 5; altered relative to pre-EBRGs) | |||
Cell junctions (upregulated) | 3 × 10−4 | ||
Intermediate filaments (downregulated) | 1 × 10−4 | ||
Fibrillar collagen (downregulated) | 7 × 10−3 | ||
NBRGs (n = 4; altered relative to pre-NBRGs) | |||
Integrin signaling pathway (upregulated) | 5 × 10−9 | ||
Extracellular matrix (upregulated) | 3 × 10−11 | ||
Collagen (upregulated) | 4 × 10−5 | ||
Integrin binding (upregulated) | 9 × 10−4 | ||
Metallopeptidase activity (upregulated) | 9.5 × 10−3 | ||
Metalloendopeptidase activity (upregulated) | 8 × 10−3 | ||
Response to wounding (upregulated) | 8 × 10−10 | ||
Regulation of cell migration (upregulated) | 8 × 10−7 | ||
Regulation of cell motility (upregulated) | 8 × 10−7 | ||
Cell-matrix adhesion (upregulated) | 8 × 10−6 | ||
Chemotaxis (upregulated) | 1 × 10−4 | ||
Wnt receptor signaling pathway (upregulated) | 8 × 10−4 | ||
M phase (downregulated) | 5 × 10−4 | ||
M phase of mitotic cell cycle (downregulated) | 2 × 10−3 |
Upregulated genes by subtype . | Log2 (fold change) . | Adjusted P . | Ref. . |
---|---|---|---|
EBRG genes (n = 5; relative to 4 NBRGs) | |||
MAPK4 | 3.3 | 0.01 | |
MAPK10 | 2.1 | 0.01 | |
NCAM-1 | 1.8 | 0.03 | |
Aquaporin 4 | 3.3 | 0.02 | |
NBRG genes (n = 4; relative to 5 EBRGs) | |||
Integrin α5 | 1.4 | 0.03 | |
Fibronectin1 | 2.2 | 0.01 | |
Laminin, alpha 4 chain | 1.2 | 0.04 | |
Platelet-derived growth factor receptor β (PDGFRβ) | 1.6 | 0.03 | |
EBRG genes (n = 5; relative to paired pre-EBRGs) | |||
Angiopoietin 2 | –1.3 | 0.0001 | (45) |
Integrin β1 | –1.0 | 0.0015 | (46) |
Integrin α8 | –1.0 | 0.0010 | (47) |
NBRG genes (n = 4; relative to paired pre-NBRGs) | |||
Fibronectin1 | 1.4 | 0.0006 | (40) |
Neutrophin 3 | 1.5 | 0.001 | (43) |
Integrin α5 | 1.0 | 0.002 | (44) |
TWIST1 | 1.0 | 0.002 | (40) |
PDGFRβ | 1.1 | 0.002 | (28) |
CXCL12 | 1.1 | 0.002 | (42) |
Interleukin 6 | 1.1 | 0.002 | (48) |
GSEA | Q | ||
EBRGs (n = 5; upregulated relative to 5 NBRGs) | |||
CNS development | 1 × 10−5 | ||
NBRGs (n = 4; upregulated relative to 4 EBRGs) | |||
Regulation of wound healing | 3 × 10−3 | ||
Collagen biosynthetic process | 9 × 10−3 | ||
Collagen metabolic process | 3 × 10−3 | ||
Collagen fibril organization | 2 × 10−3 | ||
Inflammation mediated by chemokine and cytokine signaling pathways | 2 × 10−3 | ||
EBRGs (n = 5; altered relative to pre-EBRGs) | |||
Cell junctions (upregulated) | 3 × 10−4 | ||
Intermediate filaments (downregulated) | 1 × 10−4 | ||
Fibrillar collagen (downregulated) | 7 × 10−3 | ||
NBRGs (n = 4; altered relative to pre-NBRGs) | |||
Integrin signaling pathway (upregulated) | 5 × 10−9 | ||
Extracellular matrix (upregulated) | 3 × 10−11 | ||
Collagen (upregulated) | 4 × 10−5 | ||
Integrin binding (upregulated) | 9 × 10−4 | ||
Metallopeptidase activity (upregulated) | 9.5 × 10−3 | ||
Metalloendopeptidase activity (upregulated) | 8 × 10−3 | ||
Response to wounding (upregulated) | 8 × 10−10 | ||
Regulation of cell migration (upregulated) | 8 × 10−7 | ||
Regulation of cell motility (upregulated) | 8 × 10−7 | ||
Cell-matrix adhesion (upregulated) | 8 × 10−6 | ||
Chemotaxis (upregulated) | 1 × 10−4 | ||
Wnt receptor signaling pathway (upregulated) | 8 × 10−4 | ||
M phase (downregulated) | 5 × 10−4 | ||
M phase of mitotic cell cycle (downregulated) | 2 × 10−3 |
NOTE: For comparisons of expression of specific genes in EBRGs versus NBRGs, P values adjusted for multiple testing by controlling the false discovery rate are listed for representative genes with adjusted P < 0.05. For comparisons of EBRGs and NBRGs to pretreatment specimens (pre-EBRGs and pre-NBRGs), shown are representative genes (negative log FC = downregulated; positive log FC = upregulated) with raw P < 0.002, along with supportive references. For gene subsets, Q values are adjusted for false discovery and representative subsets with Q < 0.01 are shown.
Abbreviations: Ref = Reference.
We then identified transcriptional changes after bevacizumab resistance compared with before bevacizumab treatment in EBRGs and NBRGs. Compared with their pre-bevacizumab paired specimens, in EBRGs, 58 genes were upregulated and 64 downregulated with raw P values below 0.002 (Supplementary Table S5), including downregulated angiopoietin 2 and integrins β1 and α8 (Table 1). In NBRGs, 48 upregulated and 196 downregulated genes met this criterion (Supplementary Table S6), including upregulated integrin α5, fibronectin, TWIST1, and CXCL12 (Table 1). GSEA revealed that, relative to their pre-bevacizumab–treated paired specimens, NBRGs exhibited significant (Q value adjusted for false discovery < 0.01) increased expression of 348 of 13,375 gene sets, including those describing integrin signaling, the extracellular matrix, invasion, and wound healing, with decreased expression of 40 gene sets, including those describing the mitotic phase of the cell cycle (Table 1). Relative to their pre-bevacizumab–treated specimens, EBRGs upregulated 73 gene sets, including those describing cell junctions, and downregulated 50 gene sets, including those describing intermediate filaments and collagen (Table 1).
We then analyzed the expression of the 176 genes differentially expressed between EBRGs versus NBRGs for which adjusted P values were below 0.05 in a control data set of previously published microarray data from 16 glioblastomas and their paired bevacizumab-naive recurrences (22). No paired control cases exhibited differential gene expression of the 189 genes differentially expressed between EBRGs versus NBRGs in an EBRG or NBRG pattern, suggesting that these phenotypes might be BRG specific, rather than merely reflecting enhancing or nonenhancing growth patterns (Supplementary Fig. S12).
Confirming gene expression changes using real-time RT-PCR
Real-time RT-PCR was used to confirm the differential expression of genes distinguishing EBRGs from NBRGs by microarray. Integrin β1, a downregulated invasion-mediating gene (23) in EBRGs identified by microarray, was upregulated in NBRGs and downregulated in EBRGs versus paired pretreatment specimens (P < 0.05), with no change in bevacizumab-naive recurrent GBMs relative to their paired primary specimens (P > 0.05; Fig. 5D). CXCL12, an invasion-mediating gene (24) upregulated in NBRGs by microarray was upregulated in NBRGs (P < 0.05) and unchanged in EBRGs (P > 0.05) versus paired pretreatment specimens, with no change in bevacizumab-naive recurrent GBMs relative to their paired primary specimens (P > 0.05; Fig. 5D). We then used RT-PCR to measure expression of VEGF, VEGF receptor 2 (VEGFR2), and basic fibroblast growth factor (bFGF), factors that could influence antiangiogenic therapy resistance (25), but were unaltered in microarray analysis. Expression of these factors was unchanged in BRGs compared with before bevacizumab treatment and did not differ between EBRGs and NBRGs (P = 0.5–0.8; Supplementary Fig. S13). VEGF and VEGFR2 expression were further investigated at the protein level, with IHC revealing no alteration in VEGF (P = 0.7–0.8; Supplementary Fig. S14A and B) or VEGFR2 (P = 0.8–0.9, Supplementary Fig. S15) staining in EBRGs or NBRGs versus paired pretreatment specimens and Western blot revealing no alteration in VEGF expression in EBRGs versus NBRGs versus bevacizumab-naive glioblastomas (Supplementary Fig. S14C).
Validating the gene set whose expression changes reflected bevacizumab resistance pattern with additional paired specimens
To validate whether changes in expression of the 478 genes above reflected bevacizumab resistance patterns, 3 additional paired samples from glioblastomas before and after bevacizumab resistance underwent independent microarray analysis and radiographic volumetric analysis. MRI analysis revealed that all 3 tumors after bevacizumab resistance had gadolinium-enhancing volumes above 50% of their FLAIR-bright volumes (Supplementary Fig. S16), consistent with EBRG numbers. Transcriptional profiling revealed that the differential expression of our 478 gene set in all 3 tumors matched the EBRG pattern from our primary dataset (Supplementary Fig. S16).
Identifying genes whose pretreatment expression predicted resistance patterns or survival after bevacizumab resistance
To investigate whether the pretreatment transcriptional profile predicts the subsequent pattern of bevacizumab resistance, we identified the top 100 genes whose pretreatment expression associated with the development of EBRG versus NBRG resistance patterns (Supplementary Table S7). To test these pretreatment genes as predictors of resistance patterns, we compiled a secondary dataset of 15 bevacizumab-treated glioblastomas, whose pretreatment tumors were transcriptionally profiled in a previous study (26). Five of these tumors were treated for 12 months before stopping bevacizumab electively without radiographic progression, whereas 10 progressed during bevacizumab treatment, with tissue from the bevacizumab-resistant tumors not available for these cases. Volumetric analysis of secondary dataset MRIs before and after progression during bevacizumab treatment confirmed that, in some tumors, more than 50% of the FLAIR bright volume was enhancing (n = 3), consistent with EBRGs, whereas, in the remaining tumors (n = 7), less than 10% of the FLAIR-bright volume was enhancing, consistent with NBRGs. As with our primary dataset, in our secondary dataset, pretreatment MRIs lacked a nonenhancing pattern (Supplementary Fig. S17A) and EBRGs and NBRGs exhibited similar bevacizumab treatment duration (24 weeks NBRGs vs. 19 weeks EBRGs; P = 0.5).
The top 100 genes in pretreatment tumors associating with EBRG versus NBRG development in our primary dataset (Supplementary Table S7) failed to predict resistance patterns in our secondary dataset (Supplementary Fig. S17B). Furthermore, the top 100 pretreatment genes associating with resistance pattern in the secondary dataset (Supplementary Table S8) failed to predict resistance pattern in our primary dataset (Supplementary Fig. S17C). There was one common gene between the top 100 genes predicting resistance pattern in our primary and secondary datasets, VHLL, a Von Hippel–Lindau tumor suppressor homolog, with expression elevated in pretreatment EBRGs versus pretreatment NBRGs. VHLL expression increased in each tumor developing nonenhancing resistance and decreased in each tumor developing enhancing bevacizumab resistance (Supplementary Fig. S17D). Thus, although there was a gene set whose differential expression after bevacizumab resistance versus before was associated with resistance pattern, there was no pretreatment gene set predicting which resistance pattern would develop or whether a complete treatment course could be completed without resistance.
We also investigated whether pretreatment gene expression correlated with survival after starting bevacizumab treatment by identifying the top 1,000 genes associated with OS in our primary (Supplementary Table S9) and secondary datasets (Supplementary Table S10), of which only 14 genes were in common (Supplementary Table S11), with one, TGF, beta-induced (TGFBI), promoting decreased OS with increased pretreatment expression and possessing oncologic functions as an extracellular matrix protein (27). Similar analysis identifying the top 1,000 genes associated with progression-free survival (PFS) in our primary and secondary datasets revealed 10 common genes (Supplementary Table S12), with one oncologically pertinent gene, neurofibromin (merlin; ref. 28), promoting decreased PFS with increased pretreatment expression.
Gene expression–based classification of BRGs and their paired pretreatment specimens into molecular subtypes
Classification of tumors with sufficient genetic material into the 3 glioblastoma subtypes described by Phillips and colleagues (22) revealed that, after bevacizumab resistance, 3 EBRGs transitioned from mesenchymal to proneural and 3 NBRGs transitioned from proneural or proliferative to mesenchymal. BRG classification per the 4 subtypes defined using The Cancer Genome Atlas (TCGA; ref. 28) did not change after bevacizumab resistance (Supplementary Table S2).
Discussion
Despite encouraging results in clinical trials (1, 2), the benefit of bevacizumab in glioblastoma is typically transient, followed by tumor growth resumption, a phenomenon found in other cancers treated with antiangiogenic therapy (5) and associated with a poor prognosis (13). Antiangiogenic therapy resistance has been associated with increased tumor cell invasion (10, 11), although other mechanisms, including recruitment of bone marrow–derived cells to tumors (29), have also been suggested.
In our radiographic analysis of 21 bevacizumab-resistant glioblastomas, over half, which we called enhancing bevacizumab-resistant glioblastomas (EBRG), exhibited enhancing growth on MRI that has defined glioblastoma recurrence since the 1990 advent of the Macdonald response criteria (30). The remaining cases, termed nonenhancing bevacizumab-resistant glioblastomas (NBRG), were characterized by nonenhancing FLAIR-bright growth on MRI which has been associated with glioblastoma recurrence since revised treatment response criteria were put forth by the Response Assessment in Neuro-Oncology (RANO) group in 2010 (15) to recognize changes in the presentation of recurrent glioblastoma since the advent of antiangiogenic therapy. Although the initial identification of nonenhancing recurrence after bevacizumab treatment led some to suggest that the phenomenon was common, subsequent studies (6, 8) and review of 74 glioblastomas progressing during antiangiogenic therapy at our institution (unpublished observations) showed that 60% to 75% of glioblastomas progressing during antiangiogenic therapy are enhancing local recurrences.
Other investigators have correlated glioblastoma PFS and OS after bevacizumab with the pretreatment ratio of the FLAIR bright volume to the T1 gadolinium–enhanced volume (31). Our analysis expands upon these observations by finding 2 radiographic subtypes of glioblastomas progressing during bevacizumab treatment and showing distinct transcriptional phenotypes for each subtype, with the functional correlate of increased cellular invasiveness in one subtype. Invasion, proliferation, and hypoxia have not correlated well in prior glioblastoma studies (32, 33) and did not correlate well in our bevacizumab-naive glioblastomas, so it is possible that the trend we observed of NBRGs becoming more hypoxic, more invasive, and equally proliferative relative to pretreatment may occur specifically after resistance to antiangiogenic therapy.
An unbiased gene selection approach in which differential gene expression measuring changes in expression after bevacizumab resistance compared with before was clustered on high variance probes revealed 2 clustering patterns, analysis of which confirmed one to be NBRGs and the other EBRGs. IHC showed that NBRGs, which exhibit FLAIR-bright nonenhancing tumor progression on MRI that has been associated with devascularized infiltration (6), maintained the hypoxia and reduced vascularity reported after successful bevacizumab treatment (34) and showed unchanged proliferation marker expression compared with before bevacizumab treatment. The unchanged tumor cell proliferation in NBRGs could reflect reliance on invasion, particularly perivascular invasion, to reduce angiogenesis dependence by allowing cellular migration away from devascularized areas into areas closer to blood vessels, rather than continuing the nodular enhancing exponential growth that is difficult to sustain during VEGF blockade. Conversely, EBRGs, which exhibit MRI enhancement potentially consistent with neovascularization, reacquired the reduced hypoxia and increased vascularity seen before bevacizumab treatment, which may have contributed to the observed increased cell proliferation compared with before bevacizumab treatment.
These differences in tumor cells and the tumor microenvironment between EBRGs and NBRGs were reflected in the genes found to be differentially transcribed between EBRGs versus NBRGs. Specifically, the transcriptional data showed upregulation in NBRGs of integrin α5β1 and 2 of its ligands, fibronectin, and laminin. Increased tumor cell expression of α5β1 could bind upregulated fibronectin or laminin in the vascular basement membrane, possibly promoting perivascular invasion as a mechanism of antiangiogenic therapy resistance. Conversely, EBRGs exhibited increased expression of aquaporin 4, which promotes vascular proliferation (35), possibly allowing EBRGs to reacquire their pretreatment vessel densities and obtain the blood supply needed for the increased cellular proliferation we found in EBRGs, as evidenced by increased Ki-67 staining and increased MAPK4 and MAPK10 expression. We confirmed a functional impact of these molecular differences identified in the microarray analysis by showing NBRG cells to be more invasive in Matrigel-coated Boyden chambers than EBRG cells. The hypothesis that the transcriptional reprogramming we identified is specific to bevacizumab resistance, rather than merely reflecting differences between nonenhancing versus enhancing glioblastoma, is supported by 3 lines of evidence. First, array data on glioblastomas and their paired bevacizumab-naive recurrences did not cluster toward either an EBRG or NBRG pattern using the gene set we compiled by microarray analysis. Second, confirmatory RT-PCR showed no alteration in bevacizumab-naive recurrent glioblastomas relative to their paired initial tumors of transcripts altered in NBRGs or EBRGs. Third, genes we identified to be differentially expressed between EBRGs versus NBRGs did not overlap with previously reported genes differentially expressed between enhancing versus nonenhancing glioblastoma, such as the 79 genes differentially expressed between incompletely versus completely enhancing glioblastomas (36) or the 643 genes differentially expressed between the enhancing periphery and the central necrotic core of glioblastoma (37).
It has been suggested that EMT, a phenomenon recognized in non-CNS malignancies and associated with a worsened prognosis, metastases, and chemotherapy resistance (21) can occur after antiangiogenic therapy (38). Some have proposed EMT in glioblastoma (39), but it remains unclear what the equivalent of a nonmigratory epithelial state or a mesenchymal state with metastatic potential is in glioblastoma. Although the glioblastoma subtypes of Phillips and colleagues (22) and the TCGA subtypes (28) included a mesenchymal subtype, it is unclear whether these mesenchymal subtypes embody EMT features, and the unchanging nature of TCGA subtypes implies that they do not detect EMT. Regardless, our findings that NBRGs had (i) examples that converted to the Phillips and colleagues (22) mesenchymal subtype after bevacizumab resistance; (ii) increased frequency of N-cadherin expression, a feature of EMT in non-CNS tumors (21); (iii) increased expression of TWIST1, an EMT-regulating transcription factor (40); and (iv) mesenchymal morphology with structures resembling pseudopodia (41) suggest that further work to determine what constitutes EMT in glioblastoma and under what conditions antiangiogenic therapy promotes EMT is warranted.
Although immunostaining primary tumor cells from BRGs revealed differences in N-cadherin expression between NBRGs and EBRGs, N-cadherin gene expression was unaltered by adjusted P value in microarray analysis. This discrepancy between microarray data and immunostaining suggests possible posttranscriptional or translational regulation in addition to the transcriptional differences we identified, a possible subject for future work.
There are limitations to any study using infrequently available clinical specimens. The first limitation, present in over half the cases analyzed, is that bevacizumab is usually combined with other treatments. However, there was no differential tendency for combination treatment in EBRGs versus NBRGs. Furthermore, the large number of genes changes is common between tumors, regardless of other treatments received, along with evidence supporting these genes as potential mediators of antiangiogenic therapy resistance (28, 40, 42–48), suggests that additional treatments did not confound our analysis. A second limitation is small sample size, reflecting the few cases undergoing surgery after progression during bevacizumab treatment for which pretreatment tissue was available. Although it is possible that an expanded sample size could uncover other resistance patterns with mediators distinct from those listed here, our sample size was large enough to uncover statistically significant clustering of tumors resistant to antiangiogenic therapy by differentially expressed genes and provides a platform to guide further efforts to define and disrupt resistance to antiangiogenic therapy. A third limitation is that our study selected BRGs amenable to surgical resection. Although one might hypothesize that selecting surgically resectable cases would bias our study toward enhancing local recurrences, we found 57% of resected recurrences to be EBRGs, less than the 74% of overall progressions on bevacizumab treatment that were enhancing local during this time period (unpublished observations). Regardless, by studying cases of bevacizumab resistance leading to surgery, a necessity to study this tissue, our findings may not be applicable to unresectable bevacizumab-resistant tumors. A fourth limitation is that, whereas a comparison of NBRGs to EBRGs revealed several genes with adjusted P values below 0.05, no differentially expressed oncologically pertinent transcripts in NBRGs or EBRGs relative to their paired pre-bevacizumab–treated tumors had significant adjusted P values, despite having significant raw P values, and some exhibited low upregulation (Table 1). However, these oncologically pertinent transcripts had prior evidence (28, 40, 42–48) supporting their potential involvement in antiangiogenic therapy resistance, reducing concerns about needing to adjust for multiple testing (49) and suggesting that their raw P values being statistically significant is sufficient to render them appropriate for further investigation, whereas the fold changes we detected by real-time RT-PCR were higher than those detected by microarray and have been suggested by some to better reflect changes in individual transcript levels (49).
Thus, just as prolonged treatment with temozolomide, the current standard of care for newly diagnosed glioblastoma, can create a “hypermutator” phenotype associated with recurrence (50), VEGF-targeted treatments such as bevacizumab may cause “hyperinvasive” (NBRG) or “hyperproliferative” (EBRG) phenotypes associated with glioblastoma recurrence. These findings provide important biologic insight into how tumors counteract antiangiogenic treatments, responses that have unfortunately limited the efficacy of antiangiogenic therapies used in patients to date.
Disclosure of Potential Conflicts of Interest
W.S. Carbonell has been a founder and CEO of OncoSynergy, Inc. M.K. Aghi is a Scientific Advisory Board Member and shareholder of OncoSynergy, Inc. The other authors disclosed no potential conflicts of interest.
Acknowledgments
The authors thank Joanna S. Phillips who provided neuropathology assistance, and Claudia Petritsch and Susan Chang who reviewed the manuscript.
Grant Support
This work was supported in part by funding to the laboratory of MKA from the American Brain Tumor Association, the James S. McDonnell Foundation, the American Cancer Society, the NIH (5K02NS64167-2), and the UCSF Brain Tumor SPORE CA097257. The secondary pretreatment dataset arose from work supported by Genentech (South San Francisco, CA) and the State of California through the University of California Industry–University Cooperative Research Discovery grant (Bio05-10501). A. Jahangiri is an HHMI Research Fellow.
The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.