Assessing the significance of chromosomal aberrations in cancer: Methodology and application to glioma

GISTIC Methodology.

GISTIC identifies significant aberrations through two key steps [Fig. 1 and supporting information (SI) Text]. First, the method calculates a statistic (G score) that involves both the frequency of occurrence and the amplitude of the aberration. Second, it assesses the statistical significance of each aberration by comparing the observed statistic to the results that would be expected by chance, using a permutation test that is based on the overall pattern of aberrations seen across the genome. The method accounts for multiple-hypothesis testing using the false-discovery rate (FDR) framework (9) and assigns a q value to each result, reflecting the probability that the event is due to chance fluctuation. For each significant region, the method defines a “peak region” with the greatest frequency and amplitude of aberration. Each peak is tested to determine whether the signal is due primarily to broad events, focal events, or overlapping events of both types.

Application to Glioma.

We applied the method to a collection of 141 gliomas, including 107 primary GBMs, 15 secondary GBMs, and 19 lower-grade gliomas (SI Table 2). We hybridized genomic DNA to microarrays containing probes for ≈100,000 SNPs to identify copy-number changes and loss of heterozygosity (LOH). A genome-wide view of the copy-number alterations is shown in Fig. 2a (LOH results are described in SI Note 1 in SI Text). The overall pattern is complex, with almost every region of the genome being altered in at least one tumor. Nonetheless, only 16 broad and 16 focal events are significant. Focal events are superimposed on four broad events (including two focal events on chromosome 7 and single events on chromosomes 9, 10, and 13), resulting in a total of 28 peak regions of amplification, deletion, and LOH.

The 16 broad events include six amplifications (chromosomes 7, 8q, 12p, 17q, 19p, and 20), nine deletions (6q, 9p, 10, 11p, 13, 14, 16q, 19q, and 22), and one region of copy-neutral LOH (17p) (Fig. 2b, SI Table 3, and SI Note 1 in SI Text). These events occur at high frequency (range 10–70%, median 27%). In particular, amplification of chr7 and deletion of chr10 each affect >60% of our samples, including >80% of our primary GBMs. For broad regions without superimposed focal events, the peak regions are large (median of 110 genes) (SI Table 3).

The 16 focal events tend to occur at lower frequencies than the broad aberrations (range 6–49%, median 14%) (SI Table 3). Among these, amplifications of 4q12 and 7p11.2 (18–26% of samples) and deletions of 1p36.31 and 9p21.3 (35–49%) are the most frequent. In some cases, a high degree of amplification renders amplifications highly significant even though they occur in only 6–7% of samples (for example, the regions containing CDK4 and MDM2 on chr12). Because the background rate of deletions across the genome is higher, deletions usually must occur at higher frequencies than amplifications to attain similar levels of significance (SI Table 3). The peak regions for the focal events can be localized to small regions (median of four genes).

Analysis Confirms Known Genes and Identifies New Loci.

We compared the 28 peak regions to the locations of oncogenes and tumor-suppressor genes previously implicated in the pathogenesis of glioma. A recent review (10) lists 12 such genes reported to be altered in multiple studies of glioma (TP53, RB1, CDKN2A/B, PTEN, EGFR, PDGFRA, MET, CDK4, CDK6, MDM2, MDM4, and MYC). We found that 11 of these 12 genes each correspond with one of the 28 peak regions (with one of these, MYC, lying just beyond the boundaries defined by GISTIC; as described below, this slight discrepancy is resolved with additional data) (SI Table 3). The 12th gene (CDK6) lies within the broad region of significant amplification on chr7, although it does not correspond to a peak. Interestingly, TP53 is within the single peak region of LOH that is not reflected in a peak of copy-number change, suggesting that it is primarily inactivated through copy-neutral LOH (SI Note 1 in SI Text).

An additional five peak regions contain genes that are known to play a role in other cancers (MYCN, PIK3CA, CCND2, KRAS, and CHD5) (11, 12). Our analysis suggests that chromosomal aberrations involving these genes are also relevant for glioma pathogenesis. These genes should therefore be carefully characterized in glioma.

