A method for whole genome sequencing of picogram amounts of DNA

Abstract

The present invention relates to a method for whole genome sequencing of a single cell or a group of cells to identify single base variants, chromosomal structural changes, or phasing information in the genome of a single cell or a group of cells.The present invention also describes a method for preparing an index DNA library for sequencing nucleic acid molecules; a method for preparing an index DNA library for whole genome sequencing of a single cell or a group of cells to identify single base variants, chromosomal structural changes, or phasing information in the genome of a single cell or a group of cells; and whole genome sequencing of a single cell or a group of cells to provide data for identifying single base variants (SNVs), chromosomal structural changes, or phasing information in the genome of a single cell or a group of cells.

Description

[Technical Field]

[0001] The present invention relates to a method for preparing an index DNA library for sequencing, such as whole genome sequencing of a single cell or a group of cells, to identify single nucleotide variants (SNVs) in the genome of a single cell or a group of cells, to identify chromosomal structural changes, or to identify phasing information. [Background technology]

[0002] Next-generation sequencing has revolutionized our understanding of the genetic evolution of human cells in health and disease. Bulk cancer genome sequencing allows the calculation of a tumor's clonal composition by inferring the presence of variants and the few cells that contain them. Furthermore, knowledge of clonal composition allows the construction of evolutionary trees that show how a particular tumor has evolved over time (1-3). Analysis of common mutations within individual clones can be used to infer the mutational processes that likely took place during tumor evolution. Understanding which mutational processes have occurred within a tumor and the mechanisms that drive them is highly desirable, as it may provide opportunities for therapeutic intervention or to predict the tumor's evolutionary path. However, limitations in sequencing depth mean that only frequently occurring mutations that occur early in tumor formation can be detected using standard bulk whole-genome sequencing (WGS) methods (Figures 1A and 4). As a result, our ability to model evolutionary events is limited to early events established during tumor evolution, rather than recent or current processes (1, 4). This limits the practical application of our understanding of mutational processes. The study of current or recent evolutionary events requires the identification of reliable mutations with very low incidence ( 1 , 4 ).

[0003] Sequencing single cells or small populations of spatially related cells offers promise for solving this problem (Figure 4). This provides a readout of mutational processes occurring in currently observed cells (Figure 1A). However, accurate sequencing of the small (picogram) amounts of DNA obtainable from single cells or spatially related cells is extremely challenging. When dealing with small amounts of DNA, unavoidable DNA damage due to oxidation or spontaneous deamination is particularly troublesome (5). These sources of damage generate a disproportionate number of artifactual C>A and C>T mutations, respectively (5-7). Identification of these artifactual mutations as variants during variant calling results in a large number of false positive (FP) variant calls. Therefore, whole-genome amplification introduces significant errors in single-nucleotide variant (SNV) calling, preventing accurate estimation of mutation load (6). Importantly, such mutations can also be due to biological processes such as the accumulation of DNA oxidative damage or the overactivity of members of the APOBEC family of deaminases (8-10).

[0004] Previously, when working with picogram quantities of DNA, it was not possible to distinguish between biologically driven C>A and C>T mutations and artifactual mutations arising during library preparation. In addition, the routine step of whole genome amplification (WGA) prior to sequencing increases the number of artifactual mutations and propagates errors caused by DNA damage (5). Several techniques have been proposed to reduce DNA damage during library preparation or to remove false-positive results during analysis (5, 11-13). However, to date, such techniques still retain thousands of false-positive mutations and therefore require thorough validation before definitive biological conclusions can be drawn (5, 11, 12). Because thorough validation is not possible in the majority of cases (5), robust methods for removing false-positive variants from WGA sequencing data are needed.

[0005] A long-fragment read (LFR) method for whole-genome sequencing and haplotyping from 10 to 20 human cells was previously published by Complete Genomics Inc. (Peters BA, et al. Accurate whole-genome sequencing and haplotyping from 10 to 20 human cells. Nature. 2012 Jul 11;487(7406):190-5. doi:10.1038 / nature11236. PubMed PMID:22785314; PubMed Central PMCID:PMC3397394). However, this method is highly complex, prone to bias from index cross-contamination, and generates a large number of false positives.

[0006] It is therefore an object of the present invention to provide an improved method for preparing DNA libraries for sequencing, SNV analysis, identification of chromosomal structural changes, or identification of phasing information. Summary of the Invention

[0007] In a first aspect, the present invention relates to a method for whole genome sequencing of a single cell or a group of cells to identify single base variants, chromosomal structural changes, or phasing information in the genome of a single cell or a group of cells, the method comprising: i) providing a multiwell array plate containing rows and columns of reaction wells; ii) providing genomic DNA of a single cell or group of cells, the genomic DNA being distributed into multiple reaction wells of a multi-well array plate, such that there is one single-stranded genomic DNA molecule of any given locus per reaction well; iii) performing whole genome amplification (WGA) of each genomic DNA molecule in each reaction well to obtain multiple copies of the genomic DNA molecule; iv) fragmenting the DNA molecules in each reaction well and ligating a pair of loop adapters at each end or tagging them with transposase delivery adapters to form adapted DNA fragments, wherein the loop adapters or transposase delivery adapters comprise a column index (Ci) sequence or a row index (Ri) sequence, and the Ci sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a column of the multiwell array plate, or each Ri sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a row of the multiwell array plate; vi) obtaining an index DNA library by performing index PCR on the adapted DNA fragments, wherein the adapted DNA fragments are amplified using forward and reverse index primers to form index PCR products, and the forward and reverse index primers introduce row index (Ri) sequences or column index (Ci) sequences at each end of the adapted DNA fragments, respectively, so that the resulting index PCR products contain both adjacent column index (Ci) sequence pairs common to each well in a column and adjacent row index (Ri) sequence pairs common to each well in a row; and vii) sequencing the index DNA library to obtain data for identifying single base variants, chromosomal structural changes, or phasing information in the genome of a single cell or group of cells.

[0008] Advantageously, the present invention provides an index DNA library (called DigiPico, short for Digital Sequencing of Picogram DNA) for obtaining high-quality, data-rich sequencing results from picogram amounts of DNA obtained from clinical samples using the single DNA molecule sequencing approach described herein. The present invention also provides a convenient indexing strategy to substantially eliminate cross-contamination and improve ligation efficiency. One set of indexes is first introduced into the stem-loop of a common adapter. In the ligation step, all wells in each column of the plate will receive a different indexed loop adapter or transposase delivery adapter; thus, a total of 24 different oligos will be sufficient to index all columns of the plate with the first set of indexes (column indexing). After the ligation step, all wells in each row can be pooled into separate tubes, resulting in 16 different pools. These 16 different pools can typically be purified for use in the subsequent indexing step. In the next step, the purified products from each pool can be indexed (row indexed) by PCR using only 16 different index primers. This allows for single-cell sequencing with unprecedented accuracy and represents a major technical improvement over known methods. Using the present invention, private mutations and potential neoantigens in single cells or very few cells can be identified, and these mutations or neoantigens can be used as therapeutic targets. The present invention can also be used to identify chromosomal structural changes, such as numerical or structural abnormalities, or to identify phasing information. Figure 18 herein clearly demonstrates that the method of the present invention can significantly improve the accuracy of identifying actual nucleotide variants, reducing the number of false positives and enabling more accurate discrimination between samples displaying different numbers of mutations compared to the LFR method from Complete Genomics Inc.

[0009] Cells and cell groups The cell or cells may comprise eukaryotic cells, such as mammalian cells. In one embodiment, the cells are human. In one embodiment, the cells are at least diploid cells. The cells may be cancerous or precancerous cells. The cells may comprise tumor islets. In one embodiment, the cells may be derived from a tissue biopsy of a subject.

[0010] Cells such as tumor islets can be laser capture microdissected cells. When determining the DNA sequence from multiple cells, the cells can be spatially related cells.The cells can be located in the same place in the tumor or in the region of the tumor.In another embodiment, the cells can be located close to each other or not. The SNVs to be identified may comprise single base mutations. The method may be used to identify multiple different SNVs in genomic DNA.

[0011] Preparation of nucleic acid molecules and well distribution The nucleic acid can be purified or partially purified. In another embodiment, the nucleic acid can be provided in a cell lysate. Genomic DNA can be provided as purified DNA. In another embodiment, the genomic DNA can be provided from resuspended nuclei or whole cells, such as laser capture microdissected cells.

[0012] Genomic DNA can comprise the DNA of a single cell or a group of cells (cell population), such as spatially related cells. Genomic DNA can comprise the DNA of about 1-30 cells. Genomic DNA can comprise the DNA of about 1-100 cells. In another embodiment, genomic DNA can comprise the DNA of about 1-80 cells. In another embodiment, genomic DNA can comprise the DNA of about 1-50 cells. In another embodiment, genomic DNA can comprise the DNA of about 1-40 cells. In another embodiment, genomic DNA can comprise the DNA of about 10-30 cells. In another embodiment, genomic DNA can comprise the DNA of about 20-30 cells. In another embodiment, genomic DNA can comprise the DNA of about 20-40 cells. In another embodiment, genomic DNA can comprise the DNA of about 10-40 cells.

[0013] If the nucleic acid, such as DNA, is double-stranded, the nucleic acid may be denatured before distribution to the wells. Denaturation may be achieved by heat and / or a denaturing buffer. In one embodiment, nucleic acids, such as genomic DNA, or nuclei or cells containing genomic DNA, may be denatured using a denaturing buffer, such as D2 buffer from the Repli-g single cell kit (Qiagen).

[0014] Nucleic acids, such as DNA, can be distributed into wells, resulting in one single-stranded genomic DNA molecule at any given locus per reaction well. Distribution of nucleic acids can be facilitated by dilution of the nucleic acid. Thus, in one embodiment, the nucleic acid solution can be diluted. Those skilled in the art can easily determine the dilution level and the volume of solution required to achieve one single-stranded genomic DNA molecule at any given locus per reaction well. Those skilled in the art can easily determine the dilution level and the volume of solution required to achieve one single-stranded genomic DNA molecule at any given locus per reaction well. For example, if the number of cells is known, a Poisson distribution can be used for the calculation.

[0015] In one embodiment, the DNA content of a single cell can be distributed among the wells of a single row or column.Therefore, a multi-well array plate can be used to analyze a plurality of different single cells, such as each row or column.At least one well can be used for adding cells and extracting DNA content.In another embodiment, the DNA content of a cell or a group of cells is distributed among both the wells of a row or column of a single multi-well array plate.

[0016] Those skilled in the art will understand that any standard multiwell array plate can be used in the methods of the present invention. Preferably, the multiwell array plate is compatible with any available PCR and / or sequencing instrument. The multiwell array plate can include a 384-well plate, such as a 24x16-well plate. In another embodiment, the multiwell array plate can include a 1536-well plate. Those skilled in the art will understand that a larger number of Ri and / or Ci sequences may be required to index a larger array plate.

[0017] The use of a 384 multiwell array plate can advantageously provide enough wells to distribute diluted genomic DNA strands of approximately 20-30 cells.

[0018] amplification In embodiments where the nucleic acid is genomic DNA, amplification of the genomic DNA molecule can include whole genome amplification (WGA). WGA can include adding amplification reagents for DNA amplification to genomic DNA. Amplification reagents for DNA amplification are alternatively referred to as an "amplification mixture." Those skilled in the art will appreciate that an amplification mixture can contain all the reagents necessary to amplify DNA (i.e., create multiple copies of DNA). These components can include a reaction buffer, polymerase, and dNTPs. For example, a DNA polymerization reporter molecule, such as a DNA-binding dye (e.g., Evagreen™), can be provided in the amplification mixture, allowing for monitoring of the amplification reaction using real-time PCR. The DNA-binding dye can be composed of two monomeric DNA-binding dyes linked by a flexible spacer. In the absence of DNA, the dimeric dye can transition via equilibrium to a random structure that can bind to DNA and fluoresce when DNA is available that can adopt a loop structure inactive to DNA binding. Amplification reagents may be added to each well before or after the DNA is added to the well.

[0019] One of skill in the art would be able to provide suitable conditions for the amplification reaction to occur, including suitable temperatures and incubation times. For example, the plate may be incubated at about 30°C for at least about 1 hour, followed by heat inactivation, e.g., at about 65°C for at least 5 minutes.

[0020] Fragmentation and ligation of loop adapters or transposase delivery adapters In one embodiment, loop adapters are provided such that the method includes the steps of fragmenting DNA molecules in each reaction well and then ligating the loop adapters to the fragmented DNA. The fragmented DNA may be end-repaired prior to ligation. In an alternative embodiment, transposase delivery adapters are provided such that the method includes fragmenting DNA molecules with a tagging step. Tagging can include providing a transposase, such as Tn5, with an oligonucleotide, referred to herein as a transposase delivery adapter. Those skilled in the art will be familiar with routine techniques and reagents for performing tagging and forming adapted DNA molecules.

[0021] Fragmenting the DNA molecules in each reaction well into multiple dsDNA fragments can include direct fragmentation, such as enzymatic or mechanical fragmentation. In one embodiment, fragmenting the DNA includes enzymatic fragmentation.

[0022] DNA fragmentation or tagging can be achieved by adding a fragmentation or tagging reagent to the DNA in each well. The fragmentation or tagging reagent can be added to each well simultaneously, for example, by using a multi-well dispenser such as an I-DOT (Dispendix, Germany) dispenser or similar. The fragmentation or tagging reaction is carried out for a set time to obtain fragments of the desired size. Those skilled in the art will understand that the timing of the fragmentation or tagging reaction can depend on the method used, such as the type and level of enzyme used in the reaction. Therefore, those skilled in the art can follow standard protocol timing for a given reaction, such as that of a reaction kit.

[0023] The fragmentation reagent may include a restriction enzyme or a nicking enzyme, such as DNase I. When a nicking enzyme is provided, a single-strand-specific enzyme is provided that recognizes the nick site and then cleaves the second strand. In one embodiment, a library preparation kit, such as the Lotus DNA Library Preparation Kit (IDT, USA), may be used.

[0024] After DNA fragmentation to form dsDNA fragments, the dsDNA fragments can be end-repaired and dA-tailed, allowing them to be ligated to other DNA molecules, such as loop adapters. Enzymes for end-repair and / or dA-tailing can include a DNA polymerase, such as T4 DNA polymerase, and a polynucleotide kinase (PNK), such as T4 polynucleotide kinase. T4 DNA polymerase (in the presence of dNTPs) can fill in the 5' overhang and trim the 3' overhang to the dsDNA boundary to generate a blunt end. T4 PNK can then phosphorylate the 5'-terminal nucleotide. A DNA polymerase, such as Taq DNA polymerase, which has terminal transferase activity and leaves a 3'-terminal adenine, can be used for A-tailing.

[0025] In one embodiment, dsDNA fragmentation, end repair, and dA-tail addition are all performed in a single reaction.

[0026] In one embodiment, loop adaptors can be introduced onto the fragmented DNA by ligation. Ligation of loop adaptors to dsDNA fragments can include the addition of loop adaptors and a ligase such as T4 DNA ligase.

[0027] The loop adaptor may comprise an oligonucleotide, such as DNA, having a secondary stem-loop structure. The stem-loop structure may be provided by a single oligonucleotide molecule comprising a pair of complementary sequence regions adjacent to the loop region, the complementary sequence pair being positioned to hybridize with each other to form the stem-loop structure of the loop adaptor. The loop adaptor may further encode a column index (Ci) sequence or a row index (Ri) sequence in the stem region.

