Comprehensive identification and characterization of somatic copy number alterations in triple‑negative breast cancer.

Triple‑negative breast cancer (TNBC) accounts for ~15% of all breast cancer diagnoses each year. Patients with TNBC tend to have a higher risk for early relapse and a worse prognosis. TNBC is characterized by extensive somatic copy number alterations (CNAs). However, the DNA CNA profile of TNBC remains to be extensively investigated. The present study assessed the genomic profile of CNAs in 201 TNBC samples, aiming to identify recurrent CNAs that may drive the pathogenesis of TNBC. In total, 123 regions of significant amplification and deletion were detected using the Genomic Identification of Significant Targets in Cancer algorithm, and potential driver genes for TNBC were identified. A total of 31 samples exhibited signs of chromothripsis and revealed chromosome pulverization hotspot regions. The present study further determined 199 genomic locations that were significantly enriched for breakpoints, which indicated TNBC‑specific genomic instability regions. Unsupervised hierarchical clustering of tumors resulted in three main subgroups that exhibited distinct CNA profiles, which may reveal the heterogeneity of molecular mechanisms in TNBC subgroups. These results will extend the molecular understanding of TNBC and will facilitate the discovery of therapeutic and diagnostic target candidates.

Introduction

Triple-negative breast cancer (TNBC) is defined by a lack of expression of the estrogen receptor (ER), a lack of progesterone receptor (PR) expression and the absence of ERBB2 gene amplification, which encodes human epidermal growth factor receptor 2 (HER2) (1,2). TNBC accounts for ~15% of all breast cancer cases worldwide and represents a heterogeneous group of breast tumors (3-5). Numerous patients with TNBC experience a fast relapse and commonly develop metastases, which results in a poor prognosis (2). Recent advances in genomic profiling technologies have provided significant insights into the pathogenesis of breast cancer, including TNBC (6-9). However, presently, no individualized targeted adjuvants or induction treatments for TNBC are available. Given the lack of recurrent targetable genomic alterations, functional characterization of the TNBC genome to identify genomic driver events is of utmost importance.

Somatic copy number alterations (CNAs) are a universal feature of human cancer (10-12). Compared with any other type of somatic genetic alteration, CNAs alter a greater portion of the cancer genome. In general, CNAs are associated with patient prognoses and therapeutic resistance (13,14). Different cancer types adopt copy number changes in different ways to shape their genomes (15,16). These CNAs can affect the expression of genes and/or influence the regulation of genes in their vicinity (17,18). Furthermore, CNAs serve an important role in the classification of tumor subtypes (9,19,20). However, comprehensive genomic profiling of TNBC remains to be adequately assessed. Given the prevalence of CNAs in cancer, significant progress has been made in understanding the functional impacts of CNAs, as well as the potential driver genes they contain (21,22). These studies provide a rich source of data that allows for performing meta-analysis.

Previously, chromothripsis has been described as a novel mechanism for cancer initiation and progression (23-25). In the classic model, tumorigenesis is an evolutionary process in which the transition of normal cells to neoplastic cells is mediated by the accumulation of somatic mutations. However, the phenomenon of chromothripsis reveals a new paradigm of oncogenic transformation, involving a catastrophic mutational process that is commonly observed in numerous cancer types, such as colorectal cancer, neuroblastoma and acute myeloid leukemia (26-28). Chromothripsis contrasts with the multistep model of cancer development and is characterized by the shattering of one or multiple chromosomes in a single catastrophic event. Subsequently, shattered fragments are randomly stitched together by DNA double-strand break repair to form a derivative chromosome (29-31). This process can lead to the simultaneous acquisition of multiple tumor-promoting lesions. For example, it may result in a large number of structural genome variations, including duplications, deletions, inversions and translocations. In addition, it may give rise to the amplification of oncogenes or the inactivation of tumor suppressors, which serve an important role in oncogenesis. To the best of our knowledge, a comprehensive evaluation of chromothripsis in TNBC has not yet been performed.

