{'fraction_non_alphanumeric': 0.03470579717140512, 'fraction_numerical': 0.024279577158855047, 'mean_word_length': 5.107900943396227, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 2, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Whole-genome-sequencing of triple negative breast cancers in a population-based clinical study

Summary Whole genome sequencing (WGS) brings comprehensive insights to cancer genome interpretation. To explore clinical value of WGS, we sequenced 254 triple negative breast cancers (TNBC) with associated treatment and outcome data collected between 2010-2015 via the population-based Sweden Cancerome Analysis Network-Breast (SCAN-B) project (ClinicalTrials.gov ID:NCT02306096). Applying the HRDetect mutational-signature-based algorithm to classify tumors, 59% were predicted to have Homologous-recombination-repair deficiency (HRDetect-high): 67% explained by germline/somatic mutations of BRCA1/BRCA2, BRCA1 promoter hypermethylation, RAD51C hypermethylation or biallelic loss of PALB2. A novel mechanism of BRCA1 abrogation was discovered via germline SINE-VNTR-Alu retrotransposition. HRDetect provided independent prognostic information, with HRDetect-high patients having better outcome on adjuvant chemotherapy for invasive disease-free survival (Hazard Ratio, HR=0.42, 95% confidence interval, CI=0.2-0.87), and distant relapse-free interval (HR=0.31, CI=0.13-0.76) compared to HRDetect-low, regardless of whether a genetic/epigenetic cause was identified. HRDetect-intermediate, some possessing potentially targetable biological abnormalities, had poorest outcomes. HRDetect-low cancers also had inadequate outcomes: ~4.7% were mismatch-repair-deficient - another targetable defect, not typically sought; and was enriched for (but not restricted to) PIK3CA/AKT1 pathway abnormalities. New treatment options need to be considered for now-discernible HRDetect-intermediate and HRDetect-low categories. This population-based study advocates for WGS of TNBC to better inform trial stratification and improve clinical decision-making.


Prior germline testing and classification of BRCA1 and BRCA2 germline variants
49 patients had prior clinical genetic counseling involving NGS-based screening of BRCA1 and BRCA2, or were enrolled in the SWEA research study (The Swedish BRCA1 and BRCA2 study collaborators (SWE-BRCA) Extended Analysis) for high-risk patients and screened for an extended panel of susceptibility genes. The inclusion criteria for the SWEA study were in line with the Swedish national clinical practice guidelines for breast cancer. Briefly, genetic testing was offered when there was at least a 10 % probability to detect a pathogenic germline variant in BRCA1 or BRCA2, based on the patient's age at diagnosis, histology, and family history. Detected germline variants were classified according to the ENIGMA BRCA1/2 Gene variant Classification Criteria (2017-06-29) https://enigmaconsortium.org/library/general-documents/. Only class 5 variants were considered as pathogenic, corresponding to nine BRCA1 and three BRCA2 variants.

DNA promoter methylation analysis
Bisulfite conversion of genomic DNA was performed with the column based EpiTecht Fast DNA Bisulfite kit (Qiagen GmBH, Hilden, Germany. Promoter methylation analysis was performed using a PSQ MD 96 pyrosequencing instrument (Qiagen). A fully methylated as well as an unmethylated sample was included as controls in each run. The PyroMark analysis program was used for data analysis and all electropherograms were manually checked. For BRCA1, analysis was performed as described 6 , and included analysis of two CpG island regions. We used 7% cut-off compared to previous 10% due to cleanness of the data.
For RAD51C we first performed a correlation analysis based on TCGA breast cancers used in a previous study 7 similar to the study by Polak et al. 8 . This analysis identified seven CpGs present on the Illumina 450K methylation beadchip array <1000bp upstream of the gene. Four CpGs with a Pearson correlation less than -0.44 (cg05214530, cg27221688, cg02118635, cg10487724) were selected, from which primers for RAD51C were adapted from Hansmann et al. 9 .
Final RAD51C primers were: In the final analysis, two different primer sets of four (BRCA1 and RAD51C) different CpG sites were used for each gene. CpG allele methylation percentage was averaged across each primer set and next merged to the mean of the two sets. Cut-offs were applied (BRCA1 7%, RAD51C 9%) for making a call on methylation or otherwise. The cut-off was set higher for RAD51C as this assay generated slightly higher background variability between primer sets. The cut-off was verified (and supported) by RNA sequencing data for RAD51C. Two cases failed RAD51C methylation analysis.
For PALB2 and RAD51 promoter hypermethylation analysis primers described by Wanatabe et al. 10 were used. A similar cut-off (9%) as for RAD51C was used to call methylated cases (none observed).

IHC validation of mismatch repair deficient (MMRd) cases
Suspected MMRd cases identified by whole genome sequencing were stained for MLH1, PMS2, MSH2 and MSH6 as outlined originally in Joost et al. 11 . Stained slides were evaluated for protein expression in both tumor cells and non-malignant cells.

Gene expression analyses
Gene expression data was available from Gene Expression Omnibus 4 , series GSE96058 reported elsewhere 5 . FPKM data for specific genes were extracted and log2 transformed. For TNBCtype, IC10, CIT classification, 228 cases were available for analysis based on GSE96058 5 (primary tumors only). For remaining cases, these were included separately and subtyped only using AIMS 12 (as this is a single sample predictor of molecular subtype) and analyzed for individual FPKM gene expression.
Classification according to different molecular subgroups in breast cancer was performed as follows, after i) an offset of 1 was added to all FPKM values, ii) log2 transformation: -PAM50. PAM50 subtypes were obtained using the AIMS single sample classifier 12 , based on the aims R package. All samples were classified. -TNBCtype 13,14 . For TNBCtype classification the entire GSE96058 data set was used.
Data was mean-centered across all samples for each gene, TNBC cases were extracted and uploaded as a separate data set into the web-based classifier 14 . For a few cases the web-based application called the samples as not being ER-negative. These samples were removed from the TNBC data set (inferring missing values) and remaining samples were again uploaded to the web-based application for subtyping. -IC10 15 . For IC10 classification the entire GSE96058 data set was used. Data was mean-centered across all samples for each gene. IC10 subgroups were obtained through the ic10 R package using default processing. -CIT 16 . CIT subtypes were obtained through the citbcmst R package, using pearson correlation as distance method and gene symbol as matching entity.
Calculation of six gene expression metagenes representing different biological functions in breast cancer was performed as described 17 , using the GSE96058 data set for gene-centering across samples.