[0028] The column index (Ci) sequence or row index (Ri) sequence can comprise a predetermined sequence that can label DNA by row or column, respectively. The column index (Ci) sequence or row index (Ri) sequence can be at least 3 nucleotides in length.

[0029] In one embodiment, the ends of the adaptor-attached DNA fragments can be symmetrical. In particular, the loop adaptors or transposase delivery adaptors ligated to each end of the dsDNA fragments are the same, so that each dsDNA fragment contains a pair of identical adjacent loop adaptors or transposase delivery adaptors. The Ci sequence pairs on the same adaptor-attached DNA fragments can be the same. Alternatively, if an Ri sequence is provided, the Ri sequence pairs on the same adaptor-attached DNA fragments can be the same.

[0030] Advantageously, the provision of two identical Ci or Ri sequences on the adapted DNA fragments provides a marker that avoids analysis of index DNA library sequences that are the result of cross-contamination between different columns or different rows, respectively. In particular, index DNA library sequences that do not have a matching Ci sequence at each end can be discarded from data analysis. Alternatively, if Ri sequences are provided, any index DNA library sequences that do not have a matching Ri sequence at each end can be discarded from data analysis. This can provide a first level of redundancy to eliminate index cross-contamination from subsequent data, which is a serious problem in the preparation of index libraries and their subsequent analysis.

[0031] The loop adapter can provide a 3' or 5' overhang to aid in ligation to the dsDNA fragment. When the stem regions of the loop adapter hybridize together (i.e., the loop adapter is in a secondary / stem-loop structure), a 3' or 5' overhang can be provided. The 3' or 5' overhang can correspond to a complementary overhang on the dsDNA fragment that has been end-repaired and prepared for ligation. The overhang can include a single thymine.

[0032] The loop adaptor sequence may comprise SEQ ID NO: 1, or a functional variant thereof.

[0033] After ligation of the loop adaptor, the single-stranded region of the loop DNA can be cleaved. The single-stranded region of the loop DNA can be enzymatically cleaved, for example, by a USER (Uracil-Specific Excision Reagent) enzyme, which generates a single-base gap at the position of the uracil present in the loop. Thus, in one embodiment, the loop adaptor can include a uracil in the loop region.

[0034] Pooling rows of wells If Ci sequences are provided in the adapted DNA fragments, the method may additionally include the step of pooling the adapted DNA fragments from each reaction well in a row prior to index PCR. Alternatively, if Ri sequences are provided in the adapted DNA fragments, the method may additionally include the step of pooling the adapted DNA fragments from each reaction well in a column prior to index PCR. The pooled adapted DNA fragments may then be used for index PCR for each row, or a single pooled reaction for each column, depending on which is pooled. In alternative embodiments, columns or rows may not be pooled prior to performing index PCR.

[0035] Advantageously, row or column pooling prior to index PCR greatly improves the efficiency of library preparation. For example, for a 16x24 (384) well plate, if 16 rows are pooled during the introduction of loop adapters or transposase delivery adapters and during the index PCR step, only 16 index PCR reactions are required, rather than 384 index PCR reactions if they are not pooled.

[0036] Size selection and index PCR Prior to index PCR, the size of the adapted DNA fragments can be selected; for example, if applicable, self-ligated adapters can also be removed from the reaction. An example of a desired size is approximately 300-400 bp in length. Size selection can be achieved by isolating or purifying adapted DNA fragments of the desired length using, for example, a gel or beads. Solid Phase Reversible Immobilization beads (SPRI beads) can be used for size selection. SPRI beads can contain magnetic particles coated with carboxyl groups (in the form of succinic acid) that can nonspecifically and reversibly bind to DNA.

[0037] Index PCR can include mixing an adapted DNA fragment with a set of forward and reverse index PCR primers and PCR reagents. The forward and reverse index PCR primers can include sequences configured to hybridize to sequences of the adapted DNA fragment for priming polymerization. The sequences configured to hybridize to sequences of the adapted DNA fragment for priming polymerization can be complementary sequences. The sequences for priming polymerization from the forward and reverse index PCR primers can be provided by loop adapters or transposase delivery adapters. The sequences for priming polymerization from the forward and reverse index PCR primers are adjacent to the Ci or Ri sequences of the adapted DNA fragment, thereby incorporating the Ci or Ri sequences into the index PCR product.

[0038] The complementary sequences for hybridization provided by the forward and reverse primers can each be about 15-30 nucleotides in length, for example, about 26 nucleotides in length.

[0039] In embodiments in which the adapted DNA fragments contain a Ci sequence, the forward and reverse index PCR primers may each contain an Ri sequence to obtain a pair of Ri sequences in the index PCR product. If rows are pooled, an Ri sequence is added to each adapted DNA fragment in the pool (from all wells in the row). Alternatively, if rows are not pooled, each well in the row may be provided with the same Ri sequence.

[0040] In an alternative embodiment in which the adapted DNA fragments contain an Ri sequence, the forward and reverse index PCR primers may each contain a Ci sequence to obtain a pair of Ci sequences in the index PCR product. If columns are pooled, a Ci sequence is added to each adapted DNA fragment in the pool (from all wells in the column). Alternatively, if columns are not pooled, each well in the column may be provided with the same Ci sequence.

[0041] The row index (Ri) sequences, which may be provided by the forward and reverse primers, may be the same for each adapted DNA fragment in or from a row. Alternatively, the tandem index (Ci) sequences, which may be provided by the forward and reverse primers, may be the same for each adapted DNA fragment in or from a tandem.

[0042] The row index (Ri) or row index (Ci) sequences provided by the forward and reverse primers can be at least 3 nucleotides in length, for example, about 8 nucleotides in length.

[0043] The resulting ends of the index PCR product may be symmetrical. For example, the sequences flanking the original DNA fragment sequence may be symmetrical. The index PCR product may comprise a DNA fragment sequence that is further flanked by a pair of identical Ci sequences (i.e., inner flanking) and a pair of identical Ri sequences (i.e., outer flanking). In an alternative embodiment, the index PCR product may comprise a DNA fragment sequence that is further flanked by a pair of identical Ri sequences (i.e., inner flanking) and a pair of identical Ci sequences (i.e., outer flanking).

[0044] Advantageously, providing two identical pairs of Ci or Ri sequences on the index DNA fragments, in addition to the previously provided Ri or Ci sequences (provided by the loop adapter or transposase delivery adapter), further provides a marker to avoid analysis of index DNA library sequences that are the result of cross-contamination between different columns or different rows, respectively. In particular, index DNA library sequences that do not have matching Ci sequences at each end can be discarded from data analysis. Alternatively, if Ri sequences are provided, any index DNA library sequences that do not have matching Ri sequences at each end can be discarded from data analysis. Providing both pairs of matching Ci and Ri sequences on the index DNA fragments can provide first and second levels of redundancy to eliminate index cross-contamination from subsequent data, which is a serious problem in the preparation of index libraries and their subsequent analysis.

[0045] The forward and reverse index PCR primers further comprise a sequencing adapter sequence, whereby the sequencing adapter is incorporated into the index PCR product. The sequencing adapter sequence on the primers can be 5'.

[0046] The sequencing adapters may be at the ends of the index PCR product. If sequencing adapter sequences are provided, the resulting ends of the index PCR product may not be symmetrical. For example, one end of the index PCR product may be fitted with a different sequencing adapter than the sequencing adapter at the other end. Those skilled in the art will know the sequencing adapters that may be required for a given sequencing technology. For example, in the case of dye sequencing (e.g., Illumina dye sequencing), the sequencing adapters may be P5 and P7 sequencing adapters (i.e., P5 at one end and P7 at the other end of the index PCR product). The index primer that provides the P5 sequence may comprise the sequence of SEQ ID NO:2. The index primer that provides the P7 sequence may comprise the sequence of SEQ ID NO:3.

[0047] Once formed, the index PCR products may be referred to as "index DNA library sequences" or "index DNA fragments." The pool index PCR products, index DNA sequences, or index DNA fragments may be referred to as an "index DNA library."

[0048] Index DNA Library The index DNA fragment size of the index DNA library is filtered so that only index DNA fragments of a desired or suitable length are available for sequencing. After index PCR, the index PCR fragments can be purified / isolated, for example, with beads (e.g., SPRI beads). Purification can remove undesired short fragments, primer dimers, or other PCR artifacts or reagents.

[0049] The index library can be examined for appropriate size distribution; we typically examine library size distribution using, for example, a tapestation or bioanalyzer instrument (Agilent), or similar. After preparation, the size of the index DNA library can be adjusted by dilution to, for example, about 4 nM for sequencing.

[0050] The index DNA library can be stored for subsequent use, such as sequencing. For example, the DNA library can be stored frozen or refrigerated.

[0051] Sequencing of index DNA libraries The index DNA library can be sequenced or configured to be sequenced. The sequencing can be next-generation sequencing (NGS). The sequencing can be dye sequencing (e.g., Illumina dye sequencing), nanopore sequencing, or ion torrent sequencing. Those skilled in the art will be familiar with the many different sequencing techniques / methods that can be used and the sequencing adapters required therefor.

[0052] The sequencing can be multiplex sequencing, in which multiple index DNA libraries are sequenced simultaneously.

[0053] Mutation / nucleotide change identification and data analysis The method may include identifying any actual SNVs in the genome of a single cell or a group of cells by determining whether substantially all index DNA library sequences from a single well contain the same SNV, or whether only some of the index DNA library sequences contain the same SNV. SNVs that appear in substantially all index DNA library sequences from a single well may be identified as actual SNVs in the genomic DNA. Additionally or alternatively, SNVs found in only some of the index DNA library sequences from a single well may be identified as false positive (FP) SNVs. False positive SNVs may be damage-induced errors or replication errors.

[0054] The method may further include matching an index DNA library sequence from a single well representing one strand of genomic DNA with an index DNA library sequence from another well representing the complementary strand of genomic DNA. SNVs that are substantially present in all index DNA library sequences of both complementary strands of genomic DNA may be identified as actual SNVs. SNVs that are substantially absent in all index DNA library sequences of both complementary strands of genomic DNA may be identified as false positives (i.e., not actual SNVs).

[0055] The step of determining whether substantially all index DNA library sequences from a single well contain the same SNV, or whether only some index DNA library sequences contain the same SNV, may be performed in silico, e.g., using BAM file data. Additionally or alternatively, the step of matching index DNA library sequences from a single well representing one strand of genomic DNA with index DNA library sequences from another well representing the complementary strand of genomic DNA may be performed in silico, e.g., using BAM file data.

[0056] In one embodiment, sequencing data from tumor cells, suspected tumor cells, or precancerous cells may be compared with sequencing data obtained from normal cells (i.e., non-cancerous cells) taken from normal tissue (i.e., non-cancerous tissue). Thus, in one embodiment, the method includes preparing index DNA libraries from tumor cells, suspected tumor cells, or precancerous cells, and normal (i.e., non-cancerous) cells. Index DNA libraries may be prepared for each cell type simultaneously, e.g., in different wells of the same multi-well plate, or separately. Sequencing of index DNA libraries from different cell types may be performed during the same sequencing experiment. The different cell types (e.g., cancerous or normal) may be from the same subject.

[0057] A probability score that a particular nucleotide variant is an actual SNV or a false positive can be calculated in silico, thereby identifying a given variant nucleotide as having a statistically significant probability of being an actual SNV or a false positive.

[0058] In one embodiment, sequencing the DNA library and identifying SNVs within the library includes generating multiplexed sequencing data from multiple wells and analyzing the data for SNVs.

[0059] In one embodiment, analyzing the SNV data includes demultiplexing the sequencing data, whereby data from each well is assigned to a separate group of wells. Additionally, separate index DNA libraries may be sequenced in the same sequencing experiment, whereby the method may further include demultiplexing the sequencing data to identify / group different index DNA libraries.

[0060] The provided sequence data can be in the form of a paired-read FastQ file. The sequence data, for example, the paired-read FastQ file, can be trimmed to remove adapter sequences. The sequence data, for example, the paired-read FastQ file, can also be trimmed for quality. Those skilled in the art can easily adjust the desired threshold for the quality score of each base read in the sequence, for example, using a program such as TrimGalore. The resulting data is called "trimmed data."

[0061] In one embodiment, data analysis for SNVs involves mapping the sequence data to a reference genome, such as the human hg19 reference genome, to generate aligned sequencing data in sequence alignment map (SAM) format or its binary file version format (e.g., BAM file). Mapping to the reference genome can use trimmed read data. The SAM or BAM file can be used to identify SNVs present in each well.

[0062] Mapping can be performed using a program such as Bowtie2, for example using Picard Tools, with the ignore-quals parameter enabled and duplicate reads marked.

[0063] For example, a variant caller program such as Platypus variant caller can be used to perform joint variant calling on all individual BAM files together with the merged BAM file from all wells.

[0064] Low-quality (i.e., low-confidence) variants can be filtered out of the data. For example, low-quality (i.e., low-confidence) variants can be removed from the data by applying a quality filter. Examples of quality filters in the Platypus caller include QUAL>60, FR>0.1, HP≦4, QD>10, and SbPval≦0.95. Those skilled in the art will understand that filtering out low-confidence variants is a routine procedure, and that each variant caller, depending on its algorithm, has different confidence scores for each variant that can be used to filter out low-confidence (low-quality) variants. Thus, the specific parameters can depend on the variant caller used.

[0065] The total number of wells covering each locus (Tw) and the number of wells carrying each variant (Vw) can be determined. Well count filters, e.g., Tw>5, Vw>2, and Vw / Tw>0.1, can be applied to retain only high-confidence positions for analysis.

[0066] Regions of the genome with poor mappability (i.e., known regions with a high probability of misaligning reads) can be removed from the analysis using, for example, VCFtools.

[0067] The resulting list of high-confidence variants identified in the data can then be used to perform variant recall (genotyping) on ​​WGS data from, for example, blood and bulk tumors, using, for example, Platypus. The minPosterior parameter in Platypus can be set to 0 and the minMapQual parameter to 5. Any variants not reliably retained in both standard WGS data can be extracted as UTD (Unique to DigiPico) variants. Any variants that are also reliably present in the bulk sequencing data of the blood sample (based on GATK analysis) can be extracted as TP (True Positive) variants.

[0068] Use of Artificial Neural Networks (ANN) In silico identification or matching of index DNA sequences and / or in silico calculation of probability scores may be performed according to the methods and calculations described herein. In one embodiment, in silico identification or matching of index DNA sequences and / or in silico calculation of probability scores may be performed by an artificial neural network (ANN) model, such as by a multilayer perceptron.

[0069] The multilayer perceptron may have an input layer consisting of N neurons (e.g., N=41), where N is the number of features used in each experiment. The ANN model may include at least two hidden layers with ReLU (rectified linear unit) activation. The final layer of the ANN may be a single output neuron with sigmoid activation. The loss function may be binary cross-entropy.

[0070] The ANN may be programmed, for example, using Keras in Python 3. Those skilled in the art will recognize that Keras is an open-source Python library for developing and evaluating deep learning models, although other libraries may be used. The ANN can be pre-trained on one or more datasets, for example, the ANN can be trained on a dataset containing known nucleotide variants.