The present study collected 201 TNBC samples from the National Center for Biotechnology Information (NCBI) Gene Expression Omnibus (GEO) (32) and The Cancer Genome Atlas (TCGA) (33) databases to perform a meta-analysis of genomic CNAs. Statistically significant recurrent CNAs were identified and the distribution of CNA breakpoints along the cancer genome was examined. By employing the Genomic Identification of Significant Targets In Cancer (GISTIC) algorithm, a total of 123 regions of DNA amplification and deletion were obtained (34). In these regions, a number of potential driver genes for TNBC were revealed. For deletion of the chro-mothripsis phenomenon, the typical chromothripsis pattern was identified based on CNA data, and 31 samples with signs of chromothripsis were detected (35-37). Further analysis for chromothripsis regions was performed, which revealed the chromosome pulverization hotspots in TNBC. Furthermore, unsupervised hierarchical clustering of TNBC samples demonstrated three different clusters that corresponded to the tumors with specific CNA profiles (38). The present results extended characterization of the CNA landscape of TNBC genomes and provided novel insight into the phenomenon of chromothripsis. Thus, the present findings provide information that improves understanding of the mechanisms of TNBC development and may improve targeted therapies.

Materials and methods

Sample collection and CNA calling

All samples included in the present study were obtained from the NCBI GEO (32) and TCGA databases (33). The selection criteria were as follows: i) The patient was diagnosed with primary breast cancer, and patients with distant metastasis were excluded; ii) all samples were hybridized using Affymetrix Genome-Wide Human SNP Array 6.0 platform to facilitate data integration; and iii) patients were histologically confirmed as ER-, PR- and HER2-negative, or the study clearly stated that patients were diagnosed with TNBC.

To identify somatic CNAs, raw signal intensity files were downloaded from the GEO website (https://www.ncbi.nlm.nih.gov/geo/) for re-analysis. The R package aroma. affymetrix (https://aroma-project.org/; version 3.2.0) using the CRMAv.2 method (39) was employed for data processing. To call CNA segments, HapMap (40) data were used as a control. The data annotation was based on human reference genome assembly hg19/GRCh37 (41), and the circular binary segmentation algorithm (42) was performed to segment copy number data. Next, the CNA calling cut-off values for amplifications and deletions were set to 0.2 and -0.2, respectively. To avoid gender bias, the X and Y chromosomes were excluded. TCGA data portal (https://portal.gdc.cancer.gov/) was used for downloading the genomic array data (level 3) and clinical information. Visual inspection was used for data quality control, and samples of poor quality were excluded from further analyses. A total of 201 TNBC cases were collected and are presented in Table SI. The samples were collected between October 2011 and August 2018.

Identification of significant recurrent targets

The GISTIC algorithm (34) was used for the identification of peak regions that were significantly amplified or deleted in all samples. The parameters used to run GISTIC 2.0 were as follows: i) The false discovery rate q-value was set to <0.05; ii) the arm peel method (34) was used to reduce data noise; iii) the confidence level used to calculate the region containing a driver was set to >0.95; iv) the 'Extreme’ method (34) was applied for reducing marker-level to the gene-level copy number data; and v) the log2 ratios for calling gains and losses were set to >0.2 and <-0.2, respectively. In total, 719 known cancer consensus genes were downloaded from the Catalogue of Somatic Mutations in Cancer (COSMIC) database (43).

Chromothripsis screening

