Sample-Index Misassignment Impacts Tumour Exome Sequencing

Sample pooling enabled by dedicated indexes is a common strategy for cost-effective and robust high-throughput sequencing. Index misassignment leading to mutual contamination between pooled samples has however been described as a general problem of the latest Illumina sequencing instruments utilizing exclusion amplification. Using real-life data from multiple tumour sequencing projects, we demonstrate that index misassignment can induce artefactual variant calls closely resembling true, high-quality somatic variants. These artefactual calls potentially impact cancer applications utilizing low allelic frequencies, such as in clonal analysis of tumours. We discuss the available countermeasures with an emphasis on improved library indexing methods, and provide software that can assist in the identification of variants that may be consequences of index misassignment.

Identification of somatic variants by next-generation sequencing has become an important technique in cancer research by pinpointing the genomic causes of tumour phenotypes. An increasing number of examples have further shown that genomic aberrations have prognostic value and can inform rational clinical deployment of targeted cancer drugs 1,2 . As of now, next-generation sequencing technologies enable affordable assaying of variation in the entire tumour genome within days, but despite the availability of specialized software tools, somatic variant calling continues to pose challenges. Notably, the inherent complexity of tumours, as exemplified by aneuploidy, tumour heterogeneity, and sample impurity, often leads to important somatic variants only being detectable in low allelic fractions (AFs) [3][4][5][6][7] . AFs can be further affected due to technical reasons, such as the inability to accurately represent allelic ratios at genomic loci with low coverage 8 . As a consequence of the necessarily high required sensitivity, somatic variant calling is susceptible to random noise, systematic artefacts, and sample contamination 9,10 . When unnoticed, false positive variant calls can contribute to high costs incurred by follow-up analyses and experiments, and in the worst case support inadequate therapeutic strategies.
As previously described in the context of RNA-seq applications 11 , when material from multiple samples is being pooled together before sequencing, the sequencing technology itself can be a source of noticeable sample cross-contamination. Sample pooling (also called sample multiplexing) is a standard means of dividing the throughput of a sequencing instrument among multiple samples, relying on index sequences that uniquely barcode the material of each involved sample. Sample index misassignment, a phenomenon most evident on Illumina instruments utilizing exclusion amplification chemistry (ExAmp) and patterned flow cell technology (i.e., HiSeq 3000/4000/X Ten and NovaSeq), effectively leads to transfer of individual sequencing reads between samples included into a common pool 11 .
Our analysis of high-coverage tumour sequencing data shows that index misassignment is a source of false positive somatic variant calls in a form of true variation obtained from co-multiplexed samples.

