Divergent trajectories to structural diversity impact patient survival in high grade serous ovarian cancer

Abstract

Deciphering the structural variation across tumour genomes is crucial to determine the events driving tumour progression and better understand tumour adaptation and evolution. High grade serous ovarian cancer (HGSOC) is an exemplar tumour type showing extreme, but poorly characterised structural diversity. Here, we comprehensively describe the mutational landscape driving HGSOC, exploiting a large (N = 324), deeply whole genome sequenced dataset. We reveal two divergent evolutionary trajectories, affecting patient survival and involving differing genomic environments. One involves homologous recombination repair deficiency (HRD) while the other is dominated by whole genome duplication (WGD) with frequent chromothripsis, breakage-fusion-bridges and extra-chromosomal DNA. These trajectories contribute to structural variation hotspots, containing candidate driver genes with significantly altered expression. While structural variation predominantly drives tumorigenesis, we find high mtDNA mutation loads associated with shorter patient survival. We show that a combination of mutations in the mitochondrial and nuclear genomes impact prognosis, suggesting strategies for patient stratification.

Fig. 1. Structural complexity in HGSOC.

A Structural variants are abundant across the combined tumour cohort (N = 324), reflected in the frequencies of SV and CNA calls and particularly of translocations and complex events. Samples also show frequent HRD and WGD across tumour stages and cohorts (Supplementary Data 1). Pathogenic germline variants are relatively enriched in known HGSOC susceptibility genes. B Of all protein-coding genes in the genome (N = 19,329), most are impacted by each class of variation at the cohort level (bar chart on left) but the numbers of disrupted genes per sample (forest plot on right; mean ± SD) are dominated by SVs and CNAs and impact a smaller number of genes per sample. Variants predicted to disrupt function are nonsynonymous SNVs of high/moderate impact and SVs or CNAs overlapping ≥1 protein-coding exon. C The most frequently disrupted genes overall are not enriched for known cancer genes, tend to be longer than average and to intersect recurrent CNA hotspots. The most frequently disrupted gene is TP53, where there is near ubiquitous mutation. NAALADL2 is an unusually long (1.37 Mb) gene and is located within a common fragile site frequently altered across many tumour types. (All somatic variant calls predicted to disrupt protein-coding genes are provided in Source Data.).

Fig. 2. CNAs and SVs form hotspots throughout the genome, which reveal candidate driver genes.

A GISTIC enrichment analysis reveals CNA hotspots. Duplication peaks are green, and deletion peaks are blue. B CNA hotspot activity varies across samples, with greater activity at deletion hotspots in HRD samples and greater duplication activity in WGD samples. C SVs are enriched on chromosome 19 and peak at CCNE1. In particular, there is a relative excess of inversions, translocations, breakage-fusion bridges and chromothripsis implicating this gene. D Cancer gene census genes in SV/CNA hotspots that are differently expressed across samples in the presence of the event type driving the hotspot. Genes in CNA deletion hotspots or hotspots of breakpoint enrichment have significantly lower expression, while genes in duplication (SV or CNA) or inversion hotspots have significantly higher expression. CEP89 SV deletions and SV duplications are both associated with increased expression, reflecting the high level of overlapping genomic instability occurring at this locus, which includes CCNE1. CNA deletions at LEPROTL1 are also associated with reduced expression, but this likely represents the same signal as that observed at nearby WRN. Expression fold changes are log2 transformed and are robust to the percentage of samples with the hotspot event. SV and CNA hotspot data are in Supplementary Data 7 and 8, respectively.

Fig. 3. Diverse somatic mutation classes underlie HGSOC candidate driver genes.

A Combined oncoplot indicates predicted driver genes (Methods) subject to recurrent pathogenic SNV/SV/CNA mutations predicted to impact function (SNVs annotated are nonsynonymous variants of HIGH or MODERATE impact by VEP; SVs or CNAs overlapping ≥1 exon). The unbracketed percentage is the percentage of patients with a predicted pathogenic SNV/SV or CNA. For SNVs/SVs, this is either an SNV or a deletion, or in the case of the asterisked rows, representing genes associated with gain of function, it is an SNV or duplication. The bracketed percentage is less conservative and is the percentage of patients with any annotated event. The vertical bar plot represents the mutational burden of predicted deleterious SNVs/SVs. B The horizontal bar plot (left) represents the proportion of patients with different types of pathogenic driver mutations. The forest plot (right) represents the mean number of each type of driver mutation across tumours with at least one event and the standard deviation (whiskers), based on N = 324 patients. C Total mutation frequencies in candidate driver genes. The proportions of patients with somatic alterations of any kind in each gene, whether predicted to be pathogenic or not. Clustered SVs are members of SV breakpoint clusters from ClusterSV; the Multi-Hit category represents patients with >1 somatic alteration in the gene. (All somatic variant calls predicted to disrupt protein-coding genes are provided in Source Data.).