Unsupervised clustering
All unsupervised analysis was performed in R 18 using the ConsensusClusterPlus R-package 19 . Two input formats were used, FPKM data and PCA components. When using FPKM data as input this was first offset by +0.1, log2 transformed, and mean-centered across samples for each RefSeq associated gene. A filter step based on standard deviation of expression was used as defined in result presentations. In the clustering we used Pearson correlation as distance metric and ward.D2 linkage. Additional parameters were pItem=0.8, pFeature=0.8, number of iterations = 2000.
When using PCA components as input, these were derived from a PCA analysis of all 19000 RefSeq genes available using the prcomp() function in R. All components were then used in consensus clustering with pItem=0.8, pFeature=0.98, number of iterations = 2000.

Machine learning
All machine learning was performed in R 18 using the Caret package. Model performance was assessed using ROC analysis and the area under the curve (AUC) estimate. We used 70% of samples for training and 30% for internal validation. Division into training and internal validation was performed using the createDataPartition() Caret function, balancing the sample splitting for grade (1,2,3), lymph node status (node-negative, node-positive), and age (<50, >=50 years). For FPKM data as input this was offset by +0.1, log2 transformed and mean-centered across samples in respective cohort (training and test). For machine learning using PCA components (see above for extracting these) all components were used directly in training.
We evaluated seven machine-learning methods listed below by their names in Caret: -svmLinear (linear support vector machine) In the training we used the trainControl() Caret function with parameters: method=repeatedcv number=4 repeats=10 classProbs=TRUE summaryFunction=twoClassSummary savePredictions=TRUE For the actual training we used the train() Caret function with ROC as metric, the defined trainControl parameters above, and tune.length=10. Trained models were saved as R data objects. We applied the models to the internal validation set for each evaluation group using the predict() Caret function. ROC values were calculated using the roc() function in the pROC R package.
The entire process described above was repeated 10 times with different splits of training and internal validation sets to reduce the possibility of bias in the results based on sample splitting. This generated 10 different AUC estimates from 10 potentially different predictors for each model.

Survival analyses
Definition of clinical endpoints: --Overall survival was obtained from national registries, calculated as the time from diagnosis to death of any cause.
--Invasive disease-free survival (IDFS) was defined according to STEEP guidelines 20 , as the time from diagnosis to either death of any cause or invasive breast-cancer related events (loco-regional and distant recurrence). --Distant relapse-free interval (DRFI) was defined according to STEEP guidelines as the time from surgery to diagnosis a distant relapse (event) or to last day of follow-up (censoring). Events include patients that first developed a loco-regional relapse, and then a distant relapse. For these patients the day of the distant relapse was used.
Exclusion criteria for outcome analyses: --Neoadjuvant treatment --Metastatic disease at time of diagnosis (including microinvasive disease).
--Metastatic disease identified immediately prior to, or during adjuvant chemotherapy.
--Patients not managed in an adjuvant setting (irrespective if adjuvant treatment or not provided later). --Bilateral breast cancer.
--Lost to follow-up before start of systemic treatment.
--For DRFI, patients with a relapse or death from a malignancy of uncertain origin were excluded. These patients were however included in OS and IDFS analyses.