[0071] Other Aspects In another aspect of the present invention, there is provided a method for preparing an index DNA library for performing sequencing of a nucleic acid molecule, the method comprising: i) providing a multiwell array plate containing rows and columns of reaction wells; ii) providing nucleic acid molecules, the nucleic acid molecules being distributed into a plurality of reaction wells of a multi-well array plate, such that there is one single-stranded nucleic acid molecule from any given locus per reaction well; iii) performing amplification of the nucleic acid molecule in each reaction well to obtain multiple copies of the nucleic acid molecule; iv) fragmenting the DNA molecules in each reaction well and ligating a pair of loop adapters at each end or tagging them with transposase delivery adapters to form adapted DNA fragments, wherein the loop adapters or transposase delivery adapters comprise a column index (Ci) sequence or a row index (Ri) sequence, and the Ci sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a column of the multiwell array plate, or each Ri sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a row of the multiwell array plate; vi) obtaining an index DNA library by performing index PCR on the adapted DNA fragments, wherein the adapted DNA fragments are amplified using forward and reverse index primers to form index PCR products, and the forward and reverse index primers introduce row index (Ri) sequences or column index (Ci) sequences at each end of the adapted DNA fragments, respectively, so that the resulting index PCR products contain both adjacent column index (Ci) sequence pairs common to each well in a column and adjacent row index (Ri) sequence pairs common to each well in a row; and Optionally, the forward and reverse index primers further provide 5' and 3' sequencing adapters on the index PCR products suitable for use in a sequencing reaction.

[0072] The nucleic acid may be DNA or RNA. In one embodiment, the nucleic acid is genomic DNA. In another embodiment, the nucleic acid may be mRNA.

[0073] In another aspect of the present invention, there is provided a method for preparing an index DNA library for whole genome sequencing of a single cell or a group of cells to identify single base variants, chromosomal structural changes, or phasing information in the genome of a single cell or a group of cells, the method comprising: i) providing a multiwell array plate containing rows and columns of reaction wells; ii) providing genomic DNA of a single cell or group of cells, the genomic DNA being distributed into multiple reaction wells of a multi-well array plate, such that there is one single-stranded genomic DNA molecule of any given locus per reaction well; iii) performing whole genome amplification (WGA) of each genomic DNA molecule in each reaction well to obtain multiple copies of the genomic DNA molecule; iv) fragmenting the DNA molecules in each reaction well and ligating a pair of loop adapters at each end or tagging them with transposase delivery adapters to form adapted DNA fragments, wherein the loop adapters or transposase delivery adapters comprise a column index (Ci) sequence or a row index (Ri) sequence, and the Ci sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a column of the multiwell array plate, or each Ri sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a row of the multiwell array plate; vi) obtaining an index DNA library by performing index PCR on the adapted DNA fragments, wherein the adapted DNA fragments are amplified using forward and reverse index primers to form index PCR products, and the forward and reverse index primers introduce row index (Ri) sequences or column index (Ci) sequences at each end of the adapted DNA fragments, respectively, so that the resulting index PCR products contain both adjacent column index (Ci) sequence pairs common to each well in a column and adjacent row index (Ri) sequence pairs common to each well in a row; and Optionally, the forward and reverse index primers further provide 5' and 3' sequencing adapters on the index PCR products suitable for use in a sequencing reaction. The indexed nucleic acids can be sequenced, for example, as described herein.

[0074] Thus, in another aspect of the present invention, there is provided a method for whole genome sequencing of a single cell or group of cells to provide data for the identification of single nucleotide variants (SNVs) in the genome of a single cell or group of cells, the method comprising: i) preparing an index DNA library by carrying out the method according to the invention herein, or providing an index DNA library prepared by the method according to the invention herein; ii) sequencing the index DNA library to obtain data to identify any single nucleotide variants (SNVs) in the genome of the single cell or group of cells.

[0075] Sequencing data can be used to identify SNV, for example, as described herein. Additionally or alternatively, sequencing data can be used to identify genetic changes associated with chromosomal structural changes. Chromosomal abnormalities can include numerical and / or structural abnormalities.

[0076] Additionally or alternatively, sequencing data may be used to identify phasing information in a cell or group of cells. The present invention may also include one or more of the features, alone or in combination, as described and / or disclosed in the drawings.

[0077] definition The term "spatially associated cells" is understood to mean cells that are immediately next to each other.

[0078] The term "false positive (FP) mutation" or "false positive (FP) SNV" is understood to mean a variant nucleotide that was not present in the genome prior to DNA extraction from an untreated cell; for example, the false mutation may be a damage-induced error or a replication error.

[0079] The terms "actual mutation / SNV" or "true positive mutation / SNV" are used interchangeably and are understood to mean a variant nucleotide present in the genomic DNA of a living cell prior to DNA extraction.

[0080] The term "single nucleotide variant" (SNV) can include a single nucleotide polymorphism (SNP) or any other change, such as a mutation in a sequence. A mutation or change can include a nucleotide substitution, addition, or deletion in a given sequence.

[0081] "Chromosomal abnormality" is understood to be a missing, redundant, or irregular portion of chromosomal DNA. It can be due to the typical number of chromosomes or structural abnormalities in one or more chromosomes. They include various abnormalities, such as deletion, duplication, and insertion. Balanced abnormalities such as inversion and interchromosomal and intrachromosomal translocation may occur. In addition, mobile element insertion, segment duplication, and multi-allelic chromosomal numerical abnormality may occur. The final combination of multiple above may result in complex rearrangement.

[0082] "Phasing" is understood to be the act or process of assigning alleles (As, Cs, Ts, and Gs) to paternal and maternal chromosomes. Phasing can help determine whether a match is on the paternal or maternal side, both sides, or neither side. Phasing can also aid in the process of assigning chromosomal mapping segments to specific ancestors.

[0083] Those skilled in the art will appreciate that any feature of one embodiment or aspect of the present invention may be applied to other embodiments or aspects of the present invention, where appropriate. Embodiments of the invention are hereinafter described in more detail, by way of example only, with reference to the accompanying drawings, in which: [Brief explanation of the drawings]

[0084] [Figure 1]Rationale, workflow, and performance of DigiPico sequencing. (A) WGS techniques can only identify early mutational processes (EM) in the dominant proliferating clone in a tumor (red and blue). Currently active mutational processes (CM) generate a diverse set of subclones with various clone-specific mutations. This diversity determines the evolutionary path of the tumor. (B) Template partitioning prior to WGA, such that each compartment contains one DNA molecule from each locus, allows for the identification of artifactual mutations. Because damage-induced errors (red) and replication errors (cyan) occur stochastically during replication, artifactual mutations result in biallelic compartments. Note that actual mutations are always present in all product DNA strands within the same compartment. (C) DigiPico sequencing workflow. LCM: laser capture microdissection. (D) Endpoint relative fluorescence units (RFU) from Evagreen-labeled DNA were used to ensure uniform distribution of template and WGA processes across the plate. RFU values ​​were normalized to a median of 1 for each experiment. (E) Per-well qPCR using Illumina adapter primers (P5 and P7) measures the relative amount of adapter ligation product in each well. Ct values ​​were normalized to a median of 0 for each experiment. (F) The need to streamline the DigiPico library preparation method led to the miniaturization of WGA, which can specifically and sensitively amplify sub-picogram amounts of DNA in all wells. Values ​​represent the average RFU of nine replicates. Error bars indicate SD. (G, H, and I) Preliminary analysis of DigiPico sequencing data from individual wells in each experiment demonstrates high-quality and uniform sequencing with high mapping rate, depth, and breadth of coverage, as shown. (J) Definition of DigiPico-unique (UTD) variants. Subtraction of identifiable SNVs in standard WGS data from the corresponding DigiPico data results in UTD variants. These primarily consist of artifactual mutations as well as some clone-specific mutations.Because the templates in experiment D1110 are a subset of those actually used in standard WGS, all true variants in DigiPico experiment D1110 are expected to be similarly represented in the standard WGS data. In contrast, clone-specific variants in experiment D1111 may be absent from the standard WGS data due to depth constraints, even though DNA molecules carrying such variants may be present at extremely low frequencies in the bulk DNA sample. In all box plots, horizontal lines indicate medians. Boxes indicate interquartile ranges (25th to 75th percentiles). Whiskers indicate ranges excluding outliers. Outliers are defined as data points 1.5 times higher or lower than the interquartile range. [Figure 2]MutLX algorithm, design, and results. (A) Comparison of the number of wells carrying various variants in experiment D1110 confirms that, as hypothesized, the majority of UTDs are present in only a few wells. The horizontal line is the median. The box indicates the interquartile range. The whiskers indicate the range excluding outliers, which are defined as 1.5x outside the interquartile range. (B) Similarly, the proportion of dual-allelic compartments for UTDs appears significantly higher compared to true variants. This value was calculated by dividing the number of wells containing coexisting variant and reference alleles by the total number of wells with evidence of a variant allele. (C) Illustrates key challenges in analyzing DigiPico data using ANNs. Each circle / star represents one variant. The red line indicates the behavior of the classification model. All variants above and / or to the left of the line are predicted by the model as true variants. Analysis of samples without clone-specific variants should yield accurate separation between real and artifactual mutations. In contrast, analysis of samples containing clone-specific mutations leads to suboptimal models, which can lead to overfitting for true UTDs. This forces the model to remove all FP calls at the expense of losing nearly all clone-specific variants. (D) A diagram showing the two-stage training process in MutLX. The first training step identifies some mislabeled true mutations (gray circles) in the UTDs. All potentially mislabeled data points are temporarily removed from the second training step (black), thereby obtaining a better model for assigning probability scores to all mutations. Finally, combining the probability scores obtained from the model with uncertainty estimates for these probability scores (as described in E) allows for efficient removal of FP calls while maintaining excellent sensitivity for true clone-specific variants. (E) A diagram showing the test-time dropout analysis for calculating the uncertainty estimates for the probability scores. Black neurons indicate neurons that were stopped during the dropout analysis.Accepting only variants with high probability scores and low uncertainty scores should allow for the removal of FP variant calls. (F) ROC curves are shown for the output of MutLX analysis for experiments D1110, D1111, DE011, and GM12885. Circles represent the default cutoff values ​​determined by MutLX. (G) Bar graphs showing the number of UTDs passed in the output of SCcaller, Platypus, and MutLX. Because no true UTDs are presumed to exist in experiments D1110, DE011, and GM12885, the number of UTDs in these experiments represents the FP rate for each analysis method. Platypus values ​​are based on DigiPico-specific filtering criteria before application of MutLX. [Figure 3]Identifying active mutational processes using DigiPico / MutLX. (A) Schematic diagram of tumor evolution in HGSOC patient #11152. Standard bulk WGS of various tumor samples identified approximately 11,000 common somatic mutations across all sites. The purple dotted line indicates the time point at which the most recent common ancestor of the tested tumor samples diverged. Bulk sequencing also identified approximately 5,000, 3,000, and 2,000 subclonal mutations specific to the pre-chemotherapy omental mass, PT2R recurrence, and PALNR recurrence, respectively. However, these mutations can arise at any time during the expansion of these clones, resulting in a bias toward older mutations. This is due to limitations in identifying low-incidence somatic mutations. However, in each of these samples, DigiPico sequencing of the five pre-chemotherapy tumor islands, PT2R, and PALNR recurrence sites identified a variable number of recently emerged clonally specific mutations (represented by red numbers). The significantly larger number of clonally specific variants in PT2R indicates the presence of an active mutational process. (B) This active mutational process is highlighted by a strong clonally specific kataegis event on chromosome 17 in experiment D1111. The Y-axis represents the pairwise distance of consecutive somatic mutations on a log scale. Only mutations from chromosome 17 are shown. Mutations associated with subclonal kataegis events are highlighted in boxes; nearly all of these are in the form of strand-specific C>T or C>G mutations. This suggests the involvement of APOBEC enzymes in this hypermutation process. (C) Representative examples of some mutations involved in kataegis. The presence of all mutations on the forward strand of the genome further confirms the involvement of a hypermutation event (Figure 13). [Figure 4]Challenges in Identifying Recent Mutations. While ancient mutations can be easily investigated from bulk tumor sequencing data, investigating recent mutations from such data is hindered by the low variant allele fraction (VAF) of the involved mutations. Therefore, heuristic filtering criteria are not sufficient for identifying recent mutations. Reliable investigation of recent mutations requires the investigation of single cancer cells or tumor islands isolated by laser capture microdissection (LCM) (Figure 1). However, WGA of the limited amount of template in such samples generates a large number of false-positive variant calls, which hinders the identification of island-specific variants. Our analysis pipeline, MutLX, overcomes this issue by removing FP variant calls from DigiPico sequencing data. [Figure 5]Workflow analysis for DigiPico data. (1) Next-generation sequencing reads from normal tissues, bulk tumors, and DigiPico libraries are first mapped to the human genome, generating bam files. DigiPico reads are split into 384 FastQ files, one for each well of a 384-well plate. (2) The 384 individual bam files from DigiPico are merged into a single bam file. (3) Using the Platypus variant caller, new joint variant calling is performed on the 384 bam files and the merged bam file. The addition of the merged bam file ensures that variants with low per-well coverage are not missed during variant calling. (4) The resulting new DigiPico variants are then used as a reference for variant recall from standard WGS data of normal tissues and bulk tumors. (5) The variant recall data can then be used to extract DigiPico-specific (UTD) variants by removing all variants with reads retained in the standard WGS data. (6) Standard WGS data were also used for variant calling using GATK to obtain a list of high-confidence germline SNPs. (7) This list was then used as a guide to extract TP variant calls from the DigiPico data. For this purpose, any variants identified in the bulk blood sample using GATK and similarly in the DigiPico data using Platypus were assumed to be real. (8) The UTD and germline SNPs were then used by MutLX to train a binary classification model for the extraction of clonally specific variants from the UTD (Figure 6). [Figure 6]MutLX analysis algorithm. (1) UTD variants are identified by extracting WGS data from DigiPico data. (2) The UTDs and SNPs are used as a training set to train a primary binary classification model. (3) This model is used for primary analysis of the training set, which allows for the generation of an improved training set (4). (5) The improved training set is then used to generate a classification model (6), which can be used for the analysis of UTDs. (7) A "probability score" indicating the likelihood that the variant is a real variant, and (8) an "uncertainty score" as a measure of the uncertainty of the calculated probability score are calculated for each variant using the model. (9) A high probability score with high certainty (low uncertainty score) identifies a true variant. [Figure 7] Probability score values ​​for experiments D1110, D1111, DE011, and GM12885. A cutoff value of 0.2 removes most FP variant calls from the UTD while retaining nearly all germline SNPs in all samples. [Figure 8] Data simulations confirmed that AUC correlated negatively with the number of true UTD variants. In experiments D1110 and DE111, various numbers of somatic mutations were artificially mislabeled as UTDs (UTD*), resulting in UTD* / UTD ratios of 1%, 5%, and 10%. Both experiments were performed on 200 pg of purified DNA from bulk tumor samples, neither of which was predicted to contain true UTD variants. Independent analysis of each of these simulated datasets using MutLX confirmed the negative correlation between the number of true UTDs in the dataset and AUC. Notably, a UTD* / UTD ratio of as little as 1% (16 and 36 variants in experiments D1110 and DE11, respectively) appears to decrease AUC, suggesting that the presence of even a small proportion of true clonally specific variants can perturb the ROC curve. [Figure 9]Analysis of synthetic DigiPico datasets. In experiments (A) D1110 and (B) DE111, various numbers of high-confidence somatic mutations were artificially mislabeled as UTDs (UTD*) to synthetically increase the number of true UTDs at various UTD* / UTD ratios. Both experiments were performed on 200 pg of purified DNA from bulk tumor samples, neither of which was predicted to harbor true UTD variants. The results show that the presence of true UTD variants among the majority of artifactual UTDs does not appear to compromise the robustness of the classification model generated by MutLX. Each boxplot is derived from 10 different UTD* subsets used in the analysis. To obtain comparable FP rates across experiments, analytical cutoff values ​​were used to achieve 90% TPR and 95% TPR in germline variants for all synthetic datasets in experiments D1110 (A) and DE111 (B), respectively. Boxplots show median, interquartile range, and outlier exclusion range. Outliers are defined as 1.5 times higher or lower than the interquartile range. [Figure 10] Targeted sequencing of several clone-specific variants identified in experiment D1111. Amplicon sequencing of the target sites was performed on the MiSeq platform. Three of the 14 targets (highlighted in orange) appeared to have high noise levels in the blood sample and were therefore considered inconclusive. Only one of the remaining 11 mutations did not appear to have any evidence in the PT2R tumor bulk DNA sample (highlighted in blue). VAWF: Variant allele well frequency. [Figure 11] Targeted sequencing of several artifactual variants identified in experiment DE111. [Figure 12]Frequency of various variants in FP calls identified by MutLX in DigiPico data. Green dots (left) are obtained from the analysis of somatic variants in standard WGS data from patients #11152, #11513, and OP1036. Red dots (right) are obtained from all artifactual mutations identified by MutLX in DigiPico data from the same patients. The black line represents the median value for each set. The higher proportion of C>A mutations among mutations removed by MutLX is consistent with previous studies showing that DNA oxidative damage during library preparation creates artifactual C>A mutations. [Figure 13] IGV image of SNVs identified in subclonal Kataegis of the PT2R sample. Note that nearly all mutations on the forward strand of the genome are in the form of C>T or C>G mutations. [Figure 14] Comparison of the DigiPico and DigiPico2 workflows. (A) The DigiPico workflow takes approximately 12 hours and consists of 7 steps, 5 of which occur in 384-well plates. (B) The DigiPico2 workflow takes just under 4.5 hours and consists of 5 steps, only 3 of which occur in 384-well plates. The blue reaction occurs in a 384-well plate format, the green reaction occurs in 16 wells, and the orange reaction occurs in one tube. [Figure 15]Comparison of the indexing strategies of DigiPico and DigiPico2. (A) Asymmetric ligation of two annealed indexing oligos introduces one i5 index and one i7 index. Combinations of these indexes can generate 384 indexes in DigiPico. Each strand receives one i5 index and one i7 index, so there is no redundancy. (B) In DigiPico2, a tandem index (Ci) is first introduced by efficient ligation of loop adapters. Next, index PCR is performed using row index primers (Ri). Note that each strand receives two Ci and Ri indexes, respectively, which introduces the redundancy necessary for the removal of index cross-contamination. [Figure 16] Comparison of DigiPico and DigiPico2 results. (A) Both WGA products appear relatively uniform across the plate. Values ​​represent relative fluorescence (RFU) as measured by Evagreen. (B) The frequency of indexes from each well across the plate appears to correlate better with the RFU values ​​of the WGA products in DigiPico2. (C) This fact can be quantified using a correlation plot. (D) MutLX analysis of DigiPico2 data appears to yield better discrimination between true and artifactual mutations. Note the presence of fewer mutations in the upper right part of the plot. This region marks artifactual mutations that received a high probability score in the analysis. [Figure 17]Evaluation of the single-cell DigiPico (ScDigiPico) sequencing method. To evaluate the effectiveness of the ScDigiPico method in identifying active mutational processes, we mimicked such processes by mutagenizing cultured Kuramochi cells with N-ethyl-N-nitrosourea (ENU). ENU is an alkylating agent and an extremely potent mutagen, selectively inducing T>C, T>A, and C>T mutations. In this setting, each cell acquires a different set of mutations; however, because the underlying mutagenesis mechanism is the same, the mutations are expected to be of the same type. (A) Cultured Kuramochi cells were exposed to 0.1 g / L ENU for 48 hours. Most cells die in the presence of ENU, but surviving cells accumulate numerous mutations. Single-mutagenized Kuramochi cells were then sorted into the first column of every well in a 384-well plate. During ScDigiPico, the amount of DNA from each cell is evenly distributed among all row wells before WGA and library preparation. (B) Analysis of the ScDigiPico results from mutagenized cells clearly shows that the frequency of novel mutations matches the ENU results, as expected. (C) In only one of the analyzed cells after UV irradiation, we identified a clear kataegis event on chromosome 7. [Figure 18] DigiPico / MutLX removes false positives from whole genome amplified DNA. Dots represent individual samples (approximately 20 cancer cells) from cancer patients or from whole genome amplified blood DNA (starting with picograms of blood DNA). Germline (blood) DNA should not have thousands of unique variants when compared to standard DNA sequencing methods. However, existing methods find tens of thousands of such false positive mutations. In comparison, DigiPico / MutLX removes these false positives. DETAILED DESCRIPTION OF THE INVENTION