The remaining 12 regions (43% of the total) are not associated with known cancer-related genes. These events occur at substantial frequencies (6–37%), but nine are due to broad events and two others rarely reach high amplitude. The final region (12q14.3) undergoes high-level amplifications, but always in concert with amplifications of either of two neighboring regions containing CDK4 and MDM2 (data not shown), suggesting that it may be due to structural features required to amplify these genes [as is the case in dedifferentiated liposarcoma (13, 14)]. The fact that nearly half the regions (including regions affected by highly prevalent aberrations) are not yet associated with known cancer-related genes underscores the importance of systematic analysis of the cancer genome.

Previous studies of copy-number alterations in glioma have shown distinct patterns for certain subtypes such as primary vs. secondary GBMs (15) or astrocytic vs. oligodendroglial tumors (16). To explore whether our combined analysis of these glioma subtypes prevented the detection of alterations specific to primary GBMs, we performed a separate analysis on only the 107 primary GBMs in our sample set. No additional statistically significant alterations were identified (SI Fig. 4).

Consistency Across Independent Data Sets.

We then sought to compare our results with two previous studies of copy-number alterations in glioma (178 samples on 100K SNP arrays; 37 samples on a 16K CGH arrays) (15, 17). At first glance, there appear to be striking differences. The previous studies reported many more regions of interest (208 and 97) (Table 1), but the regions included fewer of the known glioma-associated genes and show low concordance with one another. The differences are attributable to the methodology used in these studies, in which minimal common regions of copy-number change are reported, without explicitly taking into account the degree of background noise (see SI Note 2 in SI Text). Applying a similar analysis to our own data identifies a similarly large number of regions (144) (see SI Note 2 in SI Text), but these include fewer known glioma-associated genes and show low concordance with the other studies. By contrast, applying GISTIC to the raw data from the other two studies identifies 24 and 26 significant regions each. Importantly, these regions agree closely with the 27 regions (excluding copy-neutral LOH of 17p, because it is not observed in the copy-number analyses) identified above (SI Fig. 5), and they include essentially the same glioma-associated genes. Moreover, these results are specific to glioma, as seen in a comparison to lung cancer (SI Fig. 6) (2). The strong concordance across three independent data sets and two different platforms supports the validity of both the GISTIC methodology and our results for glioma.

Given the close agreement across these data sets, we combined our initial data set with the data from the prior study performed on the same platform to obtain a pooled data set with 319 glioma samples. This analysis identified 34 significant regions, including 27 of the 28 identified in the initial data set (SI Table 3). For the additional regions, the prevalence of aberrations is similar in the initial and combined data sets (11.0% vs. 11.3%), but they now exceed the significance threshold due to the larger sample size. Increasing the sample size also leads to narrower peak regions, with the median number of genes decreasing from 12 to five per region. Inclusion of additional data sets may define these regions with even greater precision and facilitate identification of the gene targets.

Overlapping Broad and Focal Events May Have Different Consequences.

