Next-generation sequencing of BRCA1 and BRCA2 genes for rapid detection of germline mutations in hereditary breast/ovarian cancer

Patients and DNA

A retrospective cohort of 11 familial breast cancer patients, training set (TS), previously tested by Sanger sequencing for all BRCA1 and BRCA2 exons and by MLPA to detect large genomic rearrangements (LGRs), and a prospective consecutive cohort of 136 uncharacterized probands, validation set (VS), enrolled at the Hereditary Tumors section of Policlinico Umberto I, University La Sapienza, between July 2015 and September 2017, were analyzed. Selection criteria were as follows: (i) three or more breast cancer cases diagnosed at any age or two first-degree family members affected before 50; (ii) early onset breast cancer (<35 years); (iii) breast and ovarian cancer in the same individual or two breast cancer cases and at least one ovarian cancer case, or one breast cancer case and one ovarian cancer case diagnosed before 50 in first-degree family members; (iv) male breast cancer; (v) BRCAPro5 greater than 10% (Capalbo et al., 2006a, 2006b; Coppa et al., 2014).

Comprehensive pre-test counseling was offered to all probands and their family members and informed consent were obtained. The TS was characterized by 11 samples carrying 11 different pathogenic mutations including single nucleotide variants (SNVs), indels, and LGRs plus 37 benign variants (Table 1; Table S1).

Genomic DNA was extracted from peripheral blood samples using a commercial kit (QIAamp Blood Kit, Qiagen, Valencia, CA, USA) according to the manufacturer’s instructions. All investigations were approved by Ethics Committee of the University of Roma “La Sapienza” (Prot.: 88/18; RIF.CE:4903, 31-01-2018) and conducted according to the principles outlined in the declaration of Helsinki.

Ion Torrent PGM library preparation

The entire coding regions of the BRCA1 and BRCA2 genes were amplified using the Ion AmpliSeq BRCA1 and BRCA2 Panel (Life Technologies, Carlsbad, CA, USA) consisting of three primer pools, covering the target regions in 167 amplicons, including all exons and 10–20 bp of intronic flanking sequences, for both genes. Ion Torrent adapter-barcode ligated libraries were generated using Ion AmpliSeq™ Library Kit 2.0 (Applied Biosystems, Houston, TX, USA; Life Technologies, Carlsbad, CA, USA) and Ion Xpress™ Barcode Adapter 1-16 Kit (Applied Biosystems, Houston, TX, USA; Life Technologies, Carlsbad, CA, USA) according to the manufacturer’s procedures to obtain 200-bp PCR fragments flanked by adaptor and barcode sequences, allowing sequencing and sample identification, respectively. In brief, DNA concentration was measured with Qubit™ dsDNA HS Assay Kit (Thermo Fisher Scientific, Waltham, MA, USA) and Qubit® Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA). Ten ng of each genomic DNA sample was amplified, using Ion AmpliSeq™ HiFi Master Mix (Applied Biosystems, Houston, TX, USA; Life Technologies, Carlsbad, CA, USA) and Ion AmpliSeq BRCA1 and BRCA2 Community Panel primer pools (3 pools), for 2 min at 99 °C, followed by 19 two-step cycles of 99 °C for 15 s and 60 °C for 4 min, ending with a holding period at 10 °C in a GeneAmp PCR System 9700 thermocycler (Applied Biosystems, Houston, TX, USA). The amplified DNA samples were digested, then barcodes and adapters were ligated to the amplicons as previously described (Belardinilli et al., 2015) and reported in Ion AmpliSeq library preparation manual. Adaptor ligated amplicon libraries were purified with the Agencourt AMPure XP system (Beckman Coulter Genomics, Danvers, MA, USA). Following Agencourt AMPure XP purification, the concentration of the libraries was determined using the Qubit™ dsDNA HS Assay Kit (Thermo Fisher Scientific, Waltham, MA, USA) and Qubit® Fluorometer (Thermo Fisher Scientific, Waltham, MA, USA). After quantification, each amplicon library was diluted to 50 pM and the same amount of the 12 libraries for four patients (Chip 314v2) or 24 libraries for eight patients (Chip 316v2) was pooled to perform the emulsion PCR reaction and the sequence reaction.

Emulsion PCR and sequencing