[0085] Example 1 - Using DigiPico / MutLX to reveal active mutational processes in tumors with unprecedented precision overview Bulk whole-genome sequencing (WGS) enables the analysis of tumor evolution, but due to depth limitations, it can only identify ancient mutational events. Discovery of current mutational processes to predict tumor evolutionary trajectories requires dense sequencing of individual clones or single cells. However, such investigations are inherently problematic when sequencing picogram quantities of DNA, resulting in excessive false-positive mutations. Data pooling to increase the confidence in discovered mutations can trace them back to a past common ancestor. Here, we report a robust whole-genome sequencing and analysis pipeline (DigiPico / MutLX) that virtually eliminates all false-positive results while retaining a high proportion of true-positives. Using our method, we identified for the first time hypermutation events in approximately 30 cell populations derived from recurrent ovarian cancer, which cannot be identified from bulk WGS data. Overall, we propose the DigiPico / MutLX method as a powerful framework for identifying clonally specific variants with unprecedented precision.

[0086] Introduction In this study, we developed a single-DNA molecule WGA and sequencing method to obtain high-quality, data-rich sequencing results from picogram amounts of DNA obtained from clinical samples (we called it DigiPico; short for Digital Sequencing of Picogram DNA). Furthermore, we implemented a complementary analysis workflow for DigiPico data using an artificial neural network (ANN)-based algorithm (MutLX, short for Mutation Learning) to remove false-positive results while maintaining excellent sensitivity for true-positive mutations at the whole-genome scale. We validate our method using data from a tumor extensively sequenced to a cumulative depth of approximately 4,200x across 45 whole-genome sequencing experiments on DNA from three different time points, derived from a single patient. We demonstrate the versatility of our method with sequencing samples from four additional cancer patients and a lymphoblastoid cell line.

[0087] material and method Patient samples and consent Patients #11152, #11502, and #11513 provided written consent to participate in the prospective biomarker validation study, Gynaecological Oncology Targeted Therapy Study 01 (GO-Target-01), under research ethics approval number 11 / SC / 0014. Patient OP1036 participated in the prospective Oxford Ovarian Cancer Predict Chemotherapy Response Trial (OXO-PCR-01), under research ethics approval number 12 / SC / 0404. Necessary informed consent from study participants was obtained as needed. Blood samples were obtained on the day of surgery. Tumor samples were collected for biopsy during laparoscopy or bulk removal surgery and immediately frozen on dry ice. All samples were stored in clearly labeled cryovials in a -80°C freezer.

[0088] cell line The GM12885 lymphoblastoid cell line (RRID:CVCL_5F01) was obtained from the Coriell Institute and cells were maintained in culture as recommended by the provider.

[0089] Sectioning and LCM Frozen tumor samples were embedded in OCT (NEG-50, Richard-Allan Scientific) and 10-15 μm sections were taken on a CryoStar cryostat microtome (Thermofisher Scientific) using a DynaSharp microtome blade (Thermofisher Scientific). Tumor sections were then transferred to PEN membrane slides (Zeiss) and immediately stained on ice (70% ethanol for 2 minutes, 1% cresyl violet (Sigma-Aldrich) in 50% ethanol for 2 minutes), followed by rinsing in 100% ethanol. Individual tumor islets were then placed into 200 μl opaque adhesive caps (Zeiss) using a PALM laser microdissection system (Zeiss).

[0090] Standard WGS and data analysis DNA was extracted using the DNeasy Blood and Tissue Kit (Qiagen). For fragmentation, up to 1 μg of DNA was diluted in 50 μl of water using a Covaris S220 focused-ultrasonicator to obtain 250–300 bp fragments. The resulting DNA fragments were then used for library preparation using the NEBNext Uitra II Library Preparation Kit (NEB) according to the manufacturer's protocol. The resulting libraries were sequenced against the human genome at 30–40x depth using an Illumina NextSeq or HiSeq platform. Sequencing reads in FastQ format were initially trimmed using TrimGalore (14) and then mapped to the human hg19 genome using Bowtie2 (15). Germline variant calling was performed using GATK's HaplotypeCaller (16). Somatic variants were called using Strelka2 with a variant allele ratio cutoff of 0.2 (17).

[0091] DigiPico sequencing Two hundred micrograms of purified DNA, 20–30 resuspended nuclei, or laser-capture microdissected tumor islets were first denatured using 5 μl of D2 buffer from the Repli-g single cell kit (Qiagen). After a 5-minute incubation at room temperature, 95 μl of water was added to the sample, and then 200 μl of denatured template was added to each well of a 384-well reaction plate containing 800 nl of WGA mixture (0.58 μl of Sc reaction buffer, 0.04 μl of Sc polymerase (Repli-g single cell kit, Qiagen), 0.075 μl of 1 mM dUTP (Invitrogen), 0.04 μl of Evagreen 20x (Biotium), and 0.065 μl of water) using a Mosquito HTS liquid handler (TTP Labtech). The plate was incubated at 30°C for 2 hours, followed by heat inactivation at 65°C for 15 minutes. Addition of Evagreen during the reaction allowed for monitoring of the WGA reaction using a real-time PCR instrument, if desired (18). A controlled enzymatic fragmentation (19) reaction step was then performed sequentially on the whole genome amplified DNA without any purification steps. Briefly, (A) 1200 nl of UDG mixture (0.08 U / μl rSAP (NEB), 0.2 U / μl UDG (NEB), 0.4 U / μl EndoIV (NEB) in 1.8x NEBuffer 3) was added, followed by a 2-hour incubation at 37°C and a 15-minute heat inactivation at 65°C. (B) 1200 nl of Pol I mixture (0.4 U / μl DNA polymerase I (NEB), 0.25 mM dNTPs, 8 mM MgCl2, and 0.8 mM DTT) was added, incubated at 37°C for 1.5 hours, and heat-inactivated at 70°C for 20 minutes. (C) 1200 nl of Klenow mixture (0.5 U / μl Klenow exo-(NEB), 0.5 mM dATP, 8 mM MgCl2, and 0.8 mM DTT) was added, incubated at 37°C for 45 minutes, and heat-inactivated at 70°C for 20 minutes.(D) 400 nl of 20 μM full-length Illumina adapter oligos with well-specific indexes (Table S1) were added to each well, followed by 1100 nl of ligation mixture (40 U / μl T4 DNA ligase (NEB), 5 mM ATP, 11.5% PEG8000 (Qiagen), and 6.8 mM MgCl2), followed by 30 min of incubation at 20°C and 15 min of heat inactivation at 65°C.

[0092] The resulting products were then pooled, and the DNA was precipitated with an equal volume of isopropanol. The DNA was then resuspended in water and subjected to double size selection using Agencourt AMPure XP SPRI magnetic beads (Beckmann Coulter) with a 0.45x bead ratio for left selection and an additional 0.32x bead ratio for right selection. The purified DNA was then resuspended in water and immediately used for limited-cycle PCR amplification using the P5 and P7 primer mix (Table S1). PCR was performed for 12 cycles with 10 seconds of annealing at 55°C and 45 seconds of extension at 72°C. The final product was bead-purified at a 0.9x ratio. The resulting library was sequenced using an Illumina sequencing platform in 2x150 paired-end sequencing mode, achieving 30-40x coverage depth against the human genome. Currently, the additional processing steps required for DigiPico library preparation add approximately £250 to the total reagent cost.

[0093] Analysis of DigiPico sequencing data The analysis pipeline for DigiPico sequencing data is shown in Supplementary Figure 5. Briefly, 384 paired-read FastQ files for each well were obtained after demultiplexing of the Illumina sequencing data. FastQ files were trimmed for adapter sequences and quality (14). The first 12 nucleotides of each read were also removed. Trimmed reads were mapped to the human hg19 reference genome using Bowtie2 (15) with the ignore-quals parameter enabled and duplicate reads marked using Picard Tools (20). Joint variant calling was performed on all 384 individual bam files along with the merged bam file from all wells using the Platypus variant caller (21). Next, all low-quality variants were removed by applying quality filters (QUAL > 60, FR > 0.1, HP ≤ 4, QD > 10, and SbPval ≤ 0.95). Additionally, the total number of wells covering each locus (Tw) and the number of wells carrying each variant (Vw) were determined, and well-counting filters (Tw > 5, Vw > 2, and Vw / Tw > 0.1) were applied to retain only high-confidence positions for analysis. Finally, all regions of the genome with poor mappability (22) were removed from the analysis, e.g., using VCFtools (23). The resulting list of high-confidence novel DigiPico variants was then used to perform variant recall (genotyping) on ​​WGS data from blood and bulk tumors using Platypus with the minPosterior parameter set to 0 and the minMapQual parameter set to 5. Any variants not reliably retained in both standard WGS data were extracted as UTD (Unique to DigiPico) variants. Any variants that were reliably present in the bulk sequencing data of blood samples (based on GATK analysis) were also extracted as TP (True Positive) variants (Figure 5).

[0094] MutLX analysis algorithm The MutLX analysis pipeline is summarized in Figure 6. Artificial Neural Network Structure The neural network model used in this study was a multilayer perceptron with an input layer consisting of N neurons (N = 41), where N is the number of features used in each experiment. It was implemented in Python 3 using Keras (24). The model has two hidden layers with ReLU activation. We varied these numbers, but did not observe significant improvement when using a larger number of neurons. The final layer is a single-output neuron with sigmoid activation. The loss function is binary cross-entropy. For training, we applied stochastic gradient descent optimization with momentum (Adam (25)) with a learning rate of 0.001, a batch size of 8, and 10 epochs. We did not observe any additional improvement in performance after 10 epochs.

[0095] Features used for training The following features extracted from the Platypus output of the DigiPico data were used as inputs for the neural network model: Platypus quality parameters: QUAL, BRF, FR, HP, HapScore, MGOF, MMLQ, MQ, QD, SbPval, NF, NR, TCF, and TCR (21).

[0096] Sequence context complexity: F 20 [1], F 20 [2], F 20 [3], F 20 [i] is the sum of the frequencies of the i most abundant nucleotides in the 10 bp sequences on either side of the variant position.

