Systems and methods for clustering sequencing reads to identify PCR duplicates
Inexact matching of sorting features in sequencing reads addresses the inefficiencies of exact matching methods, improving PCR duplicate detection and data analysis accuracy in high error rate sequencing technologies.
Patent Information
- Application Number
- PCT/US2025/044150
- Authority / Receiving Office
- WO · WO
- Patent Type
- Applications
- Current Assignee / Owner
- Priority Date
- 2024-08-30
- Filing Date
- 2025-08-29
- Publication Date
- 2026-03-05
AI Technical Summary
Existing methods for identifying PCR duplicates in sequencing data, particularly for high error rate technologies like single-molecule sequencing, are inefficient due to reliance on exact matching, leading to undercounting of duplicates and loss of accuracy.
A method for clustering sequencing reads using inexact matching of sorting features, such as alignment positions and unique molecular identifiers (UMIs), allowing for more effective identification of PCR duplicates by tolerating a predetermined distance in these features.
Improves the detection of PCR duplicates in high error rate sequencing technologies, enhancing the accuracy and efficiency of sequencing data analysis.
Smart Images

Figure US2025044150_05032026_PF_FP_ABST
Abstract
Description
PATENT Client Reference No.: P39621-WO-1 INTERNATIONAL PATENT APPLICATION Title: SYSTEMS AND METHODS FOR CLUSTERING SEQUENCING READS TO IDENTIFY PCR DUPLICATES Inventors: Marghoob Mohiyuddin, a U.S. citizen, resident of Milpitas, CA Yunfei Guo, a citizen of China, resident of Newton, MA Assignee: Roche Sequencing Solutions, Inc. 4300 Hacienda Drive Pleasanton, CA 94588 United States of America Entity: LargePATENT Client Reference No.: P39621-WO-1 SYSTEMS AND METHODS FOR CLUSTERING SEQUENCING READS TO IDENTIFY PCR DUPLICATES CROSS-REFERENCE TO RELATED APPLICATIONS
[0001] This Application claims priority to United States Provisional Patent Application No.63 / 689,556, filed on August 30, 2024, which is herein incorporated by reference in its entirety. BACKGROUND
[0002] The development of affordable and rapid DNA sequencing technologies has enabled the development of targeted therapeutics that rely on the use of DNA biomarkers to identify patients that are suitable for receiving the targeted therapy. For example, mutations in certain genes, such as genes involved in cell proliferation, are known to lead to certain types of cancers that can be treated very effectively with specific types of drugs. Other mutations are known to confer resistance to certain therapies. Therefore, there is a need for improved systems and methods to accurately read molecules of interest and identify variants from sequencing data. SUMMARY
[0003] The present invention relates generally to systems and methods of clustering reads from sequencing data, and more specifically identifying PCR duplicate reads from the sequencing data based on inexact matching of sorting features.
[0004] In accordance with a first aspect of the present disclosure, a method is provided for clustering genetic sequence reads from a sample. The method includes: sequencing the sample to generate a plurality of unclustered sequence reads, each unclustered sequence read associated with a set of sorting features; sorting the plurality of unclustered sequence reads according to a primary sorting feature; for each subset of unclustered sequence reads with identical values of the primary sorting feature: selecting a first unclustered sequence read of the subset as an anchor read for a cluster of sequence reads, adding to the cluster any unclustered sequence reads that have values of one or more secondary sorting features other than the primary sorting feature that fall within a predetermined tolerance distance of the one or more secondary sorting features of the anchor read, and, in response to determining that one or more sequence reads in the subset remain unclustered, repeating the steps of selectingPATENT Client Reference No.: P39621-WO-1 and adding to generate one or more additional clusters until each of the plurality of unclustered sequence reads have been added to a particular cluster; and generating an output that indicates a set of clusters.
[0005] In at least some embodiments of the first aspect, the set of sorting features comprises at least two of the following: a unique molecular identifier (UMI), an alignment starting position, an alignment ending position, a number of 5' and 3' soft-clipped base pairs, a number of supplementary alignments, an edit distance to a reference sequence, a fragment length, and a mapping quality.
[0006] In at least some embodiments of the first aspect, the set of sorting features comprise the UMI, the alignment starting position, and the alignment ending position. In an embodiment, the primary sorting feature is the UMI, and the one or more secondary sorting features comprises the alignment starting position and the alignment ending position. In yet another embodiment, the primary sorting feature is the alignment starting position, and the one or more secondary sorting features comprise the alignment ending position and the UMI.
[0007] In at least some embodiments of the first aspect, the method further includes, for each cluster: determining a maximum variation for the values of the one or more secondary sorting features for the sequence reads in the given cluster; comparing the maximum variation to the predetermined tolerance distance; and, if the maximum variation is less than the predetermined tolerance distance, adding to the given cluster any additional unclustered sequence reads that have values of the one or more secondary sorting features that fall within an additional tolerance distance from the values of the one or more secondary sorting feature of the anchor read.
[0008] In at least some embodiments of the first aspect, all sequence reads in a particular cluster are identified as PCR duplicate reads.
[0009] In at least some embodiments of the first aspect, the method further including modifying at least one sequencing determination based on the identification of PCR duplicate reads.
[0010] In at least some embodiments of the first aspect, the predetermined tolerance distance varies according to the secondary sorting feature being evaluated.PATENT Client Reference No.: P39621-WO-1
[0011] In at least some embodiments of the first aspect, the predetermined tolerance distance associated with an alignment starting position is greater than the predetermined tolerance distance associated with an alignment ending position.
[0012] In accordance with a second aspect of the present disclosure, a computer- implemented method is provided for identifying PCR duplicates. The method includes: receiving sequencing data comprising a plurality of sequence reads, where each sequence read is associated with at least a unique molecular identifier (UMI) and a positional coordinate comprising an alignment start position and an alignment end position; identifying an anchor coordinate for a cluster, wherein the anchor coordinate is selected as a positional coordinate of an unclustered read in the plurality of sequence reads; determining at least one candidate read in the plurality of sequence reads to add to the cluster, wherein a distance metric for each candidate read is less than a threshold distance from the anchor coordinate; in response to determining that at least one sequence read in the plurality of sequence reads remains unclustered, identifying another anchor coordinate for an additional cluster and determine at least one candidate read in the plurality of sequence reads to add to the additional cluster; repeat the identifying and determining steps until all sequence reads in the plurality of sequence reads are added to a corresponding cluster; store an output identifying the clusters of sequence reads in a memory; and identify PCR duplicates in the plurality of sequence reads based on the output.
[0013] In at least some embodiments of the second aspect, the method further includes identifying a UMI for the cluster. The determining at least one candidate read further comprises determining that the UMI for the candidate read matches the UMI for the cluster.
[0014] In at least some embodiments of the second aspect, the method further includes identifying a UMI for the cluster. The determining at least one candidate read further comprises determining whether the UMI for the at least one candidate read is within a Hamming distance of the UMI for the cluster. The at least one candidate read is added to the cluster in response to determining that the UMI for the candidate read is within the Hamming distance of the UMI for the cluster; or the candidate read is added to the cluster in response to determining that the UMI for the candidate read is not within the Hamming distance of the UMI for the cluster, but is within a Levenshtein distance of the UMI for the cluster.PATENT Client Reference No.: P39621-WO-1
[0015] In at least some embodiments of the second aspect, the at least one candidate read is not added to the cluster in response to determining that the UMI for the at least one candidate read is not within the Hamming distance or a Levenshtein distance of the UMI for the cluster.
[0016] In at least some embodiments of the second aspect, the distance metric is calculated by determining whether an absolute difference between the alignment starting position of the given candidate read and the alignment starting position of the anchor coordinate is less than a threshold distance, and determining whether an absolute difference between the alignment end position of the given candidate read and the alignment end position of the anchor coordinate is less than the threshold distance.
[0017] In at least some embodiments of the second aspect, the distance metric is calculated by determining whether a Euclidean distance between the position coordinate for the candidate read and the anchor coordinate is less than a threshold distance.
[0018] In at least some embodiments of the second aspect, identifying the PCR duplicates comprises identifying all sequence reads associated with a given cluster as PCR duplicates.
[0019] In accordance with a third aspect of the present disclosure, a system is provided for identifying PCR duplicates. The system includes at least one processor configured to: receive a plurality of unclustered sequence reads, each unclustered sequence read associated with a set of sorting features; sort the plurality of unclustered sequence reads according to a primary sorting feature; for each subset of unclustered sequence reads with identical values of the primary sorting feature: select a first unclustered sequence read of the subset as an anchor read for a cluster of sequence reads, add to the cluster any unclustered sequence reads that have values of one or more secondary sorting features other than the primary sorting feature that fall within a predetermined tolerance distance of the one or more secondary sorting features of the anchor read, and in response to determining that one or more sequence reads in the subset remain unclustered, repeat the steps of selecting and adding to generate one or more additional clusters until each of the plurality of unclustered sequence reads have been added to a particular cluster; and generate an output that indicates a set of clusters.PATENT Client Reference No.: P39621-WO-1
[0020] In at least some embodiments of the third aspect, the system further includes a sequencing instrument configured to generate the plurality of unclustered sequence reads in accordance with a sequencing by expansion process. BRIEF DESCRIPTION OF THE DRAWINGS
[0021] The detailed description below is set forth with reference to the accompanying figures, which include the following.
[0022] FIG.1 sets forth an illustrative system including a sequencing device 110 communicatively coupled to a computing system, in accordance with at least some embodiments of the present disclosure.
[0023] FIG.2 is an illustration of positional clustering with inexact matching using a square centered on a read in question, in accordance with embodiments of the present disclosure.
[0024] FIG.3 is an illustrative example of a sparse representation used to cluster reads, in accordance with embodiments of the present disclosure.
[0025] FIG.4 is an illustration positional clustering with inexact matching using an anchor read in the bottom left corner of a square representing the cluster, in accordance with embodiments of the present disclosure.
[0026] FIG.5 is an illustration of a rescue step to the positional clustering described in FIG.3, in accordance with embodiments of the present disclosure.
[0027] FIG.6 is an illustrative example of maximal clustering in accordance with embodiments of the present disclosure.
[0028] FIG.7 is a flowchart of a method of clustering to identify PCR duplicates in accordance with embodiments of the present disclosure.
[0029] FIG.8 illustrates an example computer system that may be utilized to implement techniques disclosed herein. DETAILED DESCRIPTIONPATENT Client Reference No.: P39621-WO-1
[0030] Systems and methods for identifying PCR duplicates in sequencing data, such as data generated by a sequencer instrument performing a sequencing assay, are provided herein. In particular, the sequencing data can be generated by a next generation sequencing (NGS) instrument, such as a nanopore-based sequencer instrument optimized for sequencing- by-expansion (SBX) chemistry. Definitions
[0031] Additional aspects and advantages of the present disclosure will become readily apparent to those skilled in this art from the following detailed description, wherein only illustrative embodiments of the present disclosure are shown and described. As will be realized, the present disclosure is capable of other and different embodiments, and its several details are capable of modifications in various obvious respects, all without departing from the disclosure. Accordingly, the drawings and description are to be regarded as illustrative in nature, and not as restrictive.
[0032] All of the U.S. patents, U.S. patent application publications, U.S. patent applications, foreign patents, foreign patent applications and non-patent publications referred to in this specification and / or listed in the Application Data Sheet are incorporated herein by reference, in their entirety. Aspects of the embodiments can be modified, if necessary to employ concepts of the various patents, applications and publications to provide yet further embodiments.
[0033] The section headings used herein are for organizational purposes only and are not to be construed as limiting the subject matter described.
[0034] Although various features of the disclosure can be described in the context of a single embodiment, the features can also be provided separately or in any suitable combination. Conversely, although the present disclosure can be described herein in the context of separate embodiments for clarity, different aspects of the various embodiments described in the present disclosure can also be implemented in a single embodiment. It is to be understood that the present disclosure is not limited to the particular embodiments described herein and as such can vary. Those of skill in the art will recognize that there are variations and modifications of the present disclosure, which are encompassed within its scope.PATENT Client Reference No.: P39621-WO-1
[0035] All terms are intended to be understood as they would be understood by a person skilled in the art. Unless defined otherwise, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which the disclosure pertains.
[0036] The term “about” or “approximately” means within an acceptable error range for the particular value as determined by one of ordinary skill in the art, which will depend in part on how the value is measured or determined, i.e., the limitations of the measurement system. Where particular values are described in the application and claims, unless otherwise stated the term “about” meaning within an acceptable error range for the particular value should be assumed.
[0037] When a feature or element is herein referred to as being “on” another feature or element, it can be directly on the other feature or element or intervening features and / or elements may also be present. In contrast, when a feature or element is referred to as being “directly on” another feature or element, there are no intervening features or elements present. It will also be understood that, when a feature or element is referred to as being “connected," “attached” or “coupled” to another feature or element, it can be directly connected, attached or coupled to the other feature or element or intervening features or elements may be present. In contrast, when a feature or element is referred to as being “directly connected," “directly attached” or “directly coupled” to another feature or element, there are no intervening features or elements present. Although described or shown with respect to one embodiment, the features and elements so described or shown can apply to other embodiments. It will also be appreciated by those of skill in the art that references to a structure or feature that is disposed “adjacent” another feature may have portions that overlap or underlie the adjacent feature.
[0038] Terminology used herein is for the purpose of describing particular embodiments only and is not intended to be limiting of the invention. For example, as used herein, the singular forms “a," “an” and “the” are intended to include the plural forms as well, unless the context clearly indicates otherwise. It will be further understood that the terms “comprises” and / or “comprising,” when used in this specification, specify the presence of stated features, steps, operations, elements, and / or components, but do not preclude the presence or addition of one or more other features, steps, operations, elements, components, and / or groups thereof.PATENT Client Reference No.: P39621-WO-1 As used herein, the term “and / or” includes any and all combinations of one or more of the associated listed items and may be abbreviated as “ / ."
[0039] Spatially relative terms, such as “under," “below," “lower," “over," “upper” and the like, may be used herein for ease of description to describe one element or feature’s relationship to another element(s) or feature(s) as illustrated in the figures. It will be understood that the spatially relative terms are intended to encompass different orientations of the device in use or operation in addition to the orientation depicted in the figures. For example, if a device in the figures is inverted, elements described as “under” or “beneath” other elements or features would then be oriented “over” the other elements or features. Thus, the exemplary term “under” can encompass both an orientation of over and under. The device may be otherwise oriented (rotated 90 degrees or at other orientations) and the spatially relative descriptors used herein interpreted accordingly. Similarly, the terms “upwardly," “downwardly," “vertical," “horizontal” and the like are used herein for the purpose of explanation only unless specifically indicated otherwise.
[0040] Although the terms “first” and “second” may be used herein to describe various features / elements (including steps), these features / elements should not be limited by these terms, unless the context indicates otherwise. These terms may be used to distinguish one feature / element from another feature / element. Thus, a first feature / element discussed below could be termed a second feature / element, and similarly, a second feature / element discussed below could be termed a first feature / element without departing from the teachings of the present invention. Overview
[0041] Prepared nucleic acid molecules of interest (e.g., a sequencing library) can be sequenced using a sequencing assay as part of the procedure for determining sequencing reads for a plurality of loci in the genome of an organism. Any of a number of sequencing technologies or sequencing assays can be utilized. The term "Next Generation Sequencing (NGS)" as used herein refers to sequencing methods that allow for massively parallel sequencing of clonally amplified molecules and of single nucleic acid molecules (or of nucleic acid analogues).PATENT Client Reference No.: P39621-WO-1
[0042] Non-limiting examples of sequence assays that are suitable for use with the methods disclosed herein include nanopore sequencing (US Pat. Publ. Nos.2013 / 0244340, 2013 / 0264207, 2014 / 0134616, 2015 / 0119259 and 2015 / 0337366), Sanger sequencing, capillary array sequencing, thermal cycle sequencing (Sears et al., Biotechniques, 13:626-633 (1992)), solid-phase sequencing (Zimmerman et al., Methods Mol. Cell Biol., 3:39-42 (1992)), sequencing with mass spectrometry such as matrix-assisted laser desorption / ionization time-of-flight mass spectrometry (MALDI-TOF / MS; Fu et al., Nature Biotech., 16:381-384 (1998)), sequencing by hybridization (Drmanac et al., Nature Biotech., 16:54-58 (1998), and NGS methods, including but not limited to sequencing by synthesis (e.g., HiSeq™, MiSeq™, or Genome Analyzer, each available from Illumina), sequencing by ligation (e.g., SOLiD™, Life Technologies), ion semiconductor sequencing (e.g., Ion Torrent™, Life Technologies), and SMRT® sequencing (e.g., Pacific Biosciences).
[0043] Commercially available sequencing technologies include, but are not limited to, the following: sequencing-by-hybridization platforms from Affymetrix, Inc., now ThermoFisher Scientific, (Sunnyvale, Calif.), sequencing-by-synthesis platforms from Illumina, Inc. (San Diego, Calif.) and Helicos Biosciences (Cambridge, Mass.), and sequencing-by-ligation platform from Applied Biosystems (Foster City, Calif.). Other sequencing technologies include, but are not limited to, the Ion Torrent technology from ThermoFisher Scientific (Waltham, MA); as well as nanopore sequencing from Roche Sequencing Solutions, Inc. (Santa Clara, Cal.) and Oxford Nanopore Technologies (Oxford, United Kingdom). Systems and Devices for Identifying PCR Duplicates in Sequencing Data
[0044] FIG.1 sets forth an illustrative system 100 including a sequencing device 110 (which may be alternately referred to as a sequencer instrument) communicatively coupled to a computing system 102. Sequencing device 110 can be coupled to computing system 102 either directly (e.g., through one or more communication cables) or through network 130, which may be the Internet or any other combination of wide-area, local area, wired, and / or wireless networks. In some embodiments, computing system 102 may be included in or integrated with the sequencing device 110. In some embodiments, sequencing device 110 may sequence (e.g., perform a biochemical assay) a sample containing genetic material and produce resulting sequencing data. The sequencing data can be sent to computing system 102 (e.g., through network 130) or stored on a storage device and at a later stage transferred toPATENT Client Reference No.: P39621-WO-1 computing system 102 (e.g., through network 130). In some embodiments, computing system 102 may or may not include a display 108 and one or more input devices (not illustrated) for receiving commands from a user or operator (e.g. a technician or a geneticist). In some embodiments, computing system 102 and / or sequencing device 110 can be accessed by users or other devices remotely through network 130. Thus, in some embodiments various methods discussed herein may be run remotely on computing system 102.
[0045] Computing system 102 may include one computing device or a combination of a number of computing devices of any type, such as personal computers, laptops, network servers (e.g., local servers or servers included on a public / private / hybrid cloud), mobile devices, etc., where some or all of the devices can be interconnected. Computing system 102 may include one or more processors (not illustrated), each of which can have one or more logic cores. In some embodiments, computing system 102 can include one or more general- purpose processors (e.g., CPUs), special-purpose processors such as graphics processors (GPUs), digital signal processors, or any combination of these and other types of processors. In some embodiments, some or all processors in computing system can be implemented using customized or customizable circuitry, such as application specific integrated circuits (ASICs) or field programmable gate arrays (FPGAs). Computing system 102 can also in some embodiments retrieve and execute non-transitory computer-readable instructions stored in one or more memories or storage devices (not illustrated) integrated into or otherwise communicatively coupled to computing system 102. The memory / storage devices can include any combination of non-transitory computer readable storage media including semiconductor memory chips of various types (e.g., dynamic random access memory (DRAM), static random access memory (SRAM), synchronous dynamic random access memory (SDRAM), flash memory, programmable read-only memory, etc.) and so on. Magnetic and / or optical disks can also be used. The memories / storage devices can also include removable storage media that can be readable and / or writeable; examples of such media include compact disc (CD), read-only digital versatile disc (e.g., DVD-ROM, dual-layer DVD-ROM), read-only and recordable Blu-ray® disks, ultra-density optical disks, flash memory cards (e.g., SD cards, mini-SD cards, micro-SD cards, etc.), and so on. In some embodiments, data and other information (e.g. sequencing data) can be stored in one or more remote locations, e.g., cloud storage, and synchronized with other the components of system 100.PATENT Client Reference No.: P39621-WO-1
[0046] In some embodiments, the sequencing device 110 can generate sequencing data by a sequencing-by-expansion (SBX) process. Examples of the SBX process include those described in U.S. Patent Application No.17 / 456,342 (U.S. Publication No. US20220411458A1), entitled “Translocation control elements, reporter codes, and further means for translocation control for use in nanopore sequencing,” filed November 23, 2021, which is herein incorporated by reference in its entirety. During library preparation in the SBX process, a number of surrogate molecules are derived from and characterize nucleic acid material provided in a sample.
[0047] More particularly, the SBX process may translate a sequence of DNA into a measurable surrogate molecule called an Xpandomer. Xpandomer synthesis based on the natural function of DNA replication uses expandable nucleotide triphosphates (X-NTPs) that act as substrates for template-dependent, polymerase-based replication. These Xpandomer molecules are then processed by a sequencer instrument (e.g., the sequencing device 110) to measure the sequence in the original DNA template. As the Xpandomer molecule transits through nanometer-sized openings in an electrode-resistant membrane (a “nanopore”), each nanopore corresponding to a selective channel, a distinct electrical signal is generated for each base reporter and identifiable to enable highly accurate and high throughput nanopore- based nucleic acid sequencing (also referred to generally as nanopore sequencing).
[0048] Sequencing device 110 can generate a plurality of sequence reads corresponding to a genetic sample (e.g., a sample comprising a patient’s DNA or RNA material). For example, a sequence may be identified by processing a sample, which may include (for example) a blood, saliva, or tissue biopsy collected from a subject. Sequence reads can be obtained either directly from sequencing device 110, or from one or more local or remote volatile or non-volatile memories, storage devices, or databases communicatively coupled to computing system 102. Sequence reads can be pre-processed (e.g., pre-aligned) or they can be “raw,” in which case a downstream method may include a preprocessing (e.g., pre- aligning) step. Also, while in some embodiments, entire sequence reads (as generated by sequencing device 110) can be obtained, in other embodiments only sections of the sequence reads can be obtained. Thus, “obtaining a sequence read,” as used herein, refers generally to obtaining one or more sections of one or more (e.g., adjacent) sequence reads. Bioinformatics Workflow OverviewPATENT Client Reference No.: P39621-WO-1
[0049] The output of an NGS sequencer is generally processed by a bioinformatics pipeline that processes the raw signal from the NGS sequencer and translates the raw signal into base calls, often referred to as raw reads, which are typically stored in a FASTQ file that combines the raw reads with associated quality data. This portion of the bioinformatics pipeline is often referred to as primary analysis.
[0050] The next section of the bioinformatics pipeline is called secondary analysis, and it takes the raw reads generated by the primary analysis, and performs several tasks, including alignment and variant calling.
[0051] Tertiary analysis is the final portion of the bioinformatics pipeline and uses the variant calling information to generate medical insights that health care practitioners can use to improve treatments for their patients. For example, information about observed variants in the sample can be used to identify particular cancer diagnoses or cancer risk, or to determine whether targeted therapeutics are appropriate for a particular patient. Secondary Analysis
[0052] New sequencing technologies, such as nanopore-based sequencers, generate sequencing data with different characteristics than sequencing data generated by the current market leading sequencers, such as sequencers manufactured by Illumina, Inc. For example, these differences can include differences in raw read accuracy and differences in the error profiles. Because Illumina sequencers currently dominate the market, the vast majority of the secondary analysis software tools that have been developed are custom tailored to process the type of data that is generated by the Illumina sequencers. These software tools, which typically work very well with data from Illumina sequencers, may not work well with data generated by new next generation sequencing technologies, such as nanopore sequencers. Consequently, there is a need to develop new secondary analysis tools that work well with the new sequencing technologies that are currently being developed. In addition, although the methods described herein may be particularly effective with nanopore sequencing data, the methods can also be used with other types of sequencing data, such as data from an Illumina sequencer. Method for Clustering PCR DuplicatesPATENT Client Reference No.: P39621-WO-1
[0053] Identification of PCR duplicates in sequencing data is a crucial step in sequencing data analysis. In the process of PCR, a single original molecule may yield multiple copies of said molecule, of which some can end up being sequenced. The process of PCR duplicate detection identifies duplicate copies of the same original fragment to get a more representative view of the original set of molecules present in a sample. In some applications, molecular barcodes may also be leveraged in the process of identifying these duplicates. Furthermore, depending on the application, the set of identified PCR duplicate read sequences may be replaced by a consensus sequence which is a more accurate read of the original molecule.
[0054] One key step in identifying PCR duplicates is alignment of the reads to a genomic reference. In most conventional approaches, once alignment is done, two reads are considered duplicates when their alignments start and end at the same exact locations. In cases where unique molecular identifiers (UMIs) are present, in addition to matching starting and end coordinates, UMI sequences must also match exactly for reads to be considered duplicates of one another. While existing approaches (e.g., the “MarkDuplicates” tool available in Picard, and the “GroupReadsByUmi” tool available in fgbio) that use exact matching may be suitable for low error (particularly those low in INDEL errors) sequencing technologies, such approaches would lead to a significant loss in detection of the duplicates for higher throughput and higher error technologies like single-molecule sequencing, since duplicates would be undercounted.
[0055] Accordingly, proposed herein are methods and systems for identifying PCR duplications in sequencing data by using inexact matching to both (1) identify PCR duplicates more effectively for higher error rate sequencing technologies and (2) improve efficiency in clustering potential duplicate reads. While any suitable features of a read may be used to cluster potential duplicate reads, detailed examples are provided herein of using inexact matching of both the read alignment (start and end position) coordinates as well as UMIs / molecular barcodes to identify PCR duplicates. Additional features that could be used for clustering include, but are not limited to, a number of 5' and 3' soft-clipped basepairs, a number of supplementary alignments, an edit distance to a reference sequence, a fragment length, and a mapping quality. In order to identify PCR duplicates, reads are clustered based on a set of sorting features, where the clustering allows for inexact matching in any of the set of sorting features. As will be understood in the art, UMIs used in sequencing workflowsPATENT Client Reference No.: P39621-WO-1 may be taken from a fixed set of designs (referred to hereinafter as “fixed UMIs”) or random sequences (referred to hereinafter as “random UMIs”). Accordingly, alternative related methods are provided herein for both fixed UMIs and random UMIs. In each case, the starting point is to generate unclustered sequencing reads from a sample of interest. As indicated above, any sequencing method may be used, including those described above. The resulting unclustered sequencing reads may be used as the starting point for the methods of clustering described herein. Fixed UMI Read Clustering
[0056] Once the unclustered sequencing reads have been generated, the reads may be sorted in the following order: 1) UMI, 2) alignment start position (“start_pos”), and 3) alignment end position (“end_pos”). In the case of fixed UMIs, only exact UMI matches may be clustered together, thus the reads are sorted by UMI first, and then the alignment start / end positions. This result yields a number of subsets of reads with identical UMIs sorted by both alignment start and end positions therein. Rather than clustering together only reads that share the exact UMI, start, and end position as conventional methodologies do, the proposed methodology clusters together those reads that have identical UMIs but which fall within a predetermined distance with respect to the start and end position.
[0057] In the case of fixed UMIs, the positional clustering within a predetermined distance may be thought of as clustering in two dimensions where a given read can be clustered to all reads within a square having a width centered on an anchor point associated with that read. The width may be a predetermined tolerance distance indicative of a desired tolerance for the inexact position matching in each direction. In other words, the alignment start position and alignment end position are associated with the same predetermined tolerance distance. In some embodiments, the predetermined distance can be different for the alignment start positions and the alignment end positions. In other words, a desired tolerance can be given by a width w for a predetermined tolerance distance indicative of a desired tolerance for the alignment start position and a height h for a predetermined tolerance distance indicative of a desired tolerance for the alignment end position. Thus, reads can be clustered based on whether the alignment coordinates (e.g., <start, end>) for each read fall within a window of dimensions w x h.PATENT Client Reference No.: P39621-WO-1
[0058] FIG.2 is an illustration of positional clustering with inexact matching using a square centered on a read in question, in accordance with embodiments of the present disclosure. In this example, each read may have an identical UMI and thus may be thought of as a point in two dimensions with the coordinates (start_pos, end_pos) given by the alignment information. As can be seen in plot 200 of FIG.2, the centered read 202 may be clustered to other reads within square 204. Square 204 may have a width as determined by the aforementioned predetermined tolerance distance . In other words, read 202 may be clustered together with each of reads 206a-d, all of which have a start and end position within a certain predetermined tolerant distance of the read 202. In contrast, reads 208a-c which fall outside of square 204 would not be clustered together with the read 202. It will be appreciated that the predetermined tolerance distance here is given as half of the width w of square 204 due to the anchor location being located in the center of the square 204.
[0059] For each contiguous block of reads with identical UMIs the following steps maybe taken to cluster the reads positionally. First, the reads are pre-processed to form an indexdata structure. The index data structure may be a sparse representation of the reads which allows efficient storage of the alignment start and end positions represented as coordinates [start_pos, end_pos] while also collapsing (i.e. clustering) the reads with identical start_pos and end_pos. An example of such an index data structure is provided in FIG.3. As can be seen in FIG.3, the set of reads are described and already sorted by start / end coordinates in the “Reads” array. The “Positions” array captures the unique set of start coordinates in the “Reads” array, and the “Ends” array captures the unique set of end coordinates once the start coordinates are fixed. The “Position starts” array indexes into the “Ends” array indicating where each unique combination of start, end coordinates starts. The “End starts” array indexes into the “Reads” array indicating where each unique combination of start, end coordinates can be found in the reads. As will be understood, the sparse representation enables efficient storage of the sorting features. Although used here for illustrative purposes with respect to the start and end position, it will be understood that similar index data structures may be generated to represent other combinations of sorting features.
[0060] Next, once the index is formed, the block of reads is iterated on while skipping reads which have already been clustered. Specifically, for each read i not yet clustered with coordinates (s, e), the following steps may be taken:PATENT Client Reference No.: P39621-WO-1 Iterate over the following coordinates in the 2D square [s, s+w] x [e, e+w] and cluster all reads present in this window (if not already clustered) with read i. The newly clustered reads may be marked as having been clustered so that they are not considered in further iterations. An example of such clustering is provided in FIG.4, which illustrates positional clustering with inexact matching using an anchor read in the bottom left corner of a square representing the cluster, in accordance with embodiments of the present disclosure. As can be seen in plot 400 of FIG.4, while the actual cluster square 404 may be effectively the same size as that of FIG.2, in this case, the clustering may be done by using an anchor read (i.e. the read around which the cluster is formed) 402 at the bottom left corner of the cluster square rather than clustering around a read in the center as shown in FIG.2. It will be understood that clustering with an anchor read positioned like this allows faster processing. In this instance, only unclustered reads are considered as anchor reads and the algorithm proceeds by clustering all unclustered reads in the square defined using the anchor read as the bottom left corner. Once reads are clustered they are removed from further consideration, so that each read only is assigned a cluster one time. In other embodiments, it will be appreciated that other distance metrics can be used as well, such as Euclidean distance or a Manhattan distance (sum of absolute differences for both coordinates). In other embodiments, the location of the anchor read can be positioned in a different location within the 2D square, such as the upper right corner, upper left corner, etc. It will be appreciated that the choice of position used for the anchor read can increase the complexity of the search when compared to the choice of lower left corner as shown in FIG.4. Once clustering is completed over the 2D square [s, s+w] x [e, e+w] an additional rescue step may be taken depending on the maximum variation in end coordinates within the cluster. FIG.5 is an illustration of a rescue step to the positional clustering described in FIG.4, in accordance with embodiments of the present disclosure. As can be seen in FIG.5, the maximum variation in end coordinates within the cluster may be measured by comparing read 506 to anchor read 502. In this case, the maximum variation is within the acceptable tolerance by a value represented by , implying that the anchor read 502 is not actually in the bottom left corner of cluster square 504. To remedy this, an additional rescue step isPATENT Client Reference No.: P39621-WO-1 performed to cluster more reads downward and rightward using . Specifically, is first calculated by subtracting the end coordinate distance from the anchor to the farthest read from . Then, reads that are found in the downward-rightward direction of [ , ] may be clustered with the current anchor read 502, providedthat the reads have not yet been clustered. In FIG.5, for example, read 510 would be added to the cluster by virtue of the rescue step. As above, any newly clustered reads should be marked as such to avoid duplicate clustering. 3. Once clustered using one or more of the methodologies described above, the clusters generated by each read may be output.
[0061] In some cases, there may be a desire for alternative clustering methods. For example, it may be desirable to cluster as many reads as possible to reduce the total number of clusters. In such cases, a maximal clustering algorithm may be used. In an exemplary maximal clustering algorithm with each read with start and end as [s, e] the algorithm iterates over all its left neighbors with smaller start_pos and all its right neighbors with larger start_pos, and clusters all those neighbors whose start_pos and end_pos are within s – w to s + w and e-w to e+w respectively. In such case, the algorithm may not actually track whether a read has been clustered or not, in order to minimize the number of clusters. While this may theoretically result in a cluster with an unbounded tolerance distance for either start_pos or end_pos, in actual observations for typical sequencing applications, this has rarely occurred, and the intracluster maximum tolerance has rarely exceeded a constant for high error read data
[0062] FIG.6 is an illustrative example of maximal clustering in accordance with embodiments of the present disclosure. In the illustration, read 602 represents the first read in a cluster, read 604 is the second read in the cluster which sits within the tolerance distance to first read 602, and read 606 represents the third read in the cluster that sits within the tolerance distance to second read 604. As can be seen from the figure, the distance between first read 602 and the third read 606 may have otherwise been outside the tolerable distance for a cluster depending on the anchoring used, but in this approach maximizing the number of reads in a cluster and minimizing the number of clusters is prioritized. It will be understood that maximal clustering may be used in place of and or in addition to any of the other clustering methodologies described herein, depending on the application.PATENT Client Reference No.: P39621-WO-1
[0063] Once the above steps have been performed for all remaining unclustered reads, the output clusters may be stored and saved. Based on the clusters, PCR duplicates may be identified, and further secondary analysis may proceed with a better indication of the sequencing data of the sample of interest. Random UMI Read Clustering
[0064] Having described a methodology for inexact matching in the case of fixed UMIs, the methodology for the case of random UMIs will now be described. Although there are some differences in the methodology for sequencing applications where random UMIs are used, it will be appreciated that methodologies are conceptually similar. Accordingly, it will be understood that certain steps may reference analogous steps from the fixed UMI Read Clustering described above.
[0065] Once the unclustered sequencing reads have been generated as described above, the reads may be sorted in the following order: 1) alignment start position, 2) alignment end position, and 3) UMI. As can be plainly seen, the sorting order differs from the case of fixed UMIs, since in the case of random UMI there may be instances where reads with different UMIs are nevertheless clustered together because they are close enough with respect to alignment start and end positions.
[0066] Once sorted as above, a pre-processing step may occur in which an index data structure is formed. As with the case of fixed UMIs, the data structure may be a sparse representation of the reads which allows an efficient storage, however, in this case it may store <start_pos, end_pos, UMIs> as opposed to [start_pos, end_pos]. In this case, the sparse representation may also collapse reads with identical start_pos, end_pos, and UMIs. The index allows for an efficient retrieval of reads which could be candidates for clustering to a read under consideration.
[0067] After the index data structure is formed, the method calls for iterating over the set of reads while skipping reads which have been clustered already. Specifically, for each read not yet clustered with coordinates (s, e) and UMI u, the following steps may be taken: 1. Iterate over the coordinates in the 2D square [s, s+w] x [e, e+w]. As described above, may be a predetermined tolerance distance indicative of a desired tolerance for the inexact position matching in each direction. For each UMI u’ in this window, thePATENT Client Reference No.: P39621-WO-1 Hamming distance between UMI u’ and UMI u is calculated and compared to a UMI tolerance distance d. Similar to the positional tolerance distance , the tolerance distance may be a predetermined value indicative of a desired tolerance for the inexact matching of UMIs, and may be determined by observation and / or modified based on the application. If the Hamming distance between UMI u’ and UMI u is less than d then all the reads within the given window with UMI u’ may be clustered to read i. If the Hamming distance between UMI u’ and UMI u is greater than d, then a second distance, a Levenshtein distance may be calculated between UMI u’ and UMI u. If the Levenshtein distance between UMI u’ and UMI u is less than d then all the reads within the given window with UMI u’ may be clustered to read i. Otherwise, the UMI u’ will remain unclustered. It will be appreciated that the inclusion of the calculation and comparison of the Hamming distance between UMIs may speed up the distance calculation process relative to the instance in which a Levenshtein calculation was used exclusively. This is true because it is expected that most UMIs evaluated would have a Hamming distance within d, thus obviating the need to actually calculate the more intensive Levenshtein distance between UMIs. As described previously, any newly clustered reads should be marked as such to avoid clustering again. 2. Similar to as described with respect to fixed UMIs, if the maximum variation in end coordinates in the above step is w’ and w’ < w, the step above may be repeated in the 2D rectangle with coordinates [s, s+w] x [e-(w-w’), e-1] to expand the cluster as described above. 3. In the case of a cluster of a relatively small size with identical coordinates and UMI, an additional rescue step may be performed by iterating over the reads in the 2D square [s-w, s+w] x [e-w, e+w], while only looking at already clustered reads. If a read j (j i) has UMI within the Levenshtein distance d, then examine the cluster for that read to see if other reads in that cluster have UMIs within Levenshtein distance d from read j. If so, then add this cluster to the cluster for read j.
[0068] Since the maximum distance is typically limited, the wavefront aligner (WFA) tool may be leveraged to perform the Levenshtein distance calculation whose time complexity is represented by O(LM) where L is the length of UMI and M is the maximumPATENT Client Reference No.: P39621-WO-1 distance allowed, compared with O(L^2) time complexity for conventional Levenshtein distance calculations
[0069] Although Hamming and Levenshtein distances are described above in comparing random UMIs, it will be understood that the methodology may support gap linear (which linearly penalize gaps in alignments) and gap affine (which allows for large gaps in alignments) distances offered by wavefront aligner (WFA). It will be understood that the support of these distances may make it easier to optimize the methodology for different error profiles, i.e. to weigh insertion, deletion, and mismatch errors differently.
[0070] The worst case number of UMI comparisons might be O(N^2) where N is the number of reads clustered based on start and end coordinates. To reduce the number of comparisons, minimisers may be leveraged to first group similar UMIs and then limit pairwise UMI comparisons to be within each minimiser group. As used herein, a minimiser may be defined as the smallest k-mer (a string of length k) in a UMI where k-mers are sorted by user-specified random order. The calculation of minimisers and grouping by minimisers can be done in O(N) time and space.
[0071] As described above with respect to fixed UMIs, it may be desirable to cluster as many reads as possible to reduce the total number of clusters. In such cases, a maximal clustering algorithm may be used similar to that described above with respect to fixed UMIs in place of or in addition to the clustering described herein. For example, a maximal clustering algorithm can be used to form positional clusters (based solely on start_pos and end_pos of reads) first. Then, within each positional cluster, UMIs may be clustered based on distances as described above.
[0072] If the rescue step in 3. above was not successful then output the cluster generated by read i. As described above, once the above steps have been performed for all remaining unclustered reads, the output clusters may be stored and saved. Based on the clusters, PCR duplicates may be identified, and further secondary analysis may proceed with a better indication of the sequencing data of the sample of interest.
[0073] The methodologies above are only illustrative and not intended to be limiting examples, as alternative features and processes may be used. For example, as described previously, while detailed examples are provided herein for using start position, end position,PATENT Client Reference No.: P39621-WO-1 and UMI for sorting and clustering reads, other features may also be used. For example, the number of soft clippings and / or the number of split alignments may be compared to determine whether a read belongs to the same cluster as another. As other examples, an edit distance to reference sequence, a fragment length, and a mapping quality may be compared to assist in clustering reads in order to identify PCR duplicates. In some cases, these additional features may be helpful to address situations where the presence of structural variants may create ambiguities in clustering.
[0074] Similarly, while equivalent distances are illustrated for the start and end positions, it will be understood that there may be different tolerant distances for the start and end positions, depending on the application and sample of interest, for example. Additionally, other criteria may be used to control for the tolerable variations in start and end positions. For example, instead of limiting the maximum variation in either start or end positions, we may limit the total of the variations in start and end positions. Note that when controlling for the total of the variations in start and end positions, the area of interest for clustering would be a rhombus instead of square around the anchor read.
[0075] As previously described, the implementation of the above methodologies is not only effective in identifying duplicates for higher error technologies, but is also very space efficient. In particular, the minimized physical memory utilized by the methodologies allows for a larger number of reads to be de-duplicated together which reduces the likelihood of artificially splitting read clusters as a result of insufficient memory. Thus, the methodologies also improve the storage and efficiency of de-duplication. Method for Clustering to Identify PCR Duplicates
[0076] FIG.7 is a flowchart of a method 700 of clustering to identify PCR duplicates, in accordance with embodiments of the present disclosure. This generalized method of clustering may apply to either of the UMI types (fixed or random) as will be evident to one of skill in the art.
[0077] At step 710, a sample may be sequenced to generate a plurality of unclustered sequence reads. The sample may be sequenced using any suitable sequencing technology as described herein. In some instances, the data including the unclustered sequence reads may be received by a processing module configured to cluster the reads and identify PCRPATENT Client Reference No.: P39621-WO-1 duplicates based on the clustering. Each of the unclustered sequencing reads may include a set of features that can be used for sorting the reads. For example, the sequencing reads may have features associated with the alignment to a reference sequence, including an alignment start position and an alignment end position. As another example, the sequencing reads may also have UMIs associated therewith. For instance, the UMIs may be known fixed UMIs or random UMIs. Other features used for sorting may include, but are not limited to: a number of 5' and 3' soft-clipped basepairs, a number of supplementary alignments, an edit distance to reference sequence, a length of the fragment, and a mapping quality.
[0078] At step 720, the plurality of unclustered reads may be sorted based on a primary sorting feature selected from the set of sorting features. For example, if the sequences include fixed UMIs and the set of sorting features includes UMIs, start position, and end position, the primary sorting feature may be the fixed UMIs. As another example, if the sequences include random UMIs and the set of sorting features includes UMIs, start position, and end position, then the primary sorting feature may be start position.
[0079] At step 730, the plurality of unclustered reads may be sorted based on the remaining sorting features other than the primary sorting feature. For example, if the primary sorting feature is the fixed UMI, the unclustered reads may subsequently be sorted by start position and end position. Conversely, if the primary sorting is the starting position, the unclustered reads may be subsequently sorted by end position and UMI.
[0080] At step 740, once sorted by all of the sorting features, the reads may be iteratively indexed using a sparse representation that stores values of the sorting features and clusters together reads with identical values of the sorting features. The sparse representation may be generated similar to that described above with reference to FIG.3. For example, the sparse representation of the reads may allow efficient storage of the alignment start and end positions represented as coordinates [start_pos, end_pos] while also collapsing the reads with identical start_pos and end_pos. As another example, the sparse representation may store unique <start_pos, end_pos, UMIs> while also collapsing / clustering reads with identical start_pos, end_pos, and UMIs.
[0081] At step 750, once indexed, the unclustered reads may be iterated on over predetermined tolerance distance windows to cluster together reads which have values of the sorting features within a predetermined tolerance distance of each other. For example, in thePATENT Client Reference No.: P39621-WO-1 case where fixed UMIs are used, unclustered reads with identical UMIs but with start position and end position within a predetermined tolerance distance of an initial anchor read may be clustered together, as described above with respect to FIGS.2, 4, and or 6. As another example, in the case where random UMIs are used, unclustered reads falling within a predetermined distance with respect to starting and end position, and a second predetermined distance as measured by Hamming and or Levenshtein distance with respect to the UMI sequence may be clustered together as described above.
[0082] As optional steps not depicted in flowchart, rescue steps may be performed after various iterations of clustering to expand the windows of clustering. For example, as described with respect to FIG.5, a maximum variation in a sorting feature (such as the end position) for a given cluster may be determined, and if that maximum variation is less than the predetermined tolerance distance, the window may be expanded to include an additional tolerance distance based on the difference between the maximum variation and the predetermined tolerance.
[0083] At step 760, reads which have been clustered in the steps above may be marked as clustered. Reads may be marked as clustered by any suitable mechanism, so long as it is reflected in a way to potentially avoid further clustering of an already clustered read. As described above, in some embodiments, it may be desirable to avoid clustering reads more than once, but in other embodiments, it may be desirable to re-cluster already clustered reads (for example in the case of maximal clustering described above with respect to FIG.6).
[0084] As indicated at 765, steps 740 through 760 may be iterated and performed until all reads have been clustered, or until some condition has been met with respect to clustering of the reads (for example, until a threshold of the number of clusters is met or another condition of relevance is met). Once complete, at step 770, the clustered reads may be output. In some cases, the clusters may be stored and information associated with the clusters and individual reads may be stored in memory or transmitted for further evaluation by the processing module or another processing module.
[0085] At step 780, PCR duplicate reads may be identified based on the output clusters. For example, it may be determined that reads within a given cluster are all PCR duplicates of one another. Subsequent steps may be performed that use the identification of PCR duplicates to make a sequencing determination about a sample of interest and / or modify aPATENT Client Reference No.: P39621-WO-1 previous determination made thereof. As an example, depending on the application, the set of identified PCR duplicate read sequences may be replaced by a consensus sequence which is a more accurate read of the original molecule of interest, which may aid in making further determinations regarding the molecule of interest. Exemplary Computer System
[0086] FIG.8 illustrates an example computer system that may be utilized to implement techniques disclosed herein. Any of the computer systems mentioned herein, such as for hosting the systems and implementing the processes described for calculating variants, may utilize any suitable number of subsystems. Examples of such subsystems are shown in FIG.8 as included in computer system 800. In some embodiments, a computer system includes a single computer apparatus, where the subsystems can be the components of the computer apparatus. In other embodiments, a computer system can include multiple computer apparatuses, each being a subsystem, with internal components. A computer system can include desktop and laptop computers, tablets, mobile phones, telecommunication devices or other mobile devices. In some embodiments, a cloud infrastructure (e.g., Amazon Web Services), a graphical processing unit (GPU), etc., can be used to implement the disclosed techniques.
[0087] The subsystems shown in FIG.8 are interconnected via a system bus 75. Additional subsystems such as a printer 74, keyboard 78, storage device(s) 79, monitor 76, which is coupled to display adapter 82, and others are shown. Peripherals and input / output (I / O) devices, which couple to I / O controller 71, can be connected to the computer system by any number of means known in the art such as input / output (I / O) port 77 (e.g., USB, FireWire®). For example, I / O port 77 or external interface 81 (e.g. Ethernet, Wi-Fi, etc.) can be used to connect computer system 800 to a wide area network such as the Internet, a mouse input device, or a scanner.
[0088] The interconnection via system bus 75 allows the central processor 73 to communicate with each subsystem and to control the execution of a plurality of instructions from system memory 72 or the storage device(s) 79 (e.g., a fixed disk, such as a hard drive, or optical disk), as well as the exchange of information between subsystems. The system memory 72 and / or the storage device(s) 79 may embody a computer readable medium. Another subsystem is a data collection device 85, such as a camera, microphone,PATENT Client Reference No.: P39621-WO-1 accelerometer, and the like. Any of the data mentioned herein can be output from one component to another component and can be output to the user.
[0089] A sequencer instrument 90 (e.g., a nanopore sequencer), which may sometimes be referred to as a sequencing device, is connected through external interface 81 for providing sequencing data to a data collection device 85 and / or storage devices 79.
[0090] The sequencer instrument 90 may be any sequencer instrument including, but not limited to, the following: sequencing-by-hybridization platforms from Affymetrix, Inc., now ThermoFisher Scientific, (Sunnyvale, Calif.), sequencing-by-synthesis platforms from Illumina, Inc. (San Diego, Calif.) and Helicos Biosciences (Cambridge, Mass.), and sequencing-by-ligation platform from Applied Biosystems (Foster City, Calif.). Other sequencing technologies include, but are not limited to, the Ion Torrent technology from ThermoFisher Scientific (Waltham, MA); as well as nanopore sequencing from Roche Sequencing Solutions, Inc. (Santa Clara, Cal.) and Oxford Nanopore Technologies (Oxford, United Kingdom).
[0091] A computer system 800 can include a plurality of the same components or subsystems, e.g., connected together by external interface 81 or by an internal interface such as the bus 75. In some embodiments, computer systems, subsystems, or apparatuses can communicate over a network. In such instances, one computer can be considered a client and another computer a server, where each can be part of the same computer system. A client and a server can each include multiple systems, subsystems, or components.
[0092] Aspects of embodiments can be implemented in the form of control logic using hardware (e.g. an application specific integrated circuit or field programmable gate array) and / or using computer software with a generally programmable processor in a modular or integrated manner. As used herein, a processor includes a single-core processor, multi-core processor on a same integrated chip, or multiple processing units on a single circuit board or networked. Based on the disclosure and teachings provided herein, a person of ordinary skill in the art will know and appreciate other ways and / or methods to implement embodiments of the present invention using hardware and a combination of hardware and software.
[0093] Machine learning models utilized herein may include one or more of a Naive Bayes (NB) model, a logistic regression (LR) model, a random forest (RF) model, a support
Claims
PATENT Client Reference No.: P39621-WO-1 vector machine (SVM) model, an artificial neural network model, a multilayer perceptron (MLP) model, a convolutional neural network (CNN), a Large Language model (LLM), and / or other machine learning or deep learning models, etc. The machine learning models can be updated / trained using a supervised learning technique, an unsupervised learning technique, etc. [0094] Any of the software components or functions described in this application may be implemented as software code to be executed by a processor using any suitable computer language such as, for example, Java, C, C++, C#, Objective-C, Swift, or scripting language such as Perl or Python using, for example, conventional or object-oriented techniques. The software code may be stored as a series of instructions or commands on a computer readable medium for storage and / or transmission. A suitable non-transitory computer readable medium can include random access memory (RAM), a read only memory (ROM), a magnetic medium such as a hard-drive or a floppy disk, or an optical medium such as a compact disk (CD) or DVD (digital versatile disk), flash memory, and the like. The computer readable medium may be any combination of such storage or transmission devices. [0095] Such programs may also be encoded and transmitted using carrier signals adapted for transmission via wired, optical, and / or wireless networks conforming to a variety of protocols, including the Internet. As such, a computer readable medium may be created using a data signal encoded with such programs. Computer readable media encoded with the program code may be packaged with a compatible device or provided separately from other devices (e.g., via Internet download). Any such computer readable medium may reside on or within a single computer product (e.g. a hard drive, a CD, or an entire computer system), and may be present on or within different computer products within a system or network. A computer system may include a monitor, printer, or other suitable display for providing any of the results mentioned herein to a user. [0096] Any of the methods described herein may be totally or partially performed with a computer system including one or more processors, which can be configured to perform the steps. Thus, embodiments can be directed to computer systems configured to perform the steps of any of the methods described herein, potentially with different components performing a respective step or a respective group of steps. Although presented as numbered steps, steps of methods herein can be performed at the same time or in a different order. Additionally, portions of these steps may be used with portions of other steps from otherPATENT Client Reference No.: P39621-WO-1 methods. Also, all or portions of a step may be optional. Additionally, any of the steps of any of the methods can be performed with modules, units, circuits, or other means for performing these steps.PATENT Client Reference No.: P39621-WO-1 CLAIMS 1. A method for clustering genetic sequence reads from a sample, the method comprising: sequencing the sample to generate a plurality of unclustered sequence reads, each unclustered sequence read associated with a set of sorting features; sorting the plurality of unclustered sequence reads according to a primary sorting feature; for each subset of unclustered sequence reads with identical values of the primary sorting feature: selecting a first unclustered sequence read of the subset as an anchor read for a cluster of sequence reads, adding to the cluster any unclustered sequence reads that have values of one or more secondary sorting features other than the primary sorting feature that fall within a predetermined tolerance distance of the one or more secondary sorting features of the anchor read, and in response to determining that one or more sequence reads in the subset remain unclustered, repeating the steps of selecting and adding to generate one or more additional clusters until each of the plurality of unclustered sequence reads have been added to a particular cluster; and generating an output that indicates a set of clusters.
2. The method of claim 1, wherein the set of sorting features comprises at least two of the following: a unique molecular identifier (UMI), an alignment starting position, an alignment ending position, a number of 5' and 3' soft-clipped base pairs, a number of supplementary alignments, an edit distance to a reference sequence, a fragment length, and a mapping quality.
3. The method of claim 2, wherein the set of sorting features comprise the UMI, the alignment starting position, and the alignment ending position, further wherein the primary sorting feature is the UMI, and the one or more secondary sorting features comprises the alignment starting position and the alignment ending position.PATENT Client Reference No.: P39621-WO-1 4. The method of claim 2, wherein the set of sorting features comprise the UMI, the alignment starting position, and the alignment ending position, further wherein the primary sorting feature is the alignment starting position, and the one or more secondary sorting features comprise the alignment ending position and the UMI.
5. The method of claim 1, further comprising, for each cluster: determining a maximum variation for the values of the one or more secondary sorting features for the sequence reads in the given cluster; comparing the maximum variation to the predetermined tolerance distance; and if the maximum variation is less than the predetermined tolerance distance, adding to the given cluster any additional unclustered sequence reads that have values of the one or more secondary sorting features that fall within an additional tolerance distance from the values of the one or more secondary sorting feature of the anchor read.
6. The method of claim 1, wherein all sequence reads in a particular cluster are identified as PCR duplicate reads.
7. The method of claim 6, further comprising modifying at least one sequencing determination based on the identification of PCR duplicate reads.
8. The method of claim 1, wherein the predetermined tolerance distance varies according to the secondary sorting feature being evaluated.
9. The method of claim 8, wherein the predetermined tolerance distance associated with an alignment starting position is greater than the predetermined tolerance distance associated with an alignment ending position.
10. A computer-implemented method of identifying PCR duplicates, the method comprising: receiving sequencing data comprising a plurality of sequence reads, wherein each sequence read is associated with at least a unique molecular identifier (UMI) and a positional coordinate comprising an alignment start position and an alignment end position;PATENT Client Reference No.: P39621-WO-1 identifying an anchor coordinate for a cluster, wherein the anchor coordinate is selected as a positional coordinate of an unclustered read in the plurality of sequence reads; determining at least one candidate read in the plurality of sequence reads to add to the cluster, wherein a distance metric for each candidate read is less than a threshold distance from the anchor coordinate; in response to determining that at least one sequence read in the plurality of sequence reads remains unclustered, identifying another anchor coordinate for an additional cluster and determine at least one candidate read in the plurality of sequence reads to add to the additional cluster; repeat the identifying and determining steps until all sequence reads in the plurality of sequence reads are added to a corresponding cluster; store an output identifying the clusters of sequence reads in a memory; and identify PCR duplicates in the plurality of sequence reads based on the output.
11. The method of claim 10, further comprising identifying a UMI for the cluster, wherein determining at least one candidate read further comprises determining that the UMI for the candidate read matches the UMI for the cluster.
12. The method of claim 10, further comprising identifying a UMI for the cluster, wherein determining at least one candidate read further comprises determining whether the UMI for the at least one candidate read is within a Hamming distance of the UMI for the cluster.
13. The method of claim 12, wherein: the at least one candidate read is added to the cluster in response to determining that the UMI for the candidate read is within the Hamming distance of the UMI for the cluster; or the candidate read is added to the cluster in response to determining that the UMI for the candidate read is not within the Hamming distance of the UMI for the cluster, but is within a Levenshtein distance of the UMI for the cluster.
14. The method of claim 12, wherein the at least one candidate read is not added to the cluster in response to determining that the UMI for the at least one candidate read is not within the Hamming distance or a Levenshtein distance of the UMI for the cluster.
Citation Information
Patent Citations
Nanopore Based Molecular Detection and Sequencing
US20130244340A1
DNA sequencing by synthesis using modified nucleotides and nanopore detection
US20130264207A1
Nucleic acid sequencing using tags
US20140134616A1
Nucleic acid sequencing by nanopore detection of tag molecules
US20150119259A1
Methods for creating bilayers for use with nanopore sensors
US20150337366A1