The emulsion PCR was carried out with the Ion OneTouch 2 System and Ion Torrent PGM Hi-Q OT2 Kit (Thermo Fisher Scientific, Waltham, MA, USA) according to the manufacturer’s procedure. Template-positive Ion Sphere™ Particles were enriched with the Dynabeads MyOne Streptavidin C1 Beads (Thermo Fisher Scientific, Waltham, MA, USA) and washed with Ion OneTouch Wash Solution. This process was performed using an Ion OneTouch ES system (Thermo Fisher Scientific, Waltham, MA, USA). After the Ion Sphere Particle preparation, Massive Parallel Sequencing was carried out with Ion Torrent PGM sequencer system (Thermo Fisher Scientific, Waltham, MA, USA) using the Ion Torrent PGM Hi-Q Sequencing Kit and Ion 314v2 or 316v2 Chip (Thermo Fisher Scientific, Waltham, MA, USA), according to the established procedures.

Data analysis

Raw sequence data analysis, including base calling, demultiplexing, alignment to the hg19 human reference genome (Genome Reference Consortium GRCh37), was performed using the Torrent Suite Software version 3.6 and subsequent versions up to 5.2 (Thermo Fisher Scientific, Waltham, MA, USA). The average depth of total coverage was set at minimum 500× and for variant calls at minimum of 100×.

The TS data were processed using Coverage Analysis and Variant Caller (germline low stringency) plugins available within the Torrent Suite Software version 3.6 using default analysis settings and subsequently using the version 3.6 of json parameters file associated with the BRCA1 and BRCA2 Ion Ampliseq Panel. The VS data were processed using Coverage Analysis and Variant Caller plugins available within the Torrent Suite Software versions from 4.0 to 5.2 and the json parameters file associated with the BRCA1 and BRCA2 Ion Ampliseq Panel.