Having identified various instances of overlapping broad and focal events, we explored whether they have distinct functional consequences, using chr7 as an example (SI Fig. 7). We first compared copy-number profiles to gene expression among a group of 43 primary GBMs for which we had sufficient material for a combined analysis. EGFR was overexpressed in most GBMs with focal EGFR amplification (7^gainEGFR^amp) but in none of the GBMs with broad amplification of chromosome 7 in the absence of focal EGFR amplification (7^gain) (Fig. 3a). We found three additional sources of evidence supporting the biologic distinction between these two classes of tumors. First, recognizing that EGFR-amplified tumors are known to carry a low rate of mutations in the TP53 gene (18), we sequenced TP53 and found these mutations more frequently associated with 7^gain compared with 7^gainEGFR^amp (two-sided Fisher's exact test, P = 0.03). Second, we found that 7^gain is less frequently associated with EGFR point mutations or expression of the EGFRvIII deletion mutant [determined for many of our tumors in a prior study (19)] (P = 0.001). Third, we found that 7^gain but not 7^gainEGFR^amp occurs frequently in secondary GBMs (P = 0.16). These three findings are consistent with earlier observations regarding 7^gainEGFR^amp (20, 21) and further suggest that 7^gain has distinct functional consequences.

To explore the function of 7^gain, we identified genes on the chromosome that show extreme outliers in expression in at least 10% of tumors with 7^gain, compared with 7^normal (SI Note 3 in SI Text). The notion behind this “comparative outlier analysis” is that broad events such as 7^gain may have heterogeneous effects across various tumors. We are interested in identifying genes that are strongly up-regulated in even a subset of the samples.

Strikingly, two of the top four genes in this analysis are those encoding the receptor MET and its ligand HGF (SI Table 4). Approximately one-third of 7^gain events are associated with either MET or HGF overexpression (Fig. 3a), and tumors that overexpress one tend to overexpress both (P = 0.06) (data not shown). The overexpression of MET and HGF appears to be functionally relevant: we studied glioma cell lines with 7^gain and increased expression of MET and HGF (Fig. 3b) and found phosphorylation and activation of the MET receptor even under serum-starved conditions (Fig. 3c). These cell lines showed enhanced responsiveness to the MET kinase inhibitor SU11274 (22) (Fig. 3d) at drug concentrations that inhibit the MET signaling pathway (SI Fig. 8 A and B). None of the GBM cell lines with 7^gain showed constitutive activation of EGFR or were responsive to the EGFR kinase inhibitor erlotinib (SI Fig. 8 C and D). Compared with the relatively rare (5–7%) focal amplification of the MET gene locus (SI Table 3), 7^gain with overexpression of MET and HGF may provide a more common mechanism for cell-autonomous activation of the MET signaling pathway in glioma. This finding may be relevant for the clinical deployment of inhibitors targeting this network in glioma and other cancers (23, 24).

Clinical Samples and Cell Lines.

Genomic DNA was extracted from fresh-frozen tumors samples using DNeasy (Qiagen). Non-tumor tissue, including paired normal brain corresponding to 10 gliomas, was used for germ-line control DNA. Collection and analysis of all clinical samples were approved by the University of California Institutional Review Board. RNA and DNA were also obtained from the GBM cell lines 8-MG-BA, A172, DK-MG, GAMG, HS683, LN-18, SF-268, SF-295, SNB-75, T98G, and U251. Gastric and lung cancer cell lines MKN-45 (high level MET amplified) and H3255 (L858R EGFR mutant) were included as positive controls in experiments with the MET kinase inhibitor SU11274 and the EGFR kinase inhibitor erlotinib, respectively.

SNP Arrays.

Genomic DNA was applied according to the manufacturer's instructions to oligonucleotide arrays (Affymetrix) interrogating 116,204 SNP loci on all chromosomes except Y. Arrays were scanned by using the GeneChip Scanner 3000, and genotyping was performed by using Affymetrix Genotyping Tools Version 2.0. Probe-level signal intensities were normalized to a baseline array with median intensity using invariant set normalization (27). SNP-level signal intensities were obtained by using a model-based (PM/MM) method (28). Further analytic steps are described in SI Text. SNP, gene, and cytogenetic band locations are based on the hg16 (July 2003) genome build (http://genome.ucsc.edu). Data from SF-268, SF-295, SNB-75, and six gliomas (along with paired normals) were previously published (29, 30). Data are available at www.broad.mit.edu/cancer/pub/GISTIC and through the Gene Expression Omnibus (www.ncbi.nlm.nih.gov/projects/geo, accession no. GSE9635).

Mutation Detection.

PTEN and TP53 were sequenced in 134 of the 141 samples undergoing SNP analysis as previously described (19). All exons were covered in >70% of samples except exons 1, 8, and 9 in PTEN and exons 7 and 10 in TP53. EGFR point mutations and vIII expression were determined for 133 and 58 samples, respectively, in a prior study (19).

Expression Arrays.

Expression data were obtained by using Affymetrix U133A/B and plus 2 arrays from 43 primary GBMs undergoing SNP array analysis. CEL files from U133A and plus 2 arrays were preprocessed separately by using RMA (31). Probe sets common to both arrays were used after equalizing the mean and standard deviation of the U133A and plus 2 arrays. Expression data from GBM cell lines were generated by Affymetrix cartridge arrays, except SF-268, SF-295, and SNB-75, where available U133A data were used (http://wombat.gnf.org/index.html).

EGFR and MET Inhibition.