Chromothripsis-like patterns were detected using the CTLPScanner web server (http://cgma.scu.edu.cn/CTLPScanner/) (35). The pattern of oscillating copy number changes and localized clustering of breakpoints was screened based on segmented array data. The parameters and thresholds for the screening were as follows: i) Copy number status switch times ≥20; ii) log10 of likelihood ratio ≥10; iii) minimum segment size of 10 Kb; and iv) signal value difference between two adjacent segments ≥0.4. For visualization of results, signal intensity for calling genomic gains and losses were set to 0.2 and -0.2, respectively.

Detection of chromosomal breakpoints

In the present study, both boundaries of each CNA were defined as chromosomal breakpoints. A stringent definition of CNA breakpoints was used to reduce the bias caused by technical noise. If the alteration of log2 signal intensity between two adjacent genomic segments was >0.4, the related genomic position was considered to be a breakpoint (44,45). CNAs <10 Kb were ignored. Breakpoints located in chromosomal telomeres and centromeres were excluded for further analysis. To generate a simulated distribution of breakpoints, in-house Perl scripts were used to bin the genome and randomly shuffle the positions of breakpoints 10,000 times. Common fragile sites (CFSs) and non-fragile regions of the human genome were extracted from previous studies (46-48). All the genome coordinates were converted from genome assembly hg18/NCBI36 to hg19/GRCh37 by the LiftOver tool (41). Furthermore, the copy number based unsupervised hierarchical clustering was performed based on the Euclidean distance following Ward's method.

Results

CNA profile in TNBC

In total, 98 and 103 TNBC samples were collected from the GEO and TCGA databases, respectively. The genome-wide frequency plot of CNAs is presented in Fig. 1. The CNAs were of varying sizes, ranging from whole chromosome arms to small focal amplifications and deletions. The mean number of CNAs per tumor was 414, with deletions outnumbering amplifications at ~1.7:1. From the high-resolution CNA data, it was observed that the most common gains of the entire chromosome arms included 1q, 8q, 10p and 12p, and losses of 5q, 8p and 17p. Other less frequently affected chromosomes included 9p and 20q gains and 13q, 14q and 15q losses. The most frequent focal gains were narrowed down to 3q and 19q. Focal losses were identified most often in 3q and 12q. Notably, this profile may reveal TNBC-specific CNAs that have rarely been reported, such as 9p gain and 5q, 12q loss.

Figure 1

Frequency plot of copy number alterations for each chromosome. The histogram presents the percentage of samples with specific alterations. Copy number gains and losses are depicted in red and blue, respectively.

The present study identified statistically significant recurrent focal CNAs and potential cancer driver genes by GISTIC 2.0 (34). Using this algorithm, 49 and 74 amplification and deletion peaks were identified, respectively (q<0.05). The annotation of these peaks revealed 993 targeted coding genes. The significantly amplified or deleted genomic regions as well as identified genes are presented in Table SII. There were 33 regions of interest that contained only one candidate driver gene. The aberration score from GISTIC is presented in Fig. 2, in addition to some of the peaks containing only one candidate gene.

Figure 2

Amplifications and deletion peaks identified by GISTIC algorithm. The x-axis represents the aberration score of GISTIC. Some of the peaks contain only one candidate driver gene, which are shown across the genome. GISTIC, Genomic Identification of Significant Targets In Cancer; Chr, chromosome.

Identification of chromothripsis events

Using the CTLPScanner algorithm (35), a total of 31 chromothripsis cases were identified from 201 TNBC samples, with an incidence of ~15% (Table SIII). Since the overall chromothripsis incidence in various cancer types is ~5%, TNBC has a relatively high chromothripsis incidence compared with most other tumor types, including other subtypes of breast cancer (23,36). Pulverization regions were identified in various sizes ranging from dozens of Mb to the entire chromosome arm. Fig. 3 presents examples of identified chromothripsis cases and pulverized genomic regions. The present study further investigated the number of affected chromosomes per tumor sample, and observed that ~22% (7/31) of chromothripsis cases carried two shattered chromosomes. In other cases, only one chromosome was affected. The most frequently pulverized focal genomic regions in TNBC included 6p, 11p and 17q. At the chromosome level, chromosomes 1, 5 and 12 demonstrated relatively high rates of pulverization. Fig. 4 illustrates the hotspots of chromothripsis across the cancer genome and may reveal a TNBC-specific pattern of chromothripsis. Some of the GISTIC-identified candidate driver genes that were located in the chromothripsis hotspot regions are also presented. Alterations in these genes may further contribute to chromosomal instability and the overall chromothripsis phenotype in TNBC.

Figure 3

Examples of detected chromothripsis events. (A) A case of pulverization of the short arm of chromosome 3. (B) A case of pulverization of the short arm of chromosome 12. The x-axis represents genomic locations, and the y-axis denotes the log2 ratio of the copy number. Red and blue lines represent genomic gains and losses, respectively. chr, chromosome.

Figure 4

Circos plot of chromothripsis hotspots and somatic copy number alteration breakpoints. The outermost circle presents the chromosome number and bands. The next circle represents potential cancer driver genes located in chromosome pulverization hotspot regions. The third circle represents hotspots of chromothripsis region. The depth of red represents the frequency of chromothripsis occurrence. The fourth circle shows the identified breakpoint-prone genomic regions in triple-negative breast cancer. The innermost circle represents common fragile sites and non-fragile regions in red and green, respectively.

Characterization of chromosomal breakpoints

DNA breakage is an important type of cancer-associated genomic aberration, and may cause amplifications, deletions, inversions and translocations. Since array platforms have reduced the ability for detecting inversion and translocation (16,39), the landscape of chromosomal breakage was investigated based on CNA data. In the present study, the genomic start and end of CNAs were defined as breakpoints. These breakpoints may contribute to TNBC initiation and progression. A total of 44,384 CNA breakpoints were identified in 201 TNBC samples. The number of chromosomal breaks per sample ranged between 30 and 709, with a median value of 193. To investigate CNA breakpoint hotspots across the genome, each chromosome was binned into continuous 1 Mb windows, and the density of breakpoints per bin was calculated. Next, the position of CNAs was shuffled 10,000 times and the background distribution of DNA breaks was obtained. The breakpoint-prone genomic regions were identified by comparing the actual value with the background distribution of breakpoints. In total, 199 genomic regions that were significantly enriched for breakpoints of chromosomal rearrangements were identified (Bonferroni corrected P<0.01; Table SIV). To compare these regions with published common fragile sites and non-fragile regions, related data were downloaded from existing literature (46). Among the identified hotspots, only 20 regions (~10%) overlapped with CFS, and the majority of them were located on chromosome 1 (Fig. 4). A total of 25 (~12%) hotspots were located within NFS and were relatively evenly distributed across the genome. Combined, these results possibly reveal the TNBC-specific genomic instability regions.

Subgroup analysis based on CNAs

The full dataset underwent copy number based unsupervised hierarchical clustering. Classification was based on the Euclidean distance following Ward's method, and three main clusters were identified. Results of the cluster analysis are presented as a dendrogram in Fig. 5. Cluster 1 was composed of tumors with extensive arm-level CNAs. The most common alterations include gains of 1q and 8p, and losses of 5q and 8p. These tumors were dominated by frequent CNAs, suggesting an important role of genome instability in the tumorigenesis of this type of TNBC. By contrast, Cluster 2 was characterized by few CNAs, and this cluster may represent M class cancers that were identified in previous studies (38). These DNA copy number stable tumors may be primarily driven by a mutation rather than by CNA. Tumors in Cluster 3 are extraordinarily influenced by small focal alterations rather than by arm-level events, and especially copy number losses. Furthermore, each of these three main clusters can be further decomposed into several smaller sub-clusters.

Figure 5

Unsupervised hierarchical clustering of genome-wide copy number alteration data. The tumor samples are arranged along the x-axis and ordered according to their copy number profiles. The y-axis numbered with 1-22 denotes chromosome numbers. The intensity of red and blue matches the frequency of copy number gain and loss, respectively.

Discussion

TNBC is characterized by its unique molecular profile, aggressive behavior and lack of targeted therapies. Currently, treatment options for TNBC are limited when compared to that of other types of breast cancer. Previous genomic profiling studies have reported that TNBC is a heterogeneous malignancy involving diverse genomic alterations (5). Genome-wide meta-analysis of CNAs will improve the understanding of the mechanisms responsible for TNBC initiation and progression, and facilitate the detection of more reliable tumor markers and the development of targeted therapy. The present study characterized DNA CNAs in TNBC based on published high-resolution genomic array data. From the overall CNA profile of all samples, several frequent recurrent arm-level alterations were observed, which may contribute to the malignant transformation of TNBC. For the analysis of focal CNA events, the GISTIC algorithm was utilized to discern significant CNAs, and several candidate genes were identified, including both known cancer genes and genes that have not previously been associated with any type of cancer.

Among the identified driver genes, some were known cancer genes and were recorded in the COSMIC database (43). For example CCNE1, CD274, epidermal growth factor receptor (EGFR) and JAK2 were located in genomic gain regions, while APOBEC3B, CDKN2A, FAT1 and NOTCH1 were identified in loss regions. Some novel or recently described genes were also detected from the datasets. For example, the VOPP1 WW domain binding protein (VOPP1) gene was identified to be located at the amplified region 7p11.2. It has been reported that VOPP1 is often co-amplified with EGFR or is the breakpoint location for EGFR amplification (49,50). The overexpression of VOPP1 can increase the transcriptional activity of nuclear factor κ B subunit 1 by facilitating its nuclear translocation and associated apoptotic response. VOPP1 overexpression has been observed in several tumor types, such as gastric cancer, head and neck squamous carcinoma and glioma (51-53). A recent study reported that VOPP1 can promote breast cancer development by interacting with the tumor suppressor WW domain containing oxidoreductase (54). These results suggest that VOPP1 may serve an important role in TNBC carcinogenesis and could be exploited to develop therapies for patients with TNBC. Another candidate gene dehydrogenase/reductase 4 like 2 (DHRS4L2) is located in the deletion region of 14q11.2 and is a member of the short-chain dehydrogenases/reductases family. DHRS4L2 produces multiple transcript variants through alternative splicing (55,56). The encoded protein may be an NADPH dependent retinol oxidoreductase. Genomic loss of DHRS4L2 may lead to low expression and has been demonstrated to be associated with risk of diseases (56,57). The DHRS4L2 gene has not yet been reported to be involved in cancer and could be considered a novel candidate gene for TNBC.

Furthermore, the high degree of genomic instability is a hallmark of BRCA1-deficient TNBC. This instability is a prerequisite for the development of large numbers of CNAs, which can affect tumor suppressor genes and oncogenes. Although the BRCA1 status of tumors in the present cohort was not available, the significant numbers of CNAs in these samples were indicative of genomic instability and presented a specific pattern of TNBC.

The present study extensively analyzed the distribution of the chromothripsis phenomenon in TNBC. Chromothripsis is a single catastrophic event that generates dozens of mutations sufficient to produce a malignancy and is distinct from the progressive accumulation of mutations model of cancer development (23-25). The incidence of chromothripsis is heterogeneous across different cancer types. A total of 31 chromothripsis cases were identified out of 201 TNBC samples, with an incidence rate of ~15%. We previously determined that the incidence of chromothripsis in breast cancer is ~11% (36). Thus, TNBC has a relatively higher incidence of chromothripsis and exhibits increased evidence of genomic instability compared with other breast cancer subtypes. Currently, the underlying mechanisms leading to chromothripsis remain largely unknown, although several hypotheses have been proposed, including the formation of micronuclei (58,59), premature chromosome condensation (60), abortive apoptosis (61,62) and breakage-fusion-bridge cycles (63-65). The present study further evaluated the patterns of chromosomal pulverization based on the results, in order to provide clues to identify the underlying mechanisms of chromothripsis in TNBC. In the present cohort, 22% of chromothripsis cases affected two chromosomes. Chromothripsis involving more than one chromosome can result from several chromosomes in a micronucleus or is the consequence of a process of aborted apoptosis (59). More than one mechanism may be responsible for the chromothripsis events in TNBC. Several chromosomal pulverization hotspots were identified across the genome, which may contain critical genes for genomic instability. The most frequently affected chromosome regions were identified, and these results provide a foundation for further analysis. For example, as presented in Fig. 3, there were two chromothripsis events in chromosome 3 and chromosome 12, respectively. Distinct cancer genes were located in the different pulverization regions, which may reflect different mechanisms that trigger chromothripsis events in these two samples.

The present results of the chromosomal breaks analysis provides detailed characterization of recurrent chromosomal breakpoints and affected genes in the TNBC genome. The points of copy number level shift in somatic CNA profiles indicate underlying chromosomal breaks and genomic locations affected by somatic structural aberrations. Thus, the DNA CNA data generated by high-resolution genomic arrays enabled a systematic search for regions and genes that were affected by CNA-associated breakpoints. These breakpoints may silence tumor suppressor genes or create novel gene fusions with oncogenic potential. Next, simulation experiments were performed in which CNA locations were randomly assigned throughout the genome. The simulation was conducted 10,000 times and 199 recurrent breakpoints were identified that were clustered more than would be expected. These breakpoint-prone regions were compared with known common fragile sites and non-fragile regions of the human genome to reveal TNBC-specific genome instability regions (46).

Unsupervised hierarchical clustering analysis was performed to generate a comprehensive view of genome-wide copy number changes in TNBC. Three main clusters were identified. Cluster 1 was predominantly represented by tumors with large chromosomal CNAs, while cluster 3 was primarily influenced by small focal alterations. These findings may indicate the different underlying mechanisms that drive tumors in both groups. For example, the large chromosome or chromosome arm-level CNAs are usually caused by aneuploidy or abnormal numbers of chromosome, while focal alterations often cause gene-level mutagenesis (13). Consistent with previous studies, the present study identified another cluster with few CNAs, which may represent the M class cancer that is predominantly driven by mutations rather than by CNAs (38). Since somatic point mutations cannot be detected by genomic arrays, this cluster presented a different CNA pattern compared with the other two clusters. These observations reveal biologically heterogeneous groups of TNBC, which may represent different molecular mechanisms that underlie tumor development. The clustering of these tumor samples will contribute to the classification and clinical decision making of TNBC.

Previously, several studies have focused on genomic alterations in breast cancer. Banerji et al (66) analyzed 103 breast cancer samples of diverse subtypes using whole exome sequencing, which revealed several recurrent somatic mutations and fusion genes that contributed to breast cancer progression, In addition, Kim et al (67) provided profiles of longitudinal TNBC samples during neoadjuvant chemotherapy, and revealed that resistant genotypes were pre-existing and adaptively selected by therapy. Using single-cell sequencing, Gao et al (68) identified clonal subpopulations in individual tumors that shared a common evolutionary lineage, and demonstrated that most CNAs were acquired at the earliest stages of tumor evolution. These studies provide valuable knowledge about genomic alterations in TNBC. However, the number of tumor samples, particularly TNBC samples, involved in these studies were limited. The present study focused on TNBC and included a larger numbers of samples that resulted in the generation of significant results.

The present study provides valuable new insights into the mechanisms of genomic instability in TNBC. However, there remain a number of questions that need to be studied for TNBC. For example, whether tumors at different stages or metastatic breast cancers demonstrate significantly distinct CNA profiles, and the biological or clinical significance of the observed differences is unknown. Several other types of breast cancer exist, including inflammatory tumors and HER2-positive tumors. The differences in CNA profiles between the TNBC-subtype and other tumor types needs to be further elucidated. Furthermore, the breast cancer genome evolution during disease progression is not yet fully understood yet. In addition, the existence of intra-tumor heterogeneity makes things more complicated. Therefore, in-depth analysis of these questions will help improve understanding of the etiology of TNBC.

In conclusion, a comprehensive characterization of somatic genomic alterations was performed based on a large cohort of TNBC samples. The current study presented several novel findings for TNBC, including: i) A total of 123 regions of significant amplification and deletion were determined; ii) the incidence of chromothripsis in TNBC was identified as ~15%; ii) the distribution and hotspots of CNA breakpoints were revealed; and iii) three tumor clusters and their CNA patterns were identified. The present findings contribute to an increasingly detailed portrait of genomic features of TNBC and may accelerate the rate of driver gene discovery.