All the variants were visually examined and verified using integrative genomics viewer (IGV) version 5.01 (http://www.broadinstitute.org/igv) as well as by filtering out possible strand-specific errors (i.e. mutations detected exclusively in the “plus” or in the “minus” strand, but not in both DNA strands). All mutations are reported following the human genome variation society guidelines (http://www.hgvs.org/mutnomen/) on the basis of the coding sequences NM_007294.3 and NM_000059.3 for BRCA1 and BRCA2, respectively. Common polymorphisms, which represent 5% in the general population, and variant of uncertain significance (VUS) or pathogenic variants were classified referring to the following databases: Breast Cancer Information Core BIC (https://research.nhgri.nih.gov/), Clinical Variants (https://www.ncbi.nlm.nih.gov/pubmed), dbSNP138 (https://www.ncbi.nlm.nih.gov/projects/SNP/), ARUP (http://arup.utah.edu/database/BRCA/), Universal Mutation Database (http://www.umd.be/BRCA1/, http://www.umd.be/BRCA2/), Leiden Open (source) Variation Database (LOVD) (http://www.lovd.nl/3.0/home). Unclassified variants were evaluated by the following in silico predictors: SIFT-PROVEAN (http://provean.jcvi.org/index.php), Polyphen (http://genetics.bwh.harvard.edu/pph2/), GVGD (http://agvgd.hci.utah.edu/agvgd_input.php), Mutation taster (http://www.mutationtaster.org/), Fruit Fly Splice Predictor, NNSPLICE (http://www.fruitfly.org/seq_tools/splice.html), Human Splicing Finder (http://www.umd.be/HSF3/), Splice Predictor (DK), NetGene2 (http://www.cbs.dtu.dk/services/NetGene2/), MaxEntScan (http://genes.mit.edu/burgelab/maxent/Xmaxentscan_scoreseq.html), ESEfinder3.0 (http://rulai.cshl.edu/cgi-bin/tools/ESE3/esefinder.cgi?process=home). The meaning of the variants not yet reported in the literature and in the databases was established as described by Giannini et al. (2008), according to ENIGMA criteria (https://enigmaconsortium.org/).

Concordance of the NGS results with respect to those obtained with Sanger sequencing was calculated using classifications as true positive (TP), true negative (TN), false negative (FN), or false positive (FP). Regarding specificity and negative predictive value (NPV), since these parameters are calculated considering TNs, we considered the number of wild-type bases that were called by NGS and aligned to the reference sequence.

Sanger sequencing

A subset of clinical samples (11 samples of TS and 60 samples of VS) were sequenced for the entire coding regions by Sanger sequencing. Sequencing was performed using an ABI PRISM DyeDeoxy Terminator Cycle Sequencing Kit and an ABI 3100 Genetic Analyzer (Applied Biosystems, Warrington, UK) according to Coppa et al. (2018).

The reference sequence for BRCA1 was NM_007294.3, and the reference sequence for BRCA2 was NM_000059.3.

MLPA assay

BRCA1 and BRCA2 genomic rearrangements were searched by the MLPA methodology (Schouten et al., 2002) according to the manufacturer’s instructions (MRC-Holland, Amsterdam, The Netherlands) and as described (Buffone et al., 2007; Coppa et al., 2018). Variations in peak height were evaluated comparing each sample with a normal control and by a cumulative comparison.

Coverage data analysis

The BRCA1 and BRCA2 panel used in this study generates 167 amplicons that cover all targeted coding exons and exon-intron boundaries. The data coverage in 147 patients (TS and VS) showed a mean amplicon reading depth per sample ranging from 244× to 2,236×, loading 2/4 samples on chips 314v2 or eight samples on chip 316v2. An adequate number of reads mapped onto target regions has been obtained both with chip 314v2 or chip 316v2, with a mean sequencing depth and an average uniformity of coverage of 726× and 97% and 1,504× and 97%, respectively (Table 2; Table S2 and S3). The large increase in sequencing efficiency observed with the 316v2 chips makes the analysis of a better quality and led to a lower number of false-positive indel calls. For variant calls, we chose the threshold of a minimum average depth of 100×, although a minimum average of 80× is normally required for germline mutations (Trujillano et al., 2015) and a minimum average depth in strand bias of 30× for each strand end to end. By a critical analysis of the BRCA1 and BRCA2 coverage data performed on 147 samples (TS and VS) emerged that some amplicons (herein defined as lower-performance amplicons) often did not reach a minimum coverage of 100× or showed fwd/rev end to end unbalance (Table 3; Fig. 1). In particular, the coverage of the AMPL225438570 BRCA1 amplicon was <100× in 61% of samples (Fig. 1A). This amplicon maps in the exon 2 of the BRCA1 gene, containing a large A-T reach region, likely to be responsible for lower amplification efficiency. Concerning BRCA2, we observed five amplicons with a coverage less than 100×, whose GC content, self-annealing/hairpin formation and GC clamp issues could be responsible for low performance (Table 3; Fig. 1B). The same issues may also be responsible for fwd/rev end to end unbalance observed in the following amplicons: AMPL223735053, AMPL223413081, AMPL223712774 in BRCA1 (Figs. 1C and 1D; Table 3); AMPL224626553, AMPL223938117, AMPL223512592, AMPL225349438, AMPL223730984, AMPL223515418, AMPL225441321, AMPL223959719 in BRCA2 (Figs. 1E and 1F; Table 3). Of note, the read depths of these amplicons did not increase proportionally with an increase in mean amplicon read depths (across 167 amplicons) per sample. No variants have been identified in lower-performance amplicons in our TS and VS samples, using Ion AmpliSeq BRCA1 and BRCA2 Panel (Life Technologies, Carlsbad, CA, USA). However, all amplicons with strand bias were carefully and visually verified with IGV and re-analyzed by Sanger sequencing that confirmed the absence of variants.

Detection of BRCA1 and BRCA2 variants

The TS, used to test the performance of the Ion AmpliSeq BRCA1 and BRCA2 Panel (Life Technologies, Carlsbad, CA, USA) on Ion Torrent PGM platform, contained 6 BRCA1 and 5 BRCA2 pathogenic variants previously identified (Table 1) plus 37 benign variants known (Table S1). Data obtained by Ion Torrent PGM sequencing were blindly analyzed by the Torrent Suite Software version 3.6 using the Variant Caller plugin (germline low stringency) and subsequently the specific json file. Straightforward default analysis only identified six out of 11 BRCA1 and BRCA2 pathogenic variants, four indel and two SNV (Fig. 2; Fig. S1). Setting the kmer length parameter value to 11 instead of the default value 19, allowed detection of the otherwise missed BRCA1 indel c.4964_4982del19 (p.Ser1655Tyrfs) localized at the end of exon 16 (Fig. 3A), without generating FN.

To detect the BRCA2 complex variant c.7921_7926delGAATTTinsAG (p.Glu2641Argfs) the generation of “complex variant candidates parameter” in the Torrent Suite Software has been activated (Fig. 3B). In this way, the NGS analysis correctly identified a total of 179 variants (68 in BRCA1 and 111 in BRCA2), reaching the maximum sensitivity of 98.3% (95.3% for BRCA1 and 100% for BRCA2) identifying eight out 11 pathogenic variants, missing only the three BRCA1 LGRs (Table 4). NGS called eight FPs variants (seven in BRCA1 and one in BRCA2) in the 11 samples of TS in the homopolymeric regions of the genes. The positive predictive value, calculated as (TPs/TPs+ FPs), was 95.5% (89.7% for BRCA1 and 99.1% for BRCA2). The panel covers 16,246 nucleotides (5,989 for BRCA1 and 10,257 for BRCA2) for patient, then we covered 178,706 nucleotides (65,879 for BRCA1 and 112,827 for BRCA2) in the TS. If we consider that TNs could be calculated as the total analyzed nucleotides-TPs-FPs-FNs, we could calculate a total of 178,524 TNs, that is 65,808 for BRCA1 and 112,716 for BRCA2. The specificity, calculated as TN/(TN+FP), and the NPV, calculated as TN/(TN+FN), were 100% (Table 4). Therefore, considering LGRs as FNs, the concordance with Sanger sequencing analysis was 100%.

All these data indicated that an appropriate setting of the analytical pipeline can affect the accuracy of a variant call.

This optimized pipeline was then blindly applied on the VS of 136 samples, whose genotype was unknown. We detected pathogenic variants in 30 (22%) cases, 20 in BRCA1 (15%), 10 in BRCA2 (7%), 1 BRCA1 VUS, and 1 BRCA2 VUS. All variants were confirmed by Sanger sequencing in their correct zygosity status (Table 5). The remain 104 samples, negative for pathogenic variants, subjected to MLPA analysis were also negative for the presence of LGRs. Of note, one of the BRCA2 deleterious mutations (c.72delA) was a novel disease-causing variant at exon 3, predicted to code for an early truncated protein (p.Gly25Aspfs), thus falling into class V (Plon et al., 2008; ENIGMA consortium: https://enigmaconsortium.org/library/general-documents/enigma-classification-criteria/). The novel mutation has been identified in one HBC family having a clear dominant inheritance pattern and with three cases of breast cancer before 40 years of age (Fig. S2). It has not yet been possible to perform segregation analysis in this family.

A complete BRCA1 and BRCA2 Sanger sequencing has been performed in 60 negative samples of VS (included in the 104 negative samples), because their higher expectation to be carrier of pathogenic mutations. As expected, Sanger sequencing not only confirmed the absence of BRCA1 and BRCA2 pathogenic variants, but also identified all benign variants in their correct zygosity status as detected by NGS (Table S4). The variants located in regions not covered by the respective primers have been excluded from this evaluation. Altogether, our data strongly suggest that Ion Torrent PGM and Ion AmpliSeq BRCA1 and BRCA2 Panel (Life Technologies, Carlsbad, CA, USA) represent an excellent system to be applied in the field of diagnostics with 100% of concordance and sensitivity compared to Sanger method.

Our study enabled us to generate a standardized diagnostic workflow useful in clinical practice to quickly screen the BRCA1 and BRCA2 genes (Fig. 4).

Cost efficiency and turn-around time

In comparison to traditional Sanger-based sequencing, the NGS workflow led to a reduction in cost, analysis time, and TAT for the screening of BRCA1 and BRCA2 genes. Considering the routine of our laboratory, which include PCR and Sanger sequencing, we calculated that we analyzed eight samples in 30–40 days. While, with the advent of NGS methodology we can analyze eight patients in 10 days, making the results available for the clinician in 15 working days following blood sampling. Costs entailed are obviously site-dependent, we calculated a cost of about 1,500–3,000 euro/sample for BRCA1 and BRCA2 Sanger-based sequencing, and a consistent reduction to about 300 euro/sample with NGS.