[0097] Lead distribution data: R merge [ref+var], R merge[var], W[R[ref]>0s & R[var]=0], W[R[ref]>0][0 / 0], W[R[ref]>0][0 / 1], W[R[var]>0], W[R[var]>0][1 / 1], W[R[var]>0][0 / 1], W[R[var]>0 & R[ref]=0][0 / 1], W[R[ref]>0 & R[var]>0], W[R[ref]=0 & R[var]>0], W[R[ref]=0 & R[var]>1], W[R[ref]=0 & R[var]>2], W[R[ref]=0 & R[var]>3], W[R[ref]=0 & R[var]>4], W[R[ref]=0 & R[var]>5], R max [1][var], R max [2][var], R max [3][var], R max [1][ref+var], R max [2][ref+var], R max [3][ref+var], Max c +Max r , W[R[var]>0]-(Max c +Max r ). In the formula, R merge [x] denotes the total number of reads in the merged bam file that carry allele x (ref denotes the reference allele, var denotes the variant allele). w[i][j] denotes the number of wells matching criterion i with reported genotype j (if indicated). For criterion i, R[x] denotes the number of reads in a specific well that carries allele x. R max [y][x] indicates the number of reads in the well with the yth largest number of reads carrying allele x. Finally, Max c is the number of variant-carrying wells in the column with the maximum number of wells carrying the variant allele, and Max r is the number of variant-bearing wells in the row with the maximum number of wells bearing the variant allele.

[0098] Training with MutLX For each DigiPico experiment, we consider a subset of all UTD variants (labeled 0) and heterozygous germline SNPs (labeled 1) in the full training set. The number of UTD variants in this set is much smaller than the number of heterozygous germline SNPs, making the set unbalanced. Therefore, to avoid bias toward specific labels during training, we create 25 different unbiased training subsets for each DigiPico experiment. This consists of a randomly selected subset of heterozygous germline SNPs, with each training subset having a size equivalent to the number of all variants and UTD variants. As explained previously, most UTD variants are FP variant calls with an unknown proportion of true clone-specific variants within them, thus making the 0 label noisy. To perform two-stage training that accounts for these noisy labels, we adopt the following strategy. After training an initial model on each unbiased training subset, the resulting model is applied to the variants in the full training set to obtain initial probability values ​​for each variant. These probability values ​​indicate the predicted probability of a variant belonging to the 1-label category. Therefore, any 0-label variant with a predicted probability value close to 1 is likely to be a mislabeled variant. Therefore, to reduce the level of mislabeled data in the training set, all UTD variants with a probability value greater than 0.7 and all germline SNPs with a probability value less than 0.3 are considered mislabeled and removed from the training set. The cutoff values ​​in this step were determined experimentally by analyzing various simulated data sets. Finally, a new model is trained on the remaining variants in the training set, following a similar subsampling strategy as in the initial training. This model is then used to analyze all UTD variants.

[0099] Calculating probability and uncertainty scores As previously described, in MutLX, the training process is repeated 25 times with different randomly selected germline SNP subsets, yielding a different model each time, and thus 25 predicted probability values ​​for each variant. We therefore defined a "probability score" for each variant as the average of all of its predicted probability values:

number

[0100] Furthermore, to obtain uncertainty estimates for each probability value, we performed a test-time dropout analysis (26). The trained model was applied to each mutation 100 times, during which different neurons were dropped out at rates of 0.8 and 0.7 for the first and second hidden layers of the neural network, respectively. This process yielded 100 probability values ​​for each mutation. Based on these values, we defined an "uncertainty score" for each mutation as the average dropout variance from the 25 subsets.

number

[0101] An estimated receiver operating characteristic (ROC) curve was generated using the uncertainty scores of all variants with a probability score greater than 0.2 (Figure 7). This curve was generated by considering cutoff values ​​ranging from 0.0 to 0.25 for the uncertainty score. For each cutoff, the proportion of germline SNPs with uncertainty scores below the cutoff was plotted against the corresponding number of UTDs. The area under the curve (AUC) was then calculated after normalizing the number of UTDs to 0–1. Note that if no true clonally specific variants are predicted (all UTDs are FP calls), this plot represents an ROC curve, and the AUC of this plot should approach 1, assuming a perfect model. In contrast, a significantly smaller AUC indicates the presence of true clonally specific variants in the sample. This negative correlation between the number of true UTDs and AUC was verified using a mock dataset (Figure 8). Based on these observations, for example, when the ROC curve suggests the presence of true clone-specific variants (AUC<0.9), MutLX uses an "uncertainty score" cutoff value that results in a TPR of 95% to improve the recovery rate of clone-specific variants. For datasets with AUC≥0.9, the cutoff value for filtering the data was determined based on the intersection of the threshold curve and the ROC curve.

[0102] Generation and analysis of simulated DigiPico datasets Using mock data, we (a) verified the negative correlation between the number of true UTDs and AUC (Figure 8) and (b) ensured that overfitting to potential true clone-specific variants did not occur (Figure 9).

[0103] To generate a mock dataset, we first identified somatic mutations in bulk WGS data of tumor sample PT2R from patient #11152 using the Strelka2 somatic variant caller. These somatic variants were then identified in the new variant call data of experiment D1110, and all somatic variants with Tw > 6 and Vw / Tw > 0.45 were selected as high-confidence somatic variants. We then combined various numbers of randomly selected high-confidence somatic variants into UTDs (UTDs). * ) and UTDs of 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, and 0.1. * The resulting composite list of variants was then used independently for MutLX analysis, and the UTD and UTD ratios that passed MutLX filtering were calculated. * The number of somatic variants was calculated for each experiment. To ensure robust analysis, 10 different subsets of somatic variants were analyzed for each ratio. A similar analysis was also performed on DigiPico data DE111, derived from bulk DNA extracted from the ascites sample of patient #11513.

[0104] Validation of the MutLX analysis algorithm Tumor sample PT2R from patient #11152 was used to validate the MutLX algorithm. Small pieces of tumor were macro-dissected from frozen specimens and embedded in OCT medium for sectioning. The first section (15 μm) from the tumor was collected in a separate tube, and the nuclei were resuspended in 50 μl of sterile PBS solution. The total number of nuclei in the suspension was determined, and a volume containing 30 nuclei was used for direct denaturation using an equal volume of D2 buffer from the Repli-g mini WGA kit (Qiagen). The resulting crude lysate was used directly for DigiPico library preparation for experiment D1111. The remainder of the tumor sample was then used for bulk DNA extraction using the DNeasy Blood and Tissue Kit (Qiagen). 200 pg of the resulting DNA was directly used to prepare DigiPico library D1110. 1 μg of DNA was used for standard library preparation using the NEBNext Ultra DNA Library Preparation Kit (NEB). In this setting, only experiment D1111 is expected to have true clone-specific variants. Because the templates in experiment D1110 are a subset of those used in bulk WGS analysis, nearly all actual variants in experiment D1110 will also be present in the WGS data at similar frequencies and therefore will not be identified as UTDs. Similar logic can be applied to the results of DigiPico experiments DE011 and GM12885. Because both of these DigiPico experiments were performed on 200 pg of DNA from bulk DNA extracts, true UTD variants are not expected to be present in these samples. It is also worth noting that, due to the digitized nature of the data, variants with extremely low frequency (<0.05%) show increased variant allele proportions in experiments D1110, DE011, and GM12885; however, because such variants are unlikely to appear in more than one well, they will be removed from the data based on the Vw filter. Therefore, it is safe to assume that nearly all UTD variants in these experiments are FP calls.

[0105] Applying SCcaller to DigiPico data SCcaller was originally developed for the analysis of multi-displacement amplified single-cell sequencing data (11). Because DigiPico library preparation also requires multi-displacement amplification of a limited amount of template DNA, the resulting data is essentially similar to the natural input for SCcaller. Therefore, we used the merged bam file of DigiPico data as input for SCcaller. For analysis, we used the GATK HaplotypeCaller to obtain a list of heterozygous SNPs from each bulk WGS data set, with a cutoff value of α = 0.01. Next, all filtered SNVs were used for variant recall against each standard WGS data set, and all variants that were not reliably retained in the WGS data set were extracted as UTD variants.

[0106] Mutation validation Variants that passed the MutLX analysis were validated by comparison with deep sequencing data from bulk tumors from independent sequencing platforms. All DigiPico data from patient #11152 were validated by comparison with 39 deep sequencing datasets obtained from the same tumors sequenced on the Complete Genomics sequencing platform (27). This included three Complete Genomics bulk sequencing datasets and 36 LFR (long fragment read) sequencing datasets. Because the independent sequencing data from the omental tumors were not obtained from the exact same tumors used for DigiPico sequencing, the validation rate from such comparisons for these experiments is not expected to be high.

[0107] For targeting validation, primers were designed to obtain amplicons containing variants using the Primer3 tool (Table S1). Amplicons were obtained by performing 16 cycles of two-step PCR with 1 ng of template using Phusion® High-Fidelity PCR Master Mix with GC buffer. All amplicons from each sample were then pooled and purified prior to adapter ligation and indexing using the NEBNext Ultra II kit. The resulting libraries were sequenced on the MiSeq platform. Sequencing results were mapped to the human hg19 genome using Bowtie2, and the number of reads carrying each variant was counted using the Platypus variant caller.

[0108] Local hypermutation (kataegis) analysis To generate the rainfall plots, the distances between consecutive pairs of somatic mutations on chromosome 17 were plotted against the genomic location of the second mutation in each pair using a custom script in R. The presence of clusters of localized mutations indicates a kataegis event. In these plots, each dot was colored based on the mutation type of the second mutation in the pair with respect to the hg19 human reference genome.

[0109] result Implementation of the DigiPico sequencing method A key feature of amplification errors, both with and without prior DNA damage, is that they are randomly introduced during the amplification process (6, 7, 28). Therefore, we hypothesized that when a single DNA molecule is amplified and sequenced, artifactual mutations should be present in only a small fraction of the reads obtained from sequencing the original single DNA molecule. In contrast, genuine variants should be expected to be present in all such reads. Dividing the template DNA into separate compartments prior to WGA, such that each compartment contains only one DNA molecule from every locus, would yield such single DNA molecule sequencing data (Figure 1B). Because artifactual mutations can result in compartments with reads carrying multiple alleles, this approach allows for the identification of these artifacts and the removal of FP variant calls. In addition, such a partitioning approach also allows for independent running of WGA reactions for each locus, thus providing multiple internal replicates of data for the WGA process. Genuine variants are expected to be regularly distributed across replicates, whereas artifactual mutations, due to their stochastic nature, may be limited to a few partitions. Therefore, taking both of these points into consideration, WGA and sequencing techniques may produce distinctive distribution patterns across partitions of artifactual mutations compared to real mutations. ANNs have shown the ability to extract complex patterns from high-dimensional inputs, making them good candidates for identifying and removing false-positive mutations from this type of data. Previous partitioning and sequencing methods have been described to obtain haplotype information, but no such methods exist for distinguishing true mutations from artifactual mutations (19, 29, 30).

[0110] To fully benefit from the data richness of partitioning and sequencing techniques for accurate genomics testing of clinical samples, we developed the DigiPico sequencing method (Figure 1C). To perform DigiPico sequencing, we first uniformly distribute approximately 200 pg of DNA (20–30 human cells) into individual wells of a 384-well plate. This ensures that the probability of two different DNA molecules from the same locus coexisting in the same well is less than 10% (19). After WGA, each well is independently processed into an index library and receives a unique barcode sequence before pooling and sequencing (Figure 1C). In our method, key distinguishing factors for artifactual mutations are their distribution pattern, uniform distribution, and DNA molecule amplification. In addition, consistent sequencing coverage across the wells is crucial. Achieving this uniformity ensures that the differences in the distribution patterns of true and artifactual mutations are maximized. To ensure the required uniformity was achieved, during all DigiPico library preparations, we monitored the progress of the WGA reaction and quantified the final results for every well. The former was achieved by adding Evagreen dye to the WGA reaction and monitoring the fluorescence intensity in real time every 5 min. Evagreen is an intercalating dye that binds to the minor groove of DNA and therefore does not interfere with isothermal WGA reactions (18). For the latter, we introduced a per-well qPCR step to measure the relative number of adapter-ligated fragments in each well using adapter-specific primers before pooling. Only libraries that passed both of these uniformity tests were used for sequencing (Figures 1D and 1E). Importantly, we also miniaturized the WGA reaction volume to 1 μl, which could only be achieved after identifying a compatible multiple displacement amplification (MDA). Comparing six different MDA strategies, REPLI-g Single Cell Amplification was the only method that met the sensitivity and selectivity required for our purpose (Figure 1F).Reaction miniaturization allowed us to simplify the library preparation method in 384-well plates using readily available automated pipetting equipment without the need for intermediate purification steps. Finally, we aimed to optimize the DigiPico library preparation method for frozen clinical samples. This was achieved by performing WGA reactions directly on crude lysates of small, contiguous cell populations (tumor islets) isolated via LCM (laser capture microdissection). This strategy ensured minimal loss of genomic material while simultaneously minimizing handling time, thereby reducing the chance of template oxidation.