Results
We investigated tumour-normal sample pairs of three different tumour types: diffuse large B-cell lymphoma, follicular lymphoma, and sarcoma. All samples were collected from cancer patients in Norway, and deep exome sequencing (median coverages: 315X for tumour samples, 146X for controls) was carried out on local Illumina instruments: three HiSeq 2000/2500 instruments utilizing bridge amplification (generating 42 diffuse large B-cell lymphoma samples and 31 sarcoma samples) and a HiSeq 4000 employing ExAmp (generating 81 follicular lymphoma samples and 77 sarcoma samples). All samples were subject to standardized library preparation, sequencing, and bioinformatics analysis, as adopted by the Norwegian Cancer Genomics Consortium (NCGC, http:// cancergenomics.no) (Methods, Supplementary Information).
In our assessment of index misassignment, we have examined the consequences for calling of somatic single nucleotide variants (SSNVs). For increased reliability, all tests and analyses were limited to SSNVs agreed upon by two independent variant callers (MuTect and Strelka) 3,12 . Sample-wise contamination estimates, our primary measure of contamination, were generated by Conpair 13 , based on the allelic composition on several thousand genomic marker sites in a given matched sample pair. Our tests on simulated data suggest that Conpair provides consistent quantifications, despite apparent progressive underestimation dependent on the number of contributing contamination sources (Supplementary Tables 1 and 2, Supplementary Information). We attribute this underestimation to decreasing contaminant variant AFs in sample pools of increasing size (Supplementary Fig. 1) paired with the fact that Conpair's contamination model is intended for two-sample mixtures only.
A schematic overview of all conducted experiments and their data dependencies is available on Fig. 1.

Contamination rates.
During the analysis of hundreds of sequenced tumour-normal matched pairs, we noted that data generated on the ExAmp instrument showed significantly higher sample-wise contamination levels in comparison to data from bridge amplification instruments (median per-sample contamination estimates: 0.839% vs. 0.187%, p-value < 0.001, Fig. 2a, Supplementary Fig. 2). To specifically assess the role of multiplexing on the ExAmp instrument, we sequenced 16 selected sample libraries both in pools and individually (i.e., one sample library per flow cell lane), observing significantly increased contamination rates in sequencing output from the multiplexed libraries (median per-sample contamination estimates: 0.644% vs. 0.0465%, p-value < 0.001, Supplementary Table 3). Removal of free adapters/primers prior to sequencing was in Illumina's report identified as a key measure for mitigating the cross-contamination rates (https://www.illumina.com/science/education/ minimizing-index-hopping.html), but for the 16 libraries included in our testing, performing a gel-based library purification step in combination with bead purification did not provide improvements over bead purification alone (median per-sample contamination estimates: 0.671% vs. 0.6115%, p-value = 0.159). In accordance with previous  Fig. 4). We classify the remaining variants as ambiguous.
In comparison to bridge amplification, ExAmp leads to higher occurrence of SSNVs that coincide with germline variation common in the Norwegian population ( Fig. 2b). At the same time, PC-AFs of these variants are significantly higher in the ExAmp datasets (median values: 0.508 vs. 0.125, p-value < 0.001), showing a much better correspondence to the suspected contamination source. Samples sequenced on the ExAmp instrument have significantly higher counts of suspected contaminant variants than samples sequenced on bridge amplification instruments (per-sample median counts: 4 vs. 0, p-value < 0.001) and show significant correspondence between the number of suspected contaminant variants and the estimated contamination (p-values < 0.001, Fig. 3).

Discussion
We have focused on exploring the general link between tumour sample cross-contamination and suspected contaminant variants of germline origin, but the effects of index misassignment require further consideration (Table 1). Firstly, cross-sample contamination by recurrent somatic variants would represent a class of potential false positives that would be more difficult to recognize, and may have stronger clinical implications. Concern for such cases of contamination gains relevance when several samples containing identical somatic variants are being co-multiplexed and thereby jointly contribute to elevated pool complement support (e.g., pooling serial patient biopsies or screening multiple samples likely to harbour identical hotspot mutations, particularly in sensitive analyses of high-purity samples). Secondly, when co-multiplexing tumour and control samples, false negative variant rates might be increased due to somatic variation contaminating the controls. Lastly, contamination rates within a given sample can be influenced by a combination of independent factors, such as the fraction of the pool that is constituted by the pool complement 14 or a sample's copy number profile -in copy number loss regions of high-purity tumour samples, we expect both the average AF and the number of contaminant variants to be increased due to the local underrepresentation of non-contaminant reads.
Several countermeasures have been suggested for the prevention of index misassignment, with sample storage conditions being shown by Illumina to influence the problem's severity. Sequencing of a single sample per lane proved to be an effective solution, which however might be too costly for many applications. Rigorous gel-or bead-based library purification has been strongly recommended as a means of index misassignment mitigation in case of sample multiplexing. However, in our experience, both purification methods appear to be insufficient even in combination. Dual indexing has been suggested for circumventing the problem by using sample-specific index pairs rather than individual indexes, enabling recognition of reads with unexpected index combinations caused by index misassignment 14 . Dual indexes are becoming a part of best practices for multiplexed libraries sequenced with ExAmp chemistry, but their current availability may vary depending on the applied library preparation protocol. The possible solutions therefore need to be assessed in the context of each particular project.
For data that have already been generated, suitable variant post-processing should be chosen based on the type of contamination artefacts relevant for given project (Table 1). Commonly used allelic fraction thresholds and germline variant database filters are likely to remove the majority of false positive somatic variant calls caused by index misassignment. However, a more discriminative filtering approach would be preferable in settings where sensitive detection is the priority, such as in clonal analysis. Low-AF variation can, in multiple cancer types, harbour markers of prognostic and therapeutic value [15][16][17] . When available, pool complement information can help identifying variants that are unlikely to be contamination artefacts, thereby reducing the number of potentially important true positive somatic variants discarded due to low AF or germline database presence. The code used for PC-AF calculation in our analyses is available for use in other projects (Methods). If control samples have been contaminated by true somatic variation, it might be necessary to adopt variant caller settings more permissive to somatic allele evidence in the matched normal material.
We believe our findings to be of relevance to other cancer sequencing projects that utilize the ExAmp chemistry, even though the impacts of index misassignment on somatic variant calling may vary depending on the combination of employed sequencing instruments, library preparation protocols, and bioinformatics analyses. In general, we expect the effects of index misassignment to dampen as the sequencing depth decreases and the allelic fraction threshold for accepted somatic variants increases. On the other hand, we note that besides single nucleotide variants, other types of somatic variation (e.g. insertions and deletions) are likely to be affected.

Methods
Methods used for sample collection, DNA extraction, library preparation and sequencing, as well as full details of all bioinformatics analyses, are available in the methods section of Supplementary Information.

Bioinformatics analyses.
All sequenced samples were pre-processed with BWA MEM 18 , Picard (http:// broadinstitute.github.io/picard/) and GATK 19 tools, before performing somatic calling with MuTect and Strelka (only consensus calls were considered in subsequent processing and testing). All analyses were performed using human reference genome build b37 with an added decoy contig and corresponding variation databases.
The pool complement of each individual tumour sample was formed by reads of all co-multiplexed samples originating from other individuals. For each identified somatic variant in a given tumour sample, all overlapping reads from the corresponding pool complement were extracted, and the pool complement allelic fraction (PC-AF, the fraction of reads supporting the variant within the pool complement) was calculated. Tumour Accuracy of Conpair contamination estimates was assessed with the help of artificially created exome-wide sample admixtures. Reads from multiple libraries were mixed to simulate contamination by 1, 2, 4 or 7 samples from as many different individuals, with the total contamination amounting to either ~2% or ~8% in each case, depending on the experiment. Conpair's contamination estimates were compared to true contamination content, which was validated by library-specific depth of coverage calculations.
NCGC's cohort of normal samples consisted of blood samples of 789 different individuals living in Norway. The samples were processed according to "Best Practices for Germline SNP & Indel Discovery in Whole Genome and Exome Sequence" developed by the Broad Institute, utilizing allele-specific (rather than site-specific) variant calling.  Table 1. Overview of possible contamination types, their consequences and suitable filtering options. PC-AF: pool-complement allelic fraction. *Copy number loss regions of high-purity tumour samples will be especially affected. **The enrichment will increase together with given variant's recurrence, as well as with purity of tumour samples that carry the variant.