Multivariable analyses
Analysis was performed using the coxph R function from the survival R package. Covariates in multivariable Cox regression were patient age (<50, ≥50 years), lymph node status (N0/N+), tumor size (≤20, >20mm), and tumor grade (1,2,3). Data for lymph node status, and tumor size were obtained from NKBC data. Tumor grade was obtained from clinical review of individual patient's files.

Whole Genome Sequencing Analysis
Whole genome sequencing and alignment 150 base pair paired-end sequencing was performed on Illumina X10 machines following standard library preparation according to the manufacturer's protocols. Target insert size was 450 bp. The target coverage for tumour-normal pairs was 30x for patients that received adjuvant chemotherapy (based on NKBC registry data, see above), and 15x otherwise. Pairedend reads were aligned to the reference human genome (GRCh37) using Burrows-Wheeler Aligner, BWA (v0.7.15).

Identification of somatic mutations
Paired tumour-normal bam files were interrogated for somatic mutations using the following algorithms: • . ANNOVAR: Functional annotation of genetic variants from next-generation sequencing data. 23 • All rearrangements had to pass the stage of re-assembly in Brass, and were further filtered by size (at least 1kb between breakpoints). A known library-preparation artefact was removed by filtering out inversions that were shorter than 5kb and reported by 5 or fewer reads.
For seeking cancer driver mutations, we applied less stringency in specificity.

Criteria for sample exclusion from analysis
DNA was extracted from diagnostic cancer samples in the population-based study as described above (section Tissue sampling, DNA and RNA extraction). For some patients it is possible that the amount of tumor DNA present in the sample was insufficient for wholegenome sequencing and interpretation of tumour-specific somatic mutations.
Samples were excluded from analysis if: -the library preparation failed to yield sufficient output -for samples with coverage below 20x, Battenberg algorithm failed to find a copynumber solution -for samples with coverage over 20x, Battenberg failed to find a copy-number solution and the the number of somatic point mutations identified was less than 500.

Identification and curation of germline and somatic mutations in genes related to DNA repair
We interrogated germline and somatic mutations in genes related to DNA repair. The list of 165 genes included genes tested in the SWEA study (see above in section Prior germline testing and classification of BRCA1 and BRCA2 germline variants) and genes included in the INSIGNIA study (experimental study to explore mutagenesis in isogenic cell systems and in patients with inherited DNA repair defects, www.mutationsignatures.org). Screened genes are listed in the Supplementary Data file.
In order to prioritise germline mutations which are likely pathogenic, the variants were filtered to ensure the mutations were real (supported by at least 5 reads and not in a region of repetitive sequence, by curation), and further filtered according to prior information on particular mutations (protein consequence, frequency in population through the ExAC database less than 1 in 100, and if considered pathogenic in the ClinVar database, dated 20170501). Nonsense, frameshift or essential splice-site mutations were considered pathogenic, unless they were known to be exceptions (eg. p.K3326* BRCA2 nonsense).
Germline copy number in most commonly-mutated DNA-repair genes (BRCA1, BRCA2, PALB2) were assessed through visual inspection of coverage in the bam files.
Among somatic mutations in the DNA repair genes, as pathogenic we considered those that were either nonsense, frameshift or essential splice site. Somatic rearrangements from BRASS algorithm were also included if a breakpoint disrupted the gene body.
Finally, somatic copy number at the DNA-repair genes was inferred from the genome-wide copy number profiles from the ASCAT algorithm. When total copy number is zero in the tumour, the gene was considered lost (homozygous deletion), and when it was one, it was flagged as loss of heterozygosity.