[0111] The DigiPico sequencing platform generates high-quality libraries from limited clinical samples. Having optimized all the necessary aspects of the DigiPico library preparation method, we decided to evaluate the quality of DigiPico libraries obtained from clinical samples. To this end, we prepared DigiPico libraries D1110 and D1111 from a frozen recurrent tumor sample (PT2R) obtained from a patient with high-grade serous ovarian cancer (#11152). In this experiment, the D1110 library was prepared from 200 pg of template taken from a bulk DNA extract of the PT2R sample, while the D1111 library was prepared directly from a small remaining frozen section of this tumor sample (containing approximately 30 cancer cells). Each library was sequenced on the Illumina NextSeq platform, yielding approximately 400,000,000 reads in a 150x2 paired-end format. Initial evaluation of the resulting sequences revealed that both the D1110 and D1111 libraries yielded high-quality sequencing data with an overall mapping rate of 91.35% and 94.27% to the human hg19 genome, respectively (Figure 1G). Plate-wide uniformity analysis showed that, on average, each well covered approximately 4.4% and 4.6% of the genome in experiments D1110 and D1111, with an average depth of 1.7x and 2.1x for each well, demonstrating excellent uniformity across the plate (Figures 1H and 1I). This cumulatively resulted in 92.1% and 91.1% coverage widths and 30x and 43x coverage depths for each experiment, respectively. These results demonstrate that the DigiPico sequencing method can be used to generate high-quality sequencing data with excellent coverage from limited amounts of frozen clinical samples.

[0112] Finally, we evaluated whether our initial hypothesis regarding the distinctive distribution patterns of different variants holds in the actual DigiPico dataset. To this end, we assumed that any variants common between the DigiPico dataset and standard bulk sequencing data for the same tumor sample must be true variants. These should primarily consist of germline SNPs and clonal somatic variants. Consequently, by definition, all FP variant calls and the majority of clone-specific mutations (if they are present in the sample under test) will be among variants present only in the DigiPico data and absent in the bulk WGS data. These variants are hereafter referred to as UTDs (Unique to DigiPico) for simplicity. Therefore, considering that the standard bulk sequencing data for the PT2R sample was obtained from the same DNA extract used for D1110 library preparation, nearly all UTD variants in the D1110 DigiPico experiment should be artifacts (Figure 1J). In contrast, UTD in experiment D1111 likely contains some clone-specific mutations along with artifactual mutations (Figure 1J). Therefore, we used the UTD variant in experiment D1110 as a representative artifactual mutation in our analysis. Comparison of the frequency of wells with coexistence of two alleles at the same locus in experiment D1110 (Figure 2A) as well as the number of wells harboring each variant (Figure 2B) showed that UTD had a significantly higher proportion of the former and a lower number of the latter than any other category of mutations (p<2e-16 for both analyses, one-way ANOVA followed by Tukey HSD test). This clearly supports the expected distinct distribution pattern of artifactual mutations in this DigiPico dataset.

[0113] MutLX analysis pipeline for DigiPico data After obtaining high-quality data using DigiPico sequencing, we decided to implement an analysis pipeline to remove FP variant calls based on mutation distribution patterns. As previously mentioned, ANN algorithms are ideally suited to problems involving such complex patterns. Given a representative set of precisely labeled examples (training set), ANNs can learn to classify mutations without requiring any class-specific information. However, implementing ANN algorithms for the problem of removing FP mutations from sequencing data presents two major challenges: (a) the difficulty of obtaining a generalizable model, and (b) the unavailability of a representative, precisely labeled training set. First, because mutation distribution patterns depend on various experiment-specific initial conditions (e.g., genome copy number state) that cannot be easily explained, generating a generalizable model for the analysis of all DigiPico datasets is not possible. Therefore, experiment-specific models tailored to each DigiPico experiment are required. This means that a subset of experiment-specific mutations needs to be selected as the training set for each DigiPico experiment. Second, while accurately labeled examples of true mutations can be easily extracted from known SNPs in the genome, identifying a representative and accurate set of examples of artifactual mutations is not possible. To address this issue, we assumed that UTDs are primarily composed of such mutations and considered the UTDs as a reasonable approximation for a representative set of artifactual mutations. However, this assumption can create significant challenges. By definition, UTDs consist of artifactual mutations as well as true clone-specific mutations. While artifactual mutations are abundant in DigiPico experiments, true clone-specific mutations may be present at different frequencies depending on the sample (Figure 1J). Therefore, if UTDs are considered entirely as examples of artifactual mutations, samples with more clone-specific variants will have noisy training sets. If this is not taken into account, samples with true clone-specific variants may be analytically disadvantaged because a noisy training set may produce a poor classification model (Figure 2C).In particular, the presence of actual mutations in artifactual variants in samples containing true clone-specific variants can result in overfitting of the model to these variants (Figure 2C). This can impair the model's ability to accurately identify true clone-specific variants in such samples. Given that the primary goal of DigiPico sequencing is to identify clone-specific variants, it is essential to ensure that such overfitting does not occur when analyzing DigiPico datasets.

[0114] Considering all the above limitations and issues, we designed and implemented an ANN-based binary classifier, MutLX, for the analysis of the DigiPico dataset. The goal of the DigiPico analysis pipeline was to effectively remove FP calls and accurately identify true clone-specific variants from UTDs. To address the issue of training with an incomplete training set, we adopted the following approach in training MutLX. First, we considered all UTD variants as examples of artifactual mutations (labeled as 0s) and a similar number of randomly selected heterozygous germline SNPs as examples of true variants (labeled as 1s). Since the majority of UTD variants are FP calls with an unknown proportion of true clone-specific variants, 0 labels are considered "noisy" at this stage. In other words, true clone-specific variants should be labeled as 1s, but because of their anonymity at this stage, they are specifically labeled as 0s. To accommodate this type of noise in the training dataset, we adopted a two-stage training process (Figure 2D). In the first step, a model is trained on all the example data for a given label. This initial model is then used to calculate the probability of each mutation belonging to that label category. Any mutations that appear to be mislabeled based on these model predictions are temporarily removed (pruned) from the dataset. The underlying assumption is that a model trained on noisy samples may not be as robust as one trained on a hypothetical clean dataset, and even then, it will be biased toward better prediction of accurate examples due to their higher proportion relative to noise in the dataset. Thus, examples that the model predicts for their original label may have been mislabeled in the first place. In the second step, a new classification model is trained on this pruned training set. Because the second training set is likely to contain fewer mislabeled data points, the final model is expected to be more effective at identifying true mutations, regardless of whether all UTD variants were actually artifacts (31, 32). We then use this model to assign a "probability score" to each putative mutation.This score indicates the likelihood that a particular mutation belongs to the true variant category. While this two-stage training process is expected to significantly improve the classification model, the final model is still prone to errors due to incompleteness of the training set. Therefore, we added another level of analysis to further improve the accuracy of our pipeline. This was achieved by assigning an uncertainty estimate to each mutation's "probability score." This uncertainty estimate is based on the assumption that robust predictions are supported by a large proportion of activated neurons in the hidden layer of the ANN. Therefore, any subset of these neurons will also consistently produce similar probability scores, resulting in low variance among the various "probability scores" obtained from different neuron subsets (Figure 2E). In contrast, the seemingly high "probability scores" of artifactual mutations are likely supported by only a small proportion of neurons in the hidden layer of the ANN. Therefore, different neuron subsets will produce varying "probability scores," leading to a large variance in the scores obtained from various neuron subsets for artifactual mutations (Figure 2E). As a result, an “uncertainty score” can be calculated as the variance of the “probability scores” obtained from multiple different randomly selected neuron subsets of the trained ANN in MutLX. (26) The combination of the “probability score” and the “uncertainty score” for each mutation thus allows us to accurately determine whether a called variant is a real mutation in the template or the result of an artifactual change (Figure 6).

[0115] Verification of the MutLX algorithm To validate our strategy, we chose to test the MutLX analysis pipeline on experiments D1110 and D1111. This is because these DigiPico experiments were derived from a previously extensively sequenced HGSOC (patient #11152) with available data from 48 independent whole-genome sequencing datasets spanning three different time points at a total depth of approximately 4200x from two independent sequencing platforms (33). To our knowledge, this comprises the most extensively whole-genome sequenced tumor to date. This exceptionally large dataset allows for reliable cross-validation of mutations in this tumor. To this end, we used the MutLX algorithm to analyze the sequencing data from experiments D1110 and D1111. As previously described, when using bulk sequencing data from the PT2R site to compare with these DigiPico datasets, true UTD variants (clone-specific variants) are only predicted to be present in experiment D1111, whereas nearly all UTDs in experiment D1110 are predicted to be artifacts (Figure 1J). In addition, we also analyzed DigiPico sequencing data from purified DNA of a blood sample (experiment DE011) as well as cultured GM12885 lymphoblastoid cells. Both of these are also predicted to have no true UTD mutations. De novo variant calling for these DigiPico experiments, followed by initial filtering based on well counts, resulted in the identification of thousands of UTD variants, nearly all of which were predicted to be FP calls. However, application of the MutLX algorithm to the UTD variants in experiments D1110, DE011, and GM12885 effectively reduced over 99% of FP variant calls to just 4, 7, and 3 genome-wide FP mutations, respectively, while maintaining approximately 85% sensitivity for detecting true mutations (Figures 2F and 2G). In comparison, SCcaller (11) analysis of the same data yielded 713, 712, and 13,280 FP variant calls, respectively (Figure 2G).On the other hand, MutLX identified 264 putative clone-specific variants in experiment D1111, and 238 of these (90%) were validated by comparison with an independent high-depth dataset of this tumor sample (Figure 2G). These observations were further validated by performing targeted sequencing on tumor bulk DNA, which revealed that 10 of the 11 analyzed amplicons containing clone-specific variants from experiment D1111 were reliably present at low frequency in the bulk DNA of the PT2R sample (Figure S7). Furthermore, 37 ostensibly high-quality UTD variants from experiment DE111 that were labeled as artifacts by the MutLX algorithm showed no evidence of their presence in the bulk DNA sample (Figure 11). These results clearly confirm that MutLX can distinguish artifactual mutations from real variants in DigiPico data and train accurate classification models that can efficiently identify true clone-specific variants.

[0116] Furthermore, we investigated whether the presence of true clone-specific mutations could impair the model's sensitivity through overfitting. To this end, we introduced various numbers of somatic mutations into the artificial UTD variants (UTDs) during experiments D1110 and DE111. * ) to generate synthetic datasets containing various proportions of true UTDs. These synthetic datasets were then independently analyzed by MutLX to identify the various UTDs in all synthetic datasets. * FP ratio and UTD in / UTD ratio * The results showed a high UTD of 10%. * / UTD ratio, UTD * We show that this does not significantly affect the recovery of variants, indicating that overfitting does not occur with MutLX (Figure 9).

[0117] Versatility of the DigiPico / MutLX sequencing and analysis method Finally, to ensure the versatility of our proposed method, we performed DigiPico sequencing on template DNA of various origins from four different HGSOC patients, and analyzed the resulting UTDs using the MutLX algorithm. The results clearly demonstrate that MutLX can reliably identify and remove artifactual variant calls from a diverse set of DigiPico libraries (Table 1). This strongly suggests that DigiPico / MutLX can efficiently enable the testing of recently acquired mutations in solid tumors. Importantly, analysis of the frequency of different variant types in these data showed a higher level of C>A mutations among the identified artifactual mutations, consistent with the notion that such FP calls are the result of oxidative damage to the template DNA (Figure 12).

[0118] Investigating active mutation processes using DigiPico / MutLX. We next tested the feasibility of investigating the mutational process in a patient with HGSOC (#11152). In this patient, various sequencing data were available from a prechemotherapy omental mass (standard bulk sequencing at 30x and a DigiPico experiment of five tumor islands). The patient subsequently relapsed, and tumor samples were collected from the pelvis (pelvic recurrent tumor; PT2R) and para-aortic lymph nodes (PALNR) for standard bulk sequencing of tumor islands and DigiPico sequencing of tumor islands. Analysis of the bulk prechemotherapy sequencing data identified 13,721 somatic mutations. From the DigiPico data, 84.6% of these mutations were present in at least three tumor islands. Similarly, from previously published LFR data (33), 91.4% of these mutations were present in at least three additional islands. The high incidence of mutations indicates that they were early mutations that solidified in the tumor. Analysis of DigiPico data from tumor islands revealed the presence of a limited number of clonally specific mutations not present in the bulk tumor. Five pre-chemotherapy islands contained a high number of truly unique mutations (2, 6, 8, 8, and 36), indicating their recent occurrence compared to the other islands (Figure 3A). Bulk WGS data from the PT2R recurrence showed the emergence of 3,009 novel somatic mutations not present in the pre-chemotherapy bulk sequencing data, DigiPico data, or LFR data. These mutations could have arisen at any time point, as the common ancestor of the omental mass and PT2R recurrence diverged from each other (Figure 3A). Analysis of tumor islands in the recurrence sample from patient #11152 showed that the pelvic recurrence tumor (PT2R) harbored a higher load of clonally specific mutations than the para-aortic lymph node recurrence (PALNR) or the pre-chemotherapy tumor. This observation suggests that the underlying molecular mechanisms underlying SNV mutagenesis may have been recently activated in this patient. Furthermore, analysis of clonally specific mutations in PT2R samples by rainfall plot revealed the presence of a strong subclonal local hypermutation (kataegis) event ( 8 ) on chromosome 17 ( Fig. 3B , 3C , and 11 ).Comparison of the mutations that constituted this kataegis event with the bulk sequencing data, DigiPico data of the pre-chemotherapy omental mass, and LFR data revealed that they were found only in the DigiPico PT2R data, indicating that they were genuine clonally specific mutations.

[0119] Consideration In this study, we present DigiPico / MutLX as an integrated platform for identifying mutations from small cell populations at a genome-wide scale with unprecedented precision. We believe this work provides an important stepping stone for the discovery of current or recent somatic mutational processes occurring in cancer and normal tissues. Understanding current mutational processes provides clues for predicting tumor evolutionary trajectories and, potentially, for interfering with such trajectories. Mutations identified by bulk tumor sequencing must arise at a time point during the long course of a tumor, from initiation to symptoms. In contrast, cell-specific mutations must arise during the limited lifespan of that cell. Similarly, mutations in small clones derived from a single cell are also recent. The age of such mutations cannot exceed the age of the clone, defined by the number of cell divisions required to generate the clone. Studying patterns in cell-specific or small clone-specific mutations may enable the identification of recent or current mutational processes (1). Identifying such processes is highly desirable because they inevitably link to biological or chemical phenomena and thus can yield meaningful mechanistic insights. Identifying these mechanisms has important practical implications because they can potentially be applied for therapeutic intervention or prediction of future tumor behavior. Current state-of-the-art technology does not allow for the direct and precise identification of mutations from individual cells or from individual small clones derived from tumors. DigiPico / MutLX makes this effort possible for the first time.

[0120] To overcome major technical pitfalls, primarily related to false-positive mutations, current methods for single-cell WGS analysis require thorough validation testing (11) or rely on mixed data from multiple cells to obtain reliable shared mutations among cells (12, 34). These cells are then grouped into clones derived from a common ancestor. While such techniques tend to target a more recent common ancestor than bulk sequencing, they are still not ideal because the data derived from these methods do not reflect the mutational processes occurring in existing cells. Furthermore, reducing the sequencing depth per cell to enable sequencing of a large number of cells reduces coverage, which is already compromised by the loss of genetic material during the preparation phase. This increases the number of cells required to infer and identify clones, which in turn pushes the ancestry further back in time. In addition, the lack of information about physical relatedness in single-cell analysis methods reduces the opportunity to group cells that may be derived from a single clone. This widens the gap between the inferred clonal ancestry and the present, making it difficult to identify processes active in cells currently present in the tumor.

[0121] DigiPico / MutLX has the distinct advantage of preserving spatial information. Analysis of spatially related cells preserves physical associations, allowing the hypothesis that physically related cells belong to distinct clones (9). Determining distinct structures that may have arisen from tissue-resident stem cells has also been suggested to identify and analyze clones. For example, it may be reasonable to predict that cells from a single small intestinal crypt or a single endometrial gland originate from tissue-resident stem cells (35, 36). Under these circumstances, each anatomical unit defines clones that either possess or do not possess clonally specific mutations that may be related to mutation-promoting factors. Furthermore, clonal sequencing data can be used to computationally infer subclones and predict more recent events that may have occurred within clones. This is similar to achieving bulk sequencing and analysis, but at the level of single clones consisting of a limited number of cells. Preservation of spatial information is particularly interesting because of recent advances in enabling spatial transcriptomics techniques (37). It is conceivable that the combination of highly accurate DNA sequencing and spatial transcriptomics may enable detailed investigation of tissue genetic and non-genetic heterogeneity. In short, current techniques for analyzing small clones produce numerous false-positive results, making it impossible to obtain direct, accurate clone-specific information at a genome-wide scale without exhaustive validation. Combining data from multiple clones is a common solution, but it stretches the ancestry further back in time. We have previously used this approach for the analysis of small tumor cell groups (tumor islands) (33). Due to the uncertainty associated with calling mutations from individual islands, it was necessary to call only mutations common to all tumor islands, effectively identifying only truncal mutations. This was followed by independent validation of approximately 700 mutations using targeted sequencing. Although this provided important biological insights, we were unable to investigate island-specific mutations. DigiPico / MutLX now enables such investigation.We demonstrated how direct analysis of DNA from approximately 30 cancer cells resulted in the reliable identification of subclonal kataegis events.

[0122] Overall, here we demonstrate that DigiPico and MutLX, a significant improvement over existing methods, can enable ultra-high-precision identification of somatic mutations from limited numbers of cells obtained from clinical samples. Furthermore, unlike other computational methods that rely on diploid regions of the genome to calculate amplification bias, our method is also compatible with genomes plagued by extensive copy number alterations, such as those in HGSOC. We believe that the versatility of the DigiPico / MutLX method will enable the investigation of active mutational processes in tumors as well as normal tissues.

[0123] Availability The source code for MutLX is available on Github (https: / / github.com / mmdknr / DigiPico).

[0124] Accession number All sequencing data used in this study are available at EGA (EGAD00001005118).

[0125] References 1.Turajlic, S., Sottoriva, A., Graham, T. and Swanton, C. (2019) Resolving genetic heterogeneity in cancer. Nat. Rev. Genet., 10.1038 / s41576-019-0114-6. 2.Zhang,J.,Spath,S.S.,Marjani,S.L.,Zhang,W.and Pan,X.(2018)Characterization of cancer genomic heterogeneity by next-generation sequencing advances precision medicine in cancer treatment.Precis.Clin.Med.,1,29-48. 3.Gerstung,M.,Jolly,C.,Leshchiner,I.,Dentro,S.C.,Gonzalez,S.,Mitchell,T.J.,Rubanova,Y.,Anur,P.,Rosebrock,D.,Yu,K.,et al.(2017)The evolutionary history of 2,658 cancers.bioRxiv,10.1101 / 161562. 4.Barber,L.J.,Davies,M.N.and Gerlinger,M.(2015)Dissecting cancer evolution at the macro-heterogeneity and micro-heterogeneity scale.Curr.Opin.Genet.Dev.,30,1-6. 5.Bohrson,C.L.,Barton,A.R.,Lodato,M.A.,Rodin,R.E.,Luquette,L.J.,Viswanadham,V.V,Gulhan,D.C.,Cortes-Ciriano,I.,Sherman,M.A.,Kwon,M.,et al.(2019)Linked-read analysis identifies mutations in single-cell DNA-sequencing data.Nat.Genet.,10.1038 / s41588-019-0366-2. 6.Chen,L.,Liu,P.,Evans,T.C.J.and Ettwiller,L.M.(2017)DNA damage is a pervasive cause of sequencing errors,directly confounding variant identification.Science,355,752-756. 7.Costello,M.,Pugh,T.J.,Fennell,T.J.,Stewart,C.,Lichtenstein,L.,Meldrim,J.C.,Fostel,J.L.,Friedrich,D.C.,Perrin,D.,Dionne,D.,et al.(2013)Discovery and characterization of artifactual mutations in deep coverage targeted capture sequencing data due to oxidative DNA damage during sample preparation.Nucleic Acids Res.,41,e67. 8.Nik-Zainal,S.,Alexandrov,L.B.,Wedge,D.C.,Van Loo,P.,Greenman,C.D.,Raine,K.,Jones,D.,Hinton,J.,Marshall,J.,Stebbings,L.A.,et al.(2012)Mutational processes molding the genomes of 21 breast cancers.Cell,149,979-993. 9.Martincorena,I.,Fowler,J.C.,Wabik,A.,Lawson,A.R.J.,Abascal,F.,Hall,M.W.J.,Cagan,A.,Murai,K.,Mahbubani,K.,Stratton,M.R.,et al.(2018)Somatic mutant clones colonize the human esophagus with age.Science,362,911-917. 10.Tubbs,A.and Nussenzweig,A.(2017)Endogenous DNA Damage as a Source of Genomic Instability in Cancer.Cell,168,644-656. 11.Dong,X.,Zhang,L.,Milholland,B.,Lee,M.,Maslov,A.Y.,Wang,T.and Vijg,J.(2017)Accurate identification of single-nucleotide variants in whole-genome-amplified single cells.Nat.Methods,14,491-493. 12.Zafar,H.,Wang,Y.,Nakhleh,L.,Navin,N.and Chen,K.(2016)Monovar:single-nucleotide variant detection in single cells.Nat.Methods,13,505-507. 13.Chen,C.,Xing,D.,Tan,L.,Li,H.,Zhou,G.,Huang,L.and Xie,X.S.(2017)Single-cell whole-genome analyses by Linear Amplification via Transposon Insertion(LIANTI).Science,356,189-194. 14.Krueger F.(2016)Trim Galore! 15.Langmead,B.and Salzberg,S.L.(2012)Fast gapped-read alignment with Bowtie 2.Nat Meth,9,357-359. 16.McKenna,A.,Hanna,M.,Banks,E.,Sivachenko,A.,Cibulskis,K.,Kernytsky,A.,Garimella,K.,Altshuler,D.,Gabriel,S.,Daly,M.,et al.(2010)The Genome Analysis Toolkit:a MapReduce framework for analyzing next-generation DNA sequencing data.Genome Res.,20,1297-1303. 17.Kim,S.,Scheffler,K.,Halpern,A.L.,Bekritsky,M.A.,Noh,E.,Kallberg,M.,Chen,X.,Kim,Y.,Beyter,D.,Krusche,P.,et al.(2018)Strelka2:fast and accurate calling of germline and somatic variants.Nat.Methods,15,591-594. 18.Hosokawa,M.,Nishikawa,Y.,Kogawa,M.and Takeyama,H.(2017)Massively parallel whole genome amplification for single-cell sequencing using droplet microfluidics.Sci.Rep.,7,5199. 19.Peters,B.A.,Kermani,B.G.,Sparks,A.B.,Alferov,O.,Hong,P.,Alexeev,A.,Jiang,Y.,Dahl,F.,Tang,Y.T.,Haas,J.,et al.(2012)Accurate whole-genome sequencing and haplotyping from 10 to 20 human cells.Nature,487,190-195. 20.Picard Tools(2018). 21.Rimmer,A.,Phan,H.,Mathieson,I.,Iqbal,Z.,Twigg,S.R.F.,Wilkie,A.O.M.,McVean,G.and Lunter,G.(2014)Integrating mapping-,assembly- and haplotype-based approaches for calling variants in clinical sequencing applications.Nat.Genet.,46,912-918. 22.Derrien,T.,Estelle,J.,Marco Sola,S.,Knowles,D.G.,Raineri,E.,Guigo,R.and Ribeca,P.(2012)Fast computation and applications of genome mappability.PLoS One,7,e30377-e30377. 23.Danecek,P.,Auton,A.,Abecasis,G.,Albers,C.A.,Banks,E.,DePristo,M.A.,Handsaker,R.E.,Lunter,G.,Marth,G.T.,Sherry,S.T.,et al.(2011)The variant call format and VCFtools.Bioinformatics,27,2156-2158. 24.Chollet,F.and others(2015)Keras. 25.Kingma,D.P.and Ba,J.(2014)Adam:A Method for Stochastic Optimization.CoRR,abs / 1412.6. 26.Gal,Y.and Ghahramani,Z.(2015)Dropout as a Bayesian Approximation:Representing Model Uncertainty in Deep Learning.arXiv e-prints. 27.Drmanac,R.,Sparks,A.B.,Callow,M.J.,Halpern,A.L.,Burns,N.L.,Kermani,B.G.,Carnevali,P.,Nazarenko,I.,Nilsen,G.B.,Yeung,G.,et al.(2010)Human Genome Sequencing Using Unchained Base Reads on Self-Assembling DNA Nanoarrays.Science(80-.).,327. 28.Arbeithuber,B.,Makova,K.D.and Tiemann-Boege,I.(2016)Artifactual mutations resulting from DNA lesions limit detection levels in ultrasensitive sequencing applications.DNA Res.,23,547-559. 29.Amini,S.,Pushkarev,D.,Christiansen,L.,Kostem,E.,Royce,T.,Turk,C.,Pignatelli,N.,Adey,A.,Kitzman,J.O.,Vijayan,K.,et al.(2014)Haplotype-resolved whole-genome sequencing by contiguity-preserving transposition and combinatorial indexing.Nat.Genet.,46,1343-1349. 30.Zheng,G.X.Y.,Lau,B.T.,Schnall-Levin,M.,Jarosz,M.,Bell,J.M.,Hindson,C.M.,Kyriazopoulou-Panagiotopoulou,S.,Masquelier,D.A.,Merrill,L.,Terry,J.M.,et al.(2016)Haplotyping germline and cancer genomes with high-throughput linked-read sequencing.Nat.Biotechnol.,34,303-311. 31.Northcutt,C.G.,Wu,T.and Chuang,I.L.(2017)Learning with Confident Examples:Rank Pruning for Robust Classification with Noisy Labels.In Proceedings of the Thirty-Third Conference on Uncertainty in Artificial Intelligence,UAI’17.AUAI Press. 32.Natarajan,N.,Dhillon,I.S.,Ravikumar,P.K.and Tewari,A.(2013)Learning with noisy labels.In Advances in neural information processing systems.pp.1196-1204. 33.Hellner,K.,Miranda,F.,Fotso Chedom,D.,Herrero-Gonzalez,S.,Hayden,D.M.,Tearle,R.,Artibani,M.,KaramiNejadRanjbar,M.,Williams,R.,Gaitskell,K.,et al.(2016)Premalignant SOX2 overexpression in the fallopian tubes of ovarian cancer patients:Discovery and validation studies.EBioMedicine,10,137-149. 34. Laks, E., Zahn, H., Lai, D., McPherson, A., Steif, A., Brimhall, J., Biele, J., Wang, B., Masud, T., Grewal, D., et al. (2018) Resource: Scalable whole genome sequencing of 40,000 single cells identifies stochastic aneuploidies, genome replication states and clonal repertoires.bioRxiv,10.1101 / 411058. 35. Moore, L., Leongamornlert, D., Coorens, THH, Sanders, MA, Ellis, P., Dawson, K., Maura, F., Nangalia, J., Tarpey, PS, Brunner, SF, et al. (2018) The mutational landscape of normal human endometrial epithelium. bioRxiv, 10.1101 / 505685. 36. Lee-Six, H., Ellis, P., Osborne, RJ, Sanders, MA, Moore, L., Georgakopoulos, N., Torrente, F., Noorani, A., Goddard, M., Robinson, P., et al. (2018) The landscape of somatic mutation in normal colorectal epithelial cells.bioRxiv,10.1101 / 416800. 37. Burgess, DJ (2019) Spatial transcriptomics coming of age. Nat. Rev. Genet., 20, 317. All references mentioned herein are incorporated herein by reference.

[0126] DigiPico2 Example 2 - DigiPico2, a novel method for determining whole genome sequences of picogram amounts of DNA with unprecedented accuracy

[0127] Introduction We previously described the DigiPico library preparation pipeline and MutLX analysis platform as methods for the accurate identification of single-nucleotide variants (SNVs) from limited amounts of clinical material. This was an important methodological advancement, primarily because the limited amount of genetic material obtained from clinical samples must be whole-genome amplified (WGA) before sequencing. However, the WGA process introduces up to 100,000 artifactual mutations into the amplified DNA, which floods the final analysis with false-positive variant calls that prevent meaningful genetic interpretation from the entire original sample. In the DigiPico / MutLX strategy, we overcome this obstacle by separating individual DNA molecules into independent compartments before the WGA step and indexing them after processing. By doing so, we digitize information for actual mutations, meaning whether each compartment carries a mutated allele or not. However, artifactual mutations, due to the way they are generated during the WGA process, result in compartments that contain both mutated and reference allele information (Figure 1B). Based on this information, we then developed an artificial neural network (ANN)-based algorithm, MutLX, to efficiently identify and remove these artifactual mutations from our data (Figure 2). We extensively tested our strategy on mock data as well as patient samples and demonstrated that our method is indeed effective in removing FP variant calls.

[0128] However, while generating high-quality data, the DigiPico library preparation method suffers from some technical limitations. First, the fragmentation step (CoREF) of the library preparation, which was borrowed from a previously described method, is extremely complex and time-consuming. Furthermore, CoREF requires the use of dUTP during the WGA process. Because dUTP is a non-natural nucleotide, it may introduce additional artifactual mutations into the final product. Second, we found that the adapter ligation efficiency in DigiPico is very low, which can sometimes impair library quality. Finally, due to the large number of indexes and the lack of redundancy in the index information, there is a possibility of index cross-contamination, which can adversely affect the final results. Therefore, we developed the DigiPico2 library preparation method to address all these issues.

[0129] result Improved DigiPico library prep workflow As previously described, the use of dUTP in the DigiPico method is due to the requirement for the CoREF fragmentation procedure, which is an extremely complex fragmentation method (Figure 14A). Therefore, an alternative fragmentation method could simultaneously address both the complexity and the dUTP issue. For this reason, we decided to use the existing fragmentation and end-repair strategy provided by the Lotus DNA Library Preparation Kit (IDT, USA). The Lotus DNA Library Preparation Kit uses an enzyme cocktail to fragment large DNA molecules into smaller pieces in a time-dependent manner and prepare the fragment ends for the ligation step. To make this new strategy compatible with DigiPico, we first used an I-DOT (Dispendix, Germany) dispenser to ensure that all compartments received the enzyme cocktail at approximately the same time. Next, we optimized the reaction conditions to obtain the desired length for DigiPico sequencing. By doing so, we were able to reduce the library preparation time from the original 12 hours to 4.5 hours and eliminated the need for the use of dUTP in the WGA reaction (Figure 14B).

[0130] Next, we aimed to address the low ligation efficiency in DigiPico. Originally, our adapter ligation and indexing relied on the use of an asymmetric ligation approach, in which long indexing oligos with short complementary regions were used for ligation. This was extremely inefficient (Figure 15A). A more efficient approach would require looped common adapters ligated to both ends of the fragments. However, because these adapters do not contain any indexes in their standard form, the products would need to be purified individually, and then the indexes would be introduced later via PCR reactions using index primers. However, this added another challenge, as the purification of 384 independent products would be extremely complex, time-consuming, and pose a significant risk of cross-contamination. To overcome these issues, we designed a novel indexing strategy for DigiPico2. In DigiPico2, a set of indexes is first introduced into the stem-loop of the common adapter (Figure 15B). During the ligation step, all wells in each column of the plate will receive a differently indexed common loop adapter; therefore, a total of 24 different oligos will be sufficient to index all columns of the plate with a first set of indices (column indexing). After the ligation step, all wells in each row are pooled in separate tubes, resulting in 16 different pools. These 16 different pools can typically be purified for use in the next indexing step. In the next step, the purified products from each pool are indexed by PCR using 16 different index primers (row indexing). At the end of this step, each well of the plate will have received a different column-row index combination. This indexing strategy not only greatly improves ligation efficiency, but also introduces two sets of redundancy that can be used to remove index cross-contamination from the data. In the first set, column indexes are attached to both ends of each fragment.Therefore, any cross-contamination will most likely result in fragments with non-identical indices at their ends and can therefore be easily removed from the data. In the second set, the indexing oligos for each row can be doubly indexed so that both the standard index 1 (i7) and standard index 2 (i5) sequences uniquely identify a specific row. By combining the redundancy of these sets, the index cross-contamination rate can be reduced by at least two orders of magnitude.

[0131] The DigiPico2 workflow significantly improves library quality To test the effect of these improvements on the final data quality, we performed DigiPico2 sequencing on 120 pg of DNA from the blood of patient 11152. This sample was chosen because we previously extensively sequenced tumor and normal cells from this patient. As expected, WGA, similar to the previous version, produced a very uniform product distribution (Figure 16A). However, after library preparation and sequencing, unlike DigiPico, DigiPico2 revealed that the representation of each well in the final library appeared to strongly correlate with the amount of WGA product (Figure 16A-C). This is likely a result of improved ligation efficiency. This improved correlation also allows for the introduction of a QC measure based solely on the uniformity of WGA products, which was not previously possible. Furthermore, by analyzing index redundancy information, we found that nearly 5% of reads were removed using the index cross-contamination filter. Without this novel filter, contamination would have adversely affected the analytical results. Finally, we analyzed the data from DigiPico2 using the MutLX algorithm. The final results showed a clearer distinction between artifactual and real mutations, indicating that DigiPico2 performs at least as well as the DigiPico method when analyzed using the MutLX algorithm (Figure 16D).

[0132] Extending the DigiPico2 workflow to single-cell whole-genome sequencing Having established the DigiPico2 workflow, we tested its applicability to single-cell whole-genome sequencing. This is important because active mutational processes can begin within individual cells. Therefore, we implemented a workflow for single-cell DigiPico (ScDigiPico) sequencing by partitioning DNA from individual cells into entire rows of a 384-well plate (Figure 17A). To evaluate the effectiveness of the ScDigiPico method in identifying active mutational processes, we mimicked such processes by mutagenizing cultured Kuramochi cells with N-ethyl-N-nitrosourea (ENU). ENU is an alkylating agent and an extremely potent mutagen, selectively inducing T>C, T>A, and C>T mutations. In this setting, each cell acquires a different set of mutations; however, because the underlying mutagenesis mechanism is the same, the mutations are expected to be of the same type. ScDigiPico enabled the identification of enrichment of the above-mentioned mutation types (Figure 17B). We also further explored the potential of ScDigiPico by irradiating Kuramochi cells with UV light prior to single-cell sorting and ScDigiPico library preparation (Figure ​(Figure17C),17C), and identifying kataegis events within a single cell. These results collectively demonstrate that ScDigiPico is an effective strategy for the precise identification of true mutations caused by active mutational processes in individual cancer cells.

[0133] DigiPico2 protocol 200 pg of purified DNA, 20–30 resuspended nuclei, or laser-capture microdissected tumor islets were first denatured using 5 μl of D2 buffer from the Repli-g single cell kit (Qiagen). After a 5-minute incubation at room temperature, 95 μl of water was added to the sample, and then 200 μl of denatured template was added to each well of a 384-well reaction plate containing 800 nl of WGA mixture (0.58 μl of Sc reaction buffer, 0.04 μl of Sc polymerase (Repli-g single cell kit, Qiagen), 0.04 μl of Evagreen 20x (Biotium), and 0.065 μl of water) using a Mosquito HTS liquid handler (TTP Labtech). The plate was incubated at 30°C for 1.5 hours, followed by heat inactivation at 65°C for 15 minutes. The addition of Evagreen during the reaction allows for monitoring of the WGA reaction using a real-time PCR instrument, if desired. Next, 250 nl of the WGA reaction was transferred to a new plate, and 1.1 μl of NEBNext UltraII FS reaction mix (753 nl water, 270 nl UltraII FS reaction buffer, and 77 nl UltraII FS enzyme mix) was added to each well using an I-DOT dispenser (Dispendix). The plate was incubated at 30°C for 6 hours, followed by a 30-minute incubation at 65°C. Next, 150 nl of DigiPico indexed loop adapters with column indexes were added to all wells. Note that all wells within the same column receive the same indexing oligo at this stage. Next, 1.2 μl of UitraII ligation mixture (1150 nl UitraII ligation master mix, 38 nl ligation enhancer, and 12 nl water) was added to each well with 5 cycles of mixing using a Mosquito liquid handler, and the plate was incubated at 20° C. for 15 minutes, followed by heat inactivation at 65° C. for 10 minutes. All wells in the same row were then pooled together using a Mosquito liquid handler.Then, 1.5 μl of USER Enzyme (NEB) from each row was added to the 20 μl pool, and the reaction was incubated at 37°C for 15 minutes. USER Enzyme cleaves the loop adapters at uracil positions. The products were then size-selected using SPRI beads to obtain a size range of 300-400 bp. The products from each row were then amplified for four cycles using the row index primers. The final products were pooled together, and the final library was purified using SPRI beads.

[0134] ScDigiPico protocol Individual cells in the first column of a 384-well plate were sorted into wells. Each well contained 4.5 μl of MyPK buffer (). The plate was incubated at 55°C for 30 minutes. Next, 900 nl of stop solution was added to each well, and the plate was held at 95°C for 5 minutes to inactivate the proteinase K. The lysed cells were then distributed across the entire row, 200 nl per well, using a Mosquito liquid handler. Next, 800 nl of WGA reaction mixture was added to each well, and library preparation was performed as with DigiPico2.

[0135] Indexed Loop Adapter - Column Index Some portions of the sequence are derived from the library preparation instructions (NEBNext® Multiplex Oligos for Illumin® (Index Primer Set 1) - https: / / international.neb.com / - / media / nebus / files / manuals / manuale7335.pdf?rev=4bf1622b342b4d73a2b01443068ed2c5&hash=B049D91A18CDB471AB388DC6E67E06B79263E5C5). P-[index']AGATCGGAAGAGCACACGTCTGAACTUCCCTACACGACGCTCTTCCGATCT[index] * T (SEQ ID NO: 1) where P is the 5' phosphate group and* indicates a phosphorothioate linkage. Index = column index (Ci) or row index (Ri) sequence, which acts as a unique barcode for each column or row, respectively.

[0136] Row index primers (oligonucleotide sequences copyright 2007-2013 Illumina, Inc. All rights reserved) P5:AATGATACGGCGACCACCGAGATCTACAC[r-index]ACACTCTTTCCCTACACGACGCTCTTCCGATC * T (SEQ ID NO: 2) P7:CAAGCAGAAGACGGCATACGAGAT[r-index]GTGACTGGAGTTCAGACGTGTGCTCTTCCGATC * T (SEQ ID NO: 3) * indicates a phosphorothioate bond.

Claims

1. 1. A method for whole genome sequencing of a single cell or a group of cells for identifying single base variants, chromosomal structural changes, or phasing information in the genome of the single cell or a group of cells, comprising: i) providing a multiwell array plate containing rows and columns of reaction wells; ii) providing genomic DNA of a single cell or group of cells, said genomic DNA being distributed into a plurality of reaction wells of said multi-well array plate, so that there is one single-stranded genomic DNA molecule of any given locus per reaction well; iii) performing whole genome amplification (WGA) of each genomic DNA molecule in each reaction well to obtain multiple copies of said genomic DNA molecule; iv) fragmenting copies of the genomic DNA molecules in each reaction well to form DNA fragments and ligating a pair of loop adapters at each end or tagging with transposase delivery adapters to form adapted DNA fragments, wherein the loop adapters or transposase delivery adapters comprise a column index (Ci) sequence or a row index (Ri) sequence, and the Ci sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a column of the multiwell array plate, or each Ri sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a row of the multiwell array plate; v) obtaining an index DNA library by performing index PCR on the adapted DNA fragments, wherein the adapted DNA fragments are amplified using forward and reverse index primers to form index PCR products, and the forward and reverse index primers introduce row index (Ri) sequences or column index (Ci) sequences at each end of the adapted DNA fragments, respectively, so that the resulting index PCR products contain both adjacent column index (Ci) sequence pairs common to each well in a column and adjacent row index (Ri) sequence pairs common to each well in a row; and vi) sequencing the index DNA library to obtain data for identifying any single base variants, chromosomal structural changes, or phasing information in the genome of the single cell or group of cells; A method comprising:

2. 2. The method of claim 1, wherein the cell or cell group is derived from a tissue biopsy of a subject.

3. 3. The method of claim 1 or claim 2, wherein the cell or cell group comprises a cancer cell, a precancerous cell, or a suspected cancer cell, or a combination of these cells.

4. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 3, wherein the genomic DNA comprises DNA from approximately 1 to 30 cells.

5. 5. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 4, wherein the DNA amount of a single cell is distributed among wells in a single row; or the DNA amount of the cell or group of cells is distributed among both rows and columns of a single multiwell array plate.

6. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 5, wherein the multiwell array plate comprises a 384-well plate.

7. A method for determining the whole genome base sequence of a single cell or a group of cells described in any one of claims 1 to 6, wherein the DNA polymerization reporter molecule is provided in an amplification mixture.

8. a loop adaptor is provided, whereby the method comprises the step of fragmenting copies of the genomic DNA molecule in each reaction well to form DNA fragments, followed by a ligation reaction to ligate loop adaptors to the DNA fragments; or 8. The method for whole genome sequencing of a single cell or group of cells according to any one of claims 1 to 7, wherein the transposase delivery adaptors are provided such that the tagging step comprises fragmenting copies of the genomic DNA molecule.

9. 9. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 8, wherein fragmenting the DNA molecules in each reaction well into a plurality of dsDNA fragments comprises direct enzymatic fragmentation.

10. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 9, wherein the fragmentation or tagging reagents are added to each well simultaneously.

11. 11. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 10, wherein after fragmenting the genomic DNA molecule to form DNA fragments, the DNA fragments are end-repaired and dA-tailed, so that they can be ligated to the loop adapters.

12. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 11, wherein the loop adaptor comprises an oligonucleotide having a secondary stem-loop structure, the loop adaptor comprises a pair of complementary sequence regions adjacent to the loop region, and the pair of complementary sequence regions are arranged so as to hybridize with each other to form the stem-loop structure of the loop adaptor.

13. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 12, wherein both ends of the adapted DNA fragment are symmetrical.

14. 14. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 13, wherein the loop adaptor contains uracil in the loop region and after ligation of the loop adaptor to form loop DNA, and the loop DNA is cleaved at the position of the uracil.

15. Ci sequences are provided in the adapted DNA fragments, and the method may additionally comprise the step of pooling the adapted DNA fragments of each reaction well in a row prior to the index PCR; or 15. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 14, wherein Ri sequences are provided in the adapted DNA fragments, and the method further comprises the step of pooling the adapted DNA fragments of each reaction well in a column prior to the index PCR.

16. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 15, wherein the adapted DNA fragment contains a Ci sequence, and the forward and reverse index PCR primers each contain an Ri sequence to obtain a pair of Ri sequences in the index PCR product.

17. 17. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 16, wherein the forward and reverse index PCR primers further comprise sequencing adapter sequences, whereby sequencing adapters are incorporated into the index PCR products.

18. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 17, wherein the index DNA fragment sizes of the index DNA library are filtered.

19. the method comprising identifying any actual SNVs in the genome of the single cell or group of cells by determining whether substantially all index DNA library sequences from a single well contain the same SNV, or whether only some of the index DNA library sequences contain the same SNV; 19. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 18, wherein SNVs that appear in substantially all of the index DNA library sequences derived from a single well are identified as actual SNVs in the genomic DNA, and SNVs found in only some of the index DNA library sequences derived from a single well are identified as false positive (FP) SNVs.

20. the method further comprising matching an index DNA library sequence from a single well representing a first strand of the genomic DNA with an index DNA library sequence from another well representing a strand of genomic DNA complementary to the first strand; The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 19, wherein an SNV that is substantially present in all index DNA library sequences of both complementary strands of the genomic DNA is identified as an actual SNV, and an SNV that is not substantially present in all index DNA library sequences of both complementary strands of the genomic DNA is identified as a false positive.

21. 21. The method for whole genome sequencing of a single cell or a group of cells according to claim 19 or 20, wherein the identification is performed in silico using BAM file data generated from mapping sequence data to a reference genome.

22. 22. The method for whole genome sequencing of a single cell or a group of cells according to claim 21, wherein the in silico identification or matching of index DNA sequences and / or calculation of probability scores is performed by an artificial neural network (ANN) model.

23. 23. The method for whole genome sequencing of a single cell or a group of cells according to any one of claims 1 to 22, wherein the method comprises preparing index DNA libraries from both tumor cells, suspected tumor cells, or pre-cancerous cells and normal (i.e., non-cancerous) cells, and comparing the index DNA library from the tumor cells, suspected tumor cells, or pre-cancerous cells with the index DNA library obtained from the normal (i.e., non-cancerous) cells taken as a control from normal tissue.

24. 24. The method of whole genome sequencing of a single cell or group of cells according to any one of claims 1 to 23, wherein a probability score that a particular nucleotide variant is an actual SNV or a false positive is calculated in silico, whereby a given variant nucleotide is identified as having a statistically significant probability of being an actual SNV or a false positive.

25. The method for determining the whole genome base sequence of a single cell or a group of cells according to any one of claims 1 to 24, wherein the index DNA library is provided in the form of a paired-read FastQ file.

26. 26. The method of whole genome sequencing of a single cell or cell population of claim 25, wherein the sequence data of the paired read FastQ file is trimmed for adapter sequence removal and for quality to provide trimmed data.

27. 1. A method for preparing an index DNA library for sequencing a nucleic acid molecule, comprising: i) providing a multiwell array plate containing rows and columns of reaction wells; ii) providing nucleic acid molecules, said nucleic acid molecules being distributed into a plurality of reaction wells of said multiwell array plate, whereby there is one single-stranded nucleic acid molecule from any given locus per reaction well; iii) performing amplification of said nucleic acid molecule in each reaction well to obtain multiple DNA molecule copies of said nucleic acid molecule; iv) fragmenting copies of the DNA molecules in each reaction well and ligating a pair of loop adapters at each end or tagging them with transposase delivery adapters to form adapted DNA fragments, wherein the loop adapters or transposase delivery adapters comprise a column index (Ci) sequence or a row index (Ri) sequence, and the Ci sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a column of the multiwell array plate, or each Ri sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a row of the multiwell array plate; v) obtaining the index DNA library by performing index PCR on the adapted DNA fragments, wherein the adapted DNA fragments are amplified to form index PCR products using forward and reverse index primers, which respectively introduce row index (Ri) or column index (Ci) sequences at each end of the adapted DNA fragments, such that the resulting index PCR products contain both a pair of adjacent column index (Ci) sequences common to each well in a column and a pair of adjacent row index (Ri) sequences common to each well in a row; A method comprising:

28. 1. A method for preparing an index DNA library for whole genome sequencing of a single cell or a group of cells for identifying single base variants, chromosomal structural changes, or phasing information in the genome of the single cell or a group of cells, comprising: i) providing a multiwell array plate containing rows and columns of reaction wells; ii) providing genomic DNA of a single cell or group of cells, said genomic DNA being distributed into a plurality of reaction wells of said multi-well array plate, so that there is one single-stranded genomic DNA molecule of any given locus per reaction well; iii) performing whole genome amplification (WGA) of each genomic DNA molecule in each reaction well to obtain multiple copies of said genomic DNA molecule; iv) fragmenting copies of the genomic DNA molecules in each reaction well and ligating a pair of loop adapters at each end or tagging them with transposase delivery adapters to form adapted DNA fragments, wherein the loop adapters or transposase delivery adapters comprise a column index (Ci) sequence or a row index (Ri) sequence, and the Ci sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a column of the multiwell array plate, or each Ri sequence is common to each loop adapter or transposase delivery adapter of all reaction wells in a row of the multiwell array plate; v) obtaining the index DNA library by performing index PCR on the adapted DNA fragments, wherein the adapted DNA fragments are amplified to form index PCR products using forward and reverse index primers, which respectively introduce row index (Ri) or column index (Ci) sequences at each end of the adapted DNA fragments, such that the resulting index PCR products contain both a pair of adjacent column index (Ci) sequences common to each well in a column and a pair of adjacent row index (Ri) sequences common to each well in a row; A method comprising:

29. 1. A method for whole genome sequencing of a single cell or a group of cells to provide data for identifying single nucleotide variants (SNVs), identifying chromosomal structural changes, or identifying phasing information in the genome of said single cell or group of cells, comprising: i) preparing an index DNA library by carrying out the method of claim 27 or 28 or by providing an index DNA library prepared according to claim 27 or 28; ii) sequencing the index DNA library to obtain data for identifying all single nucleotide variants (SNVs), chromosomal structural changes, or phasing information in the genome of the single cell or group of cells; A method comprising:

Images