Fig. 4. Patterns of co-occurrence and mutual exclusivity of complex structural variant classes and their relationship with WGD and HRD.

A Abundance and co-occurrence of all complex events reveals two clusters of events defined by negative correlation between HRD and WGD and other cSV types. B Complex SVs are depleted in HRD samples (purple; N = 181) with the exception of chromoplexy. Breakage fusion bridge cycles, ecDNA and chromothripsis are enriched in whole genome duplicated samples (N = 158). Points represent log odds ratios, and lines represent 95% confidence intervals. C Biased co-occurrence and mutual exclusivity of cSV classes support divergent tumour evolutionary trajectories with significant association seen between HRD and chromoplexy, but exclusivity between HRD and all other features. Cell counts represent observed vs expected counts of event co-occurrence in samples. D Association between the presence of chromothripsis and breakage fusion bridges (x-axis) and the estimated fraction of the genome duplicated (y-axis) across all samples. Enrichments (two-sided Fisher’s exact tests) in WGD samples of: chromothripsis (OR (95% CI) = 2.1 (1.3, 3.3); adj. p-value = 5.12 × 10⁻³); and breakage fusion bridges (OR (95% CI) = 2.3 (1.4, 3.8); adj. p-value = 5.12 × 10⁻³). HRD, WGD and complex SV data are in Supplementary Data 1.

Fig. 5. Deleterious mtDNA mutation loads are a potent biomarker of overall survival.

A Somatic (inner ring) and germline (outer) SNV frequencies across the cohort in mitochondrial encoded genes (black: single nucleotide variants; red: indels). B Abundant somatic SNVs disproportionately impact protein-coding genes in mtDNA (758 SNVs in coding genes of 1245 SNVs in total). C SNVs categorised by VEP functional impact annotation include many protein-altering variants expected to alter mitochondrial complex functions. D Overall survival Cox proportional hazard (PH) ratio (95% confidence intervals in blue) increases with increasing heteroplasmy of deleterious mtDNA mutations (Cox PH p-value = 0.0002; N = 277 patients, stratified by cohort and adjusted for age, stage and HRD status). E Overall survival Cox PH ratio (95% confidence intervals in blue) is stable with increasing heteroplasmy of synonymous mtDNA mutations (Cox PH p-value = 0.91; N = 277 patients, stratified by cohort and adjusted for age, stage and HRD status). mtDNA SNV data are in Supplementary Data 13.

Fig. 6. Multivariable modelling of the impact of genomic features of HGSOC on overall survival adjusted for baseline clinical factors.

A Univariable modelling of 33 genomic features using a Cox Proportional Hazards (PH) model adjusted for HRD, age at diagnosis and stage at diagnosis and stratified by cohort (N = 277 samples with complete overall survival time and tumour stage data). Forest plot shows log hazard ratios (HRs) and 95% confidence interval (CI) per feature; p-values adjusted for multiple testing, the null hypothesis HR = 1. B As (A) but HRs from one multivariable model including all 33 genomic features plus adjustments. C HRs with 95% CI for selected features from elastic net penalised Cox PH regression model (round points coloured as in B), versus boxplots for HR estimates (centre = median, box = 25% and 75% percentiles, whiskers = 1.5*IQR) from 10 cross-validations of elastic net for the same features. Kaplan–Meier (K–M) plots (N = 279 patients) show overall survival with 95% CI for presence (green curve) and absence (purple curve) of D CDK12 SNV, E deleterious mtDNA mutation F severe chromothripsis (>2 chromosomes affected by chromothripsis). Cox PH HRs and 95% CI are reported for each K–M plot, with p-value rejecting HR = 1 (N = 279 patients). Survival data are in Supplementary Data 12 and 14.