Identification of cancer drivers
Among point mutations and small insertions and deletions, only coding and essential splice site mutations were considered as potential drivers. To be classified as putative drivers, missense or in-frame variants in dominant genes had to either be recurrent in literature (e.g. COSMIC database) or be deleterious to previously reported recessive cancer genes. Mutations in 740 genes identified in literature as cancer drivers were considered, together with classification regarding dominant or recessive status of each gene.
Copy number profiles of the cancer genomes were scrutinised for changes to common cancer genes. The ASCAT algorithm provided estimates of overall tumour ploidy, cellularity and copy number profiles. Copy number changes were considered putative drivers when they coincided with frequently altered genes 24 . Specifically, copy number changes were deemed as an amplifications when total copy number in the tumour exceeded (p*2+1), where p stands for overall tumour ploidy. Homozygous deletions were declared when total copy number in tumour was zero. Loss of heterozygosity was determined when total copy number of a gene was one and minor copy number was zero.
Finally, previous analysis 24 revealed increased frequencies of somatic rearrangements in gene bodies of several recessive cancer genes, compared to the surrounding chromosomal regions. We declared some somatic rearrangements as putative drivers when they disrupted gene bodies, with the exception of instances when both breakpoints of rearrangements hit the same intron. The list of recessive cancer genes included: ARID1B, CDKN2A, FBXW7, MAP2K4, MAP3K1, PTEN, RB1, TP53, and MLLT4.

Application of the HRDetect algorithm
A high-specificity somatic mutation dataset was first compiled for all cancer genomes in the study as described in above (section Identification of somatic mutations).
Cellularity and ploidy estimates were obtained from the output of ASCAT.
Point mutations were annotated with the sequence context (and reverse-complemented when necessary to put them into pyrimidine context). Signatures of point mutations were assigned with the SigFit R package using the 30 signatures from the COSMIC signature database (https://cancer.sanger.ac.uk/cosmic/signatures). Signatures of rearrangements were also analysed by SigFit using the six rearrangement signatures from Nik-Zainal et al. 24 . Small deletions were classified into the following categories: micro-homology mediated, repeatmediated and others. HRD score was calculated from ASCAT copy number profiles as described previously 24 .
Exposures of signatures, HRDetect score, and proportion of indels at microhomology were tabulated for the entire cohort. The features were normalised as described in the original report describing HRDetect 25 . The HRDetect weights were applied to the normalised features, to obtain HRDetect probability for every patient.

Principal-component analysis of components of HRDetect
HRDetect assigns a single probability of HR deficiency to every sample. However, we have previously demonstrated that BRCA1-null and BRCA2-null tumors were distinguishable. To explore the sub-structure of cancers with HR deficiency, we conducted a Principal-component analysis (PCA) of components of HRDetect.
Specifically, the products of normalised features and HRDetect weights were tabulated as matrix N x 6, where N is the number of samples in the study and 6 represents the 6 non-zero weights of HRDetect. Each entry in the matrix is a product of HRDetect weight and the respective normalised sample feature. Prior to the PCA analysis, the variables in the new matrix were not scaled, so that the dimensions with highest weights of HRDetect exhibited highest variance.
Two-dimensional visualisations of the PCA results are presented in Figure 3.

Mobile element (MELT) analysis
Mobile element (ME) insertions not present in the GRCh37 reference genome were called and genotyped using the Mobile Element Locator Tool (MELT) version 2.1.5 26 . The ME reference sequences provided with MELT were used for discovery of three major classes of MEs: Alu, SVA (SINE-VNTR-Alu) and LINE1 (long interspersed nucleotide element-1). Default parameters were used except -z (set to 50000) and read length, coverage and insert length, which were determined for each bam file. Tumor and normal samples from the same individual were analyzed together to increase sensitivity (MELT GroupAnalysis and MakeVCF). ME insertions previously discovered as part of the 1000 Genomes Projects, Phase III, were used as priors 27,28 . Supporting aligned read data for ME insertions discovered in or near (within 5 kb) 11 genes with potential impact on HR (ATM, BARD1, BRCA1, BRCA2, BRIP1, CHEK2, MRE11, NBN, PALB2, RAD51C, RAD51D) were visually inspected in Integrative Genomics Viewer (IGV 29 ).
On average, 1302 Alu, 186 LINE1 and 90 SVA insertions were detected in each pair of tumor and normal sample. Of these, 90.5% of the Alu, 79.5% of the LINE1 and 86.3% of the SVA insertions were previously seen in the 1000 Genomes Projects. Most insertions were detected in both the normal and tumor data (93.1% of Alu, 89.3% of LINE1 and 86.3% of SVA elements).
Two ME insertions were found in the 11 HR related genes. One was a known Alu insertion in BRCA2 at chr13: 32922439 found in two sample pairs (ID: ALU_umary_ALU_9673, detected in both normal and tumor). This Alu insertion has an allele frequency of 4.2% in the 1000 genomes project, phase III dataset (ref 3). The other was an SVA element inserted 1.8kb downstream of the first coding exon in BRCA1 at chr17:41274217. This ME was not seen in the 1000 genomes project, phase III dataset.