A lead-layer specific noise model for analyzing DNA data
A site-specific noise model with Bayesian inference and read layer-specific parameters addresses the challenge of distinguishing true positives from false positives in DNA sequencing, improving mutation detection reliability.
Patent Information
- Application Number
- JP2022541189
- Authority / Receiving Office
- JP · JP
- Patent Type
- Patents
- Current Assignee / Owner
- Priority Date
- 2019-09-09
- Filing Date
- 2020-09-08
- Publication Date
- 2024-12-17
- Estimated Expiration
- 2040-09-08
AI Technical Summary
Existing methods struggle to distinguish true positives from false positives in DNA sequencing data, particularly for detecting mutations associated with cancer, due to low levels of circulating tumor DNA in cell-free DNA, leading to unreliable results.
A site-specific noise model is developed, classified into multiple lead layers, incorporating Bayesian inference and read layer-specific parameters such as trinucleotide context and sequencing depth, to enhance the identification of true positives by stratifying sequencing reads and combining noise model outputs.
The model improves the sensitivity and specificity of variant calling by accurately distinguishing true positives from false positives, enhancing the reliability of mutation detection in DNA sequencing data.
Smart Images

Figure 0007601883000083 
Figure 0007601883000084 
Figure 0007601883000085
Abstract
Description
[Technical Field]
[0001] The present disclosure relates generally to noise models for determining quality scores for nucleic acid sequencing datasets. [Background technology]
[0002] Computational techniques can be used on DNA sequencing data to identify mutations or variants in DNA that may correspond to various types of cancer or other diseases. Therefore, cancer diagnosis or prediction can be performed by analyzing biological samples such as tissue biopsies or blood collected from individuals, animals, plants, etc. Detecting DNA derived from tumor cells in blood samples is challenging because circulating tumor DNA (ctDNA) is present at low levels relative to other molecules in cell-free DNA (cfDNA) extracted from blood. The inability of existing methods to distinguish true positives (e.g., indicative of cancer in a subject) from signal noise reduces the ability of known and future systems to distinguish true positives from false positives caused by noise sources, which can lead to unreliable results for variant calling or other types of analysis. [Prior art documents] [Patent documents]
[0003] [Patent Document 1] U.S. Patent Application Serial No. 16 / 153,593 Summary of the Invention
[0004] Disclosed herein are systems and methods for training and applying site-specific noise models categorized into multiple read tiers. The noise models can determine the likelihood of true positives in targeted sequencing. True positives can include single-base variants, insertions, or deletions of base pairs. Specifically, the models can use Bayesian inference to determine, for example, a noise rate or level for each position in a nucleic acid sequence, indicating the expected likelihood of several mutations. Each model can be specific to a read tier. The read tier can be determined based on whether potential variant locations are located in overlapping and / or complementary regions of processed sequencing reads. Each read tier-specific model can be a hierarchical model that accounts for covariates (e.g., trinucleotide context, mappability, or segment duplication) and various types of parameters specific to the read tier (e.g., mixture components or sequencing read depth). The models can be trained from sequencing reads of healthy subjects, also stratified by read tier. The outputs of different noise models can be combined to generate an overall quality score. The overall pipeline incorporating various read tier models can identify true positives and exclude false positives with higher sensitivity when compared to a single model that does not separate sequencing reads by read tier.
[0005] By way of example, in various embodiments, a method for processing a DNA sequencing dataset of a sample (e.g., an individual) can include accessing a DNA sequencing dataset generated by DNA sequencing, the DNA sequencing dataset including a plurality of processed sequence reads including variant locations. The method can also include stratifying the plurality of processed sequence reads into a plurality of read stratums. The method can further include determining, for each read stratum, a stratified sequencing depth at the variant locations. The method can further include determining, for each read stratum, one or more noise parameters conditioned on the stratified sequencing depth of the read stratum, the one or more noise parameters corresponding to a read stratum-specific noise model. The method can further include generating, for each read stratum, a read stratum-specific noise model output based on the one or more noise parameters conditioned on the stratified sequencing depth of the read stratum. The method can further include combining the generated noise model outputs to produce a combined result. The combined result can then represent the likelihood that the total variant count for the observed data is greater than or equal to the total variant count observed in multiple processed sequencing reads due to noise.
[0006] In one or more embodiments, the multiple lead layers include one or more of: (1) a double-stranded, sewn lead layer; (2) a double-stranded, unsewn lead layer; (3) a single-stranded, sewn lead layer; and (4) a single-stranded, unsewn lead layer.
[0007] In one or more embodiments, the mutation at the variant location is one of a single base variant, an insertion, and a deletion.
[0008] In one or more embodiments, the method may further include determining a quality score of the combined result, where the quality score is a Phred scale score.
[0009] In one or more embodiments, the method may further include indicating, in response to the quality score being higher than a predetermined threshold, that the sample has a high likelihood of harboring a mutation at the variant location.
[0010] In one or more embodiments, for each read stratum, determining one or more noise parameters conditional on the stratified sequencing depth of the read stratum can include accessing a parameter distribution specific to the read stratum, the parameter distribution describing the distribution of a set of DNA sequencing samples associated with the read stratum. The noise parameters are determined from the parameter distribution.
[0011] In one or more embodiments, for each lead stratum, the set of DNA sequencing samples associated with the lead stratum includes sequence reads stratified into the lead stratum and corresponds to one or more healthy individuals.
[0012] In one or more embodiments, for the lead layer, the lead layer specific noise model is a Bayesian hierarchical model, and the parameter distribution is based on a gamma distribution.
[0013] In one or more embodiments, a first noise parameter corresponding to a noise model specific to a first lead layer has a different value than a corresponding second noise parameter corresponding to a noise model specific to a second lead layer.
[0014] In one or more embodiments, for each lead stratum, the determined one or more noise parameters include a mean of the noise distribution conditional on the stratified sequencing depth of the lead stratum.
[0015] In one or more embodiments, each noise distribution is a negative binomial distribution conditional on the stratified sequencing depth of each read stratum.
[0016] In one or more embodiments, for each lead layer, the determined one or more noise parameters further include a variability parameter.
[0017] In one or more embodiments, the output of each noise model is one or more noise parameters conditional on the stratified sequencing depth determined for the read stratum.
[0018] In one or more embodiments, the generated output of each noise model is one or more noise parameters conditional on the stratified sequencing depth determined for the read stratum.
[0019] In one or more embodiments, the generated output of each noise model includes the likelihood that the stratified variant count for the lead stratum exceeds a threshold.
[0020] In one or more embodiments, combining the generated noise model outputs includes combining the mean variant count and variance from each noise model output to produce an overall mean variant count and an overall variability parameter that represents the overall noise distribution for the combined result.
[0021] In one or more embodiments, the overall noise distribution is modeled based on a negative binomial distribution. Determining the overall mean variant count and the overall variability parameter can include determining an average variant count for each lead stratum based on the stratified sequencing depth of the lead stratum. The determining step can also include determining a variance for each lead stratum. The determining step can further include summing the average variant counts for each lead stratum to determine the overall mean variant count. The determining step can further include combining the variances for each lead stratum to determine the overall variance. The determining step can further include determining the overall variance parameter based on the overall mean variant count and the overall variance.
[0022] In one or more embodiments, combining the outputs for each noise model to generate a combined result can include determining an observed stratified variant count for each lead stratum. The combining step can also include determining possible events within each lead stratum that are more likely than the observed stratified variant count for each lead stratum. The combining step can further include identifying combinations of possible events associated with a higher likelihood of occurrence than the observed stratified variant count for each lead stratum. The combining step can further include summing the probabilities of the identified combinations to determine a statistical complement. The combining step can further include determining a likelihood value by subtracting the statistical complement from 1.0.
[0023] In one or more embodiments, a first identified combination comprising one double-stranded read is equal to a second identified combination comprising two single-stranded reads.
[0024] In one or more embodiments, the determined likelihood values are equal to or greater than the likelihood of occurrence of the observed stratified variant counts for each lead stratum.
[0025] In one or more embodiments, the method may further include training a machine learning model to determine the likelihood value.
[0026] In one or more embodiments, the method may further include receiving a bodily fluid sample from the individual. The method may further include performing DNA sequencing on the cfDNA of the bodily fluid sample. The method may further include generating raw sequence reads based on the results of the DNA sequencing. The method may further include collapsing and stitching the raw sequence reads to generate a plurality of processed sequence reads.
[0027] In one or more embodiments, the bodily fluid sample is a sample of one of the individual's blood, whole blood, plasma, serum, urine, cerebrospinal fluid, feces, saliva, tears, tissue biopsy, pleural fluid, pericardial fluid, or peritoneal fluid.
[0028] In one or more embodiments, the plurality of processed sequencing reads are sequenced from a tumor biopsy.
[0029] In one or more embodiments, the plurality of processed sequence reads are sequenced from an isolate of cells from blood, the isolate of cells comprising at least buffy coat leukocytes or CD4+ cells.
[0030] In one or more embodiments, the DNA sequencing is a type of massively parallel DNA sequencing.
[0031] In various embodiments, a non-transitory computer-readable medium includes instructions that, when executed by one or more processors, cause the one or more processors to perform any of the steps described above and disclosed herein.
[0032] Further, in various embodiments, a system is provided having a computer processor and a memory that stores computer program instructions, whereby the instructions are executed by the computer processor to cause the processor to perform any of the steps described above and disclosed herein.
[0033] Embodiments according to the present invention are specifically disclosed in the appended claims directed to methods and computer program products, and any feature, e.g., a method, recited in one claim category may also be claimed in another claim category, e.g., a computer program product, a system, or a storage medium. Dependencies or references in the appended claims are chosen for public reasons only. However, any subject matter resulting from deliberate reference to any previous claim (especially multiple dependencies) is also claimable, and therefore any combination of claims and their features is disclosed and claimable regardless of the dependencies chosen in the appended claims. Patentable subject matter includes not only combinations of features recited in the appended claims, but also any other combination of features in the claims, and each feature recited in a claim may be combined with any other feature or other combination of features in the claims. Moreover, any of the embodiments and features described or depicted herein may be claimed in a separate claim and / or in any combination with any embodiment or feature described or depicted herein or with any feature of the appended claims. [Brief explanation of the drawings]
[0034] [Figure 1] 1 is a flowchart of a method for preparing a nucleic acid sample for sequencing according to various embodiments of the present disclosure. [Figure 2] FIG. 1 is a block diagram of a processing system for processing sequence reads according to various embodiments of the present disclosure. [Figure 3] 1 is a flowchart of a method for determining variants of a sequencing read according to various embodiments of the present disclosure. [Figure 4] FIG. 1 is a diagram of an example application of a Bayesian hierarchical model, according to various embodiments of the present disclosure. [Figure 5A] FIG. 1 illustrates dependencies between parameters and sub-models of a Bayesian hierarchical model for determining true single-base variants, according to various embodiments of the present disclosure. [Figure 5B] FIG. 1 illustrates the dependencies between parameters and sub-models of a Bayesian hierarchical model for determining true insertions or deletions, according to various embodiments of the present disclosure. [Figure 6A] FIG. 10 is a diagram of a distribution of noise rates associated with a Bayesian hierarchical model, according to various embodiments of the present disclosure. [Figure 6B] 10A-10C are diagrams relating to distributions of parameters given by variants associated with Bayesian hierarchical models, according to various embodiments of the present disclosure. [Figure 7A] FIG. 10 illustrates determining parameters by fitting a Bayesian hierarchical model, according to various embodiments of the present disclosure. [Figure 7B] FIG. 10 illustrates determining the likelihood of a false positive using parameters from a Bayesian hierarchical model, according to various embodiments of the present disclosure. [Figure 8A] FIG. 1 illustrates different lead layers of a sequence read, according to various embodiments of the present disclosure. [Figure 8B] 8B is a graph illustrating experimental results showing different qualities of the lead layer of FIG. 8A, according to various embodiments of the present disclosure. [Figure 8C] Graph illustrating experimental results of a first lead tier obtained by stratifying sequence reads into lead tiers and further into sub-read tiers based on the type of nucleotide substitution, according to various embodiments of the present disclosure. [Figure 8D]Graph illustrating experimental results of a second read tier obtained by stratifying sequence reads into read tiers and further into sub-read tiers based on the type of nucleotide substitution, according to various embodiments of the present disclosure. [Figure 8E] Graph illustrating experimental results of a third read tier obtained by stratifying sequence reads into read tiers and further into sub-read tiers based on the type of nucleotide substitution, according to various embodiments of the present disclosure. [Figure 8F] Graph illustrating experimental results of a fourth read tier obtained by stratifying sequence reads into read tiers and further into sub-read tiers based on the type of nucleotide substitution, according to various embodiments of the present disclosure. [Figure 8G] Graph illustrating experimental results of a fifth read tier obtained by stratifying sequence reads into read tiers and further into sub-read tiers based on the type of nucleotide substitution, according to various embodiments of the present disclosure. [Figure 8H] 1 is a graph illustrating average error rates across four read tiers based on alternative allele types, according to various embodiments of the present disclosure. [Figure 9] 1 is a flowchart depicting a process for analyzing DNA sequencing samples using a stratified noise model, according to various embodiments of the present disclosure. [Figure 10] 10 is a flowchart depicting a process for combining the outputs of stratified noise models for different lead strata using moment matching, according to various embodiments of the present disclosure. [Figure 11A] 10 is a flowchart depicting a process for combining the outputs of stratified noise models for different lead strata using integration, according to various embodiments of the present disclosure. [Figure 11B] FIG. 10 illustrates counting more extreme events in a multi-dimensional space according to various embodiments of the present disclosure. [Figure 12A]10 is an exemplary plot of observed quality scores versus default quality scores, according to various embodiments of the present disclosure. [Figure 12B] 10 is another exemplary plot of observed quality scores versus default quality scores, according to various embodiments of the present disclosure. [Figure 13A] 10 is a graph illustrating experimental results of quality scores using lead layers, according to various embodiments of the present disclosure. [Figure 13B] 10 is a graph illustrating experimental results of quality scores using a noise model without separating lead strata, according to various embodiments of the present disclosure. [Figure 14] 1 is a flowchart depicting an exemplary process for identifying potential mutation sites in an individual, according to various embodiments of the present disclosure. [Figure 15] FIG. 1 is a block diagram of an exemplary computing device in accordance with various embodiments of the present disclosure. DETAILED DESCRIPTION OF THE INVENTION
[0035] The figures depict embodiments of the present invention for illustrative purposes only. Those skilled in the art will readily recognize from the following description that alternative embodiments of the structures and methods shown herein may be used without departing from the principles of the present invention as described herein.
[0036] I. Definition The term "individual" refers to an individual human being. The term "healthy individual" refers to an individual who is presumed to be cancer or disease-free. The term "subject" refers to an individual being tested for cancer or disease.
[0037] The term "sequence read" refers to a nucleotide sequence read from a sample obtained from an individual. Sequence reads can be obtained through a variety of methods known in the art.
[0038] The term "read segment" or "read" refers to any nucleotide sequence, including a sequence read obtained from an individual and / or a nucleotide sequence obtained from an initial sequence read from a sample obtained from an individual. For example, a read segment can refer to an aligned sequence read, a collapsed sequence read, or a stitched read. Additionally, a read segment can refer to individual nucleotide bases, such as single-base variants.
[0039] The term "single nucleotide variant" or "SNV" refers to a substitution of one nucleotide with a different nucleotide at a position (e.g., site) in a nucleotide sequence, e.g., a sequencing read, from an individual. A substitution of a first nucleobase X with a second nucleobase Y may be represented as "X>Y." For example, a cytosine to thymine SNV may be represented as "C>T."
[0040] The term "indel" refers to an insertion or deletion of one or more base pairs having a length and a position within a sequence read (which may also be referred to as an anchor position). An insertion corresponds to a positive length, and a deletion corresponds to a negative length.
[0041] The term "variant" refers to one or more SNVs or indels. A variant location refers to a location of interest within a DNA sequence that potentially contains an SNV or indel.
[0042] The term "true positive" refers to real biology, e.g., a mutation that indicates the presence of a potential cancer, disease, or germline mutation in an individual. A true positive is not caused by naturally occurring mutations in healthy individuals (e.g., recurrent mutations) or other sources of artifact, such as process errors during assay preparation of nucleic acid samples.
[0043] The term "false positive" refers to a mutation that is erroneously determined to be a true positive. In general, false positives may be more likely to occur when processing sequence reads associated with a higher average noise rate or a higher noise rate uncertainty.
[0044] The term "cell-free DNA" or "cfDNA" refers to nucleic acid fragments that circulate within an individual's body (e.g., in the bloodstream) and originate from one or more healthy cells and / or from one or more cancer cells.
[0045] The term "circulating tumor DNA" or "ctDNA" refers to nucleic acid fragments originating from tumor cells or other types of cancer cells that can be released into an individual's bloodstream as a result of biological processes such as apoptosis or necrosis of dying cells, or that can be actively released by viable tumor cells.
[0046] The term "alternate allele" or "ALT" refers to an allele that has one or more mutations relative to a reference allele that corresponds to, for example, a known gene.
[0047] The term "sequencing depth" or "depth" refers to the total number of read segments from a sample obtained from an individual.
[0048] The term "alternative depth" or "AD" refers to the number of read segments in a sample that support ALT, e.g., contain ALT mutations.
[0049] The term "alternative frequency" or "AF" refers to the frequency of a given ALT. AF can be determined by dividing the corresponding AD of a sample by the depth of the sample for the given ALT.
[0050] II. Exemplary Assay Protocols 1 is a flowchart of a method 100 for preparing a nucleic acid sample for sequencing, according to various embodiments. Method 100 includes, but is not limited to, the following steps. For example, any step of method 100 can include a quantification substep for quality control or other laboratory assay procedures known to those skilled in the art. Method 100 can accommodate massively parallel DNA sequencing, e.g., a type of next-generation sequencing (NGS).
[0051] In step 110, a nucleic acid sample (DNA or RNA) is extracted from a subject. The subject may be an individual. The sample may be any subset of the human genome or the entire genome. The sample can be extracted from a subject known to have or suspected of having cancer. The sample can include blood, plasma, serum, urine, feces, saliva, other types of bodily fluids, or any combination thereof. In some embodiments, methods for collecting a blood sample (e.g., syringe or finger prick) may be less invasive than procedures for obtaining a tissue biopsy, which may require surgery. The extracted sample may include cfDNA and / or ctDNA. In healthy individuals, the human body can naturally remove cfDNA and other cellular debris. If the subject has cancer or disease, ctDNA in the extracted sample may be present at detectable levels for diagnosis.
[0052] In step 120, a sequencing library is prepared. During library preparation, a nucleic acid sample is randomly cleaved into thousands or millions of fragments. Unique molecular identifiers (UMIs) are added to nucleic acid fragments (e.g., DNA fragments) through adapter ligation. UMIs are short nucleic acid sequences (e.g., 4-10 base pairs) that are added to the ends of DNA fragments during adapter ligation. In some embodiments, UMIs are degenerate base pairs that serve as unique tags that can be used to identify sequence reads originating from a particular DNA fragment. During PCR amplification following adapter ligation, the UMIs are replicated along with the attached DNA fragments, providing a way to identify sequence reads originating from the same original fragment in downstream analysis.
[0053] In step 130, target DNA sequences are enriched from the library. During enrichment, hybridization probes (also referred to herein as "probes") are used to target and destroy nucleic acid fragments that provide information about the presence or absence of cancer (or disease), cancer status, or cancer classification (e.g., cancer type or tissue of origin). For a given workflow, probes can be designed to anneal (or hybridize) to a target (complementary) strand of DNA or RNA. The target strand can be the "positive" strand (e.g., the strand transcribed into mRNA and then translated into protein) or the complementary "negative" strand. Probes can be tens, hundreds, or thousands of base pairs in length. In some embodiments, probes are designed based on gene panels to analyze specific mutations or target regions of the genome (e.g., of a human or another organism) suspected to correspond to some cancer or other type of disease. Furthermore, probes can cover overlapping portions of the target regions. By using a targeted gene panel rather than sequencing every expressed gene in the genome, also known as "whole exome sequencing," method 100 can be used to increase the sequencing depth of the target region, where depth refers to the count of the number of times a given target sequence in a sample is sequenced. Increasing sequencing depth reduces the required input amount of nucleic acid sample. After the hybridization step, hybridized nucleic acid fragments can be captured and amplified using PCR.
[0054] In step 140, sequence reads are generated from the enriched DNA sequences. Sequencing data can be obtained from the enriched DNA sequences by means known in the art. For example, method 100 can include next-generation sequencing (NGS) techniques, including synthesis technology (Illumina), pyrosequencing (454 Life Sciences), ion semiconductor technology (Ion Torrent sequencing), single-molecule real-time sequencing (Pacific Biosciences), sequencing by ligation (SOLiD sequencing), nanopore sequencing (Oxford Nanopore Technologies), or paired-end sequencing. In some embodiments, massively parallel sequencing is performed using sequencing-by-synthesis with reversible dye terminators.
[0055] In some embodiments, the sequence reads can be aligned to a reference genome using methods known in the art to determine alignment position information. The alignment position information can indicate the start and end positions of a region in the reference genome corresponding to the start and end nucleotide bases of a given sequence read. The alignment position information can also include the sequence read length, which can be determined from the start and end positions. The region in the reference genome can be associated with a gene or a segment of a gene.
[0056] In various embodiments, a sequence read consists of a read pair denoted as R1 and R2. For example, the first read R1 can be sequenced from a first end of a nucleic acid fragment, while the second read R2 can be sequenced from a second end of the nucleic acid fragment. Thus, the nucleotide base pairs of the first read R1 and the second read R2 can be consistently aligned (e.g., in opposite orientations) with the nucleotide bases of a reference genome. The alignment position information obtained from the read pair R1 and R2 can include a start position in the reference genome corresponding to the end of the first read (e.g., R1) and an end position in the reference genome corresponding to the end of the second read (e.g., R2). In other words, the start and end positions in the reference genome represent possible locations in the reference genome to which the nucleic acid fragment corresponds. As described below with reference to FIG. 2, an output file having a SAM (sequence alignment map) format or a BAM (binary) format can be generated and output for further analysis, such as variant calling.
[0057] III. Exemplary Processing System FIG. 2 is a block diagram of a processing system 200 for processing sequence reads, according to various embodiments. The processing system 200 includes a sequence processor 205, a sequence database 210, a model database 215, a machine learning engine 220, a model 225 (e.g., a Bayesian hierarchical model corresponding to different read strata), a parameter database 230, a score engine 235, and a variant caller 240. FIG. 3 is a flowchart of a method 300 for determining variants of sequence reads, according to various embodiments. In some embodiments, the processing system 200 implements the method 300 to perform variant calling (e.g., for SNVs and / or indels) based on input sequencing data. Additionally, the processing system 200 can obtain the input sequencing data from an output file associated with a nucleic acid sample prepared using the method 100 described above. The method 300 includes, but is not limited to, the following steps, which are described with respect to the components of the processing system 200. In some embodiments, one or more steps of method 300 can be replaced by steps of a different process for generating variant calls using a variant calling format (VCF), such as, for example, HaplotypeCaller, VarScan, Strelka, or SomaticSniper.
[0058] In step 300, the sequence processor 205 collapses sequence reads of the input sequencing data. In some embodiments, collapsing the sequence reads includes collapsing multiple sequence reads into a consensus sequence to determine the most likely sequence of a nucleic acid fragment or portion thereof using UMIs and, optionally, alignment position information from the sequencing data in the output file (e.g., from method 100 illustrated in FIG. 1). Because UMIs are replicated along with the ligated nucleic acid fragments through enrichment and PCR, the sequence processor 205 can determine that several sequence reads originated from the same molecule within the nucleic acid sample. In some embodiments, sequence reads that have the same or similar alignment position information (e.g., start and end positions within a threshold offset) and contain a common UMI are collapsed, and the sequence processor 205 generates a collapsed read (also referred to herein as a consensus read) that represents the nucleic acid fragment. The sequence processor 205 designates a consensus read as a "duplex" if the corresponding pair of collapsed reads have a common UMI indicating that both the positive and negative strands of the originating nucleic acid molecule are captured; otherwise, the collapsed read is designated as a "non-duplex." In some embodiments, the sequence processor 205 can perform other types of error correction on the sequence reads as an alternative to, or in addition to, collapsing the sequence reads.
[0059] In step 305, the sequence processor 205 stitches the collapsed reads based on portions of overlapping nucleotide sequences between two or more sequence reads. In some embodiments, the sequence processor 205 compares the nucleotide sequences between the first read and the second read to determine whether nucleotide base pairs of the first read and the second read overlap. The two sequence reads can also be compared to a reference genome. In an exemplary use case, in response to determining that the overlap (e.g., of a given number of nucleotide bases) between the first read and the second read is greater than a threshold length (e.g., a threshold number of nucleotide bases), the sequence processor 205 designates the first read and the second read as "stitched"; otherwise, the collapsed read is designated as "not stitched." In some embodiments, the first read and the second read are stitched if the overlap is greater than the threshold length and if the overlap is not a sliding overlap. For example, the sliding overlap can include a homopolymer run (e.g., a single repeating nucleotide base), a dinucleotide run (e.g., a two-nucleotide base sequence), or a trinucleotide run (e.g., a three-nucleotide base sequence), where the homopolymer run, dinucleotide run, or trinucleotide run has at least a threshold length of base pairs.
[0060] In step 310, the sequence processor 205 assembles the reads into a path. In some embodiments, the sequence processor 205 assembles the reads to generate a directed graph, such as a de Bruijn graph, for the target region (e.g., a gene). The undirected edges of the directed graph represent sequences of k nucleotide bases (also referred to herein as "k-mers") within the target region, and the edges are connected by vertices (or nodes). The sequence processor 205 aligns the collapsed reads to the directed graph so that any of the collapsed reads can be represented in order by a subset of the edges and corresponding vertices.
[0061] In some embodiments, the sequence processor 205 determines a set of parameters that describe the directed graph and processes the directed graph. Additionally, the set of parameters can include a count of k-mers successfully aligned from the collapsed reads to k-mers represented by nodes or edges in the directed graph. The sequence processor 205 stores the directed graph and the corresponding set of parameters, e.g., in the sequence database 210, which can be retrieved to update the graph or generate a new graph. For example, the sequence processor 205 can generate a condensed version of the directed graph (e.g., or modify an existing graph) based on the set of parameters. In an exemplary use case, to filter out data in the directed graph that has a lower level of importance, the sequence processor 205 removes (e.g., "cuts" or "prunes") nodes or edges with counts below a threshold and retains nodes or edges with counts greater than or equal to the threshold.
[0062] In step 315, variant caller 240 generates variant candidate reads from the pathways assembled by sequence processor 205. The variants can correspond to SNVs or indels. In some embodiments, variant caller 240 can generate variant candidate reads by comparing the directed graph (which may have been compressed by pruning edges or nodes in step 310) to a reference sequence for the target region of the genome. Variant caller 240 can align the edges of the directed graph to the reference sequence and record the genomic locations of incorrect edges and incorrect nucleotide bases adjacent to the edges as candidate variant locations. Additionally, variant caller 240 can generate variant candidate reads based on the sequencing depth of the target region. Specifically, for example, variant caller 240 can be more confident in identifying variants in target regions with greater sequencing depth because a greater number of sequence reads helps resolve variable mismatches or other base pairs between sequences (e.g., using redundancy).
[0063] In some embodiments, variant reads can be classified into different read tiers based on the quality of the variant read. The quality of a variant read can correspond to the location of potential variant locations compared to overlapping and / or complementary locations of the collapsed sequence. During sample preparation (e.g., the library preparation process) in massively parallel sequencing, an individual subject's nucleic acid sample can be randomly cleaved before parallel sequencing is performed. The same copy of a nucleic acid sequence can be randomly cleaved in different ways. Thus, some enriched fragments can have overlapping regions that can be stitched with other enriched fragments, while other enriched fragments do not. Some enriched fragments also have complementary sequences that are also enriched, thus generating double-stranded fragments during sequencing. As a result, variant reads for different sequence locations can correspond to different qualities. For example, a variant read in a location where both complementary strands of a fragment are enriched often has better quality than another variant read in another location that finds support only from single-stranded fragments. Details of the read tiers of variant reads are further discussed in Figures 8A-8B.
[0064] In some embodiments, the variant caller 240 generates variant candidate reads using models 225 to determine expected noise rates for sequence reads from a subject. Each of the models 225 may be a Bayesian hierarchical model. A Bayesian hierarchical model can be one of many possible model architectures that can be used to generate variant candidates, all of which are related to each other in that they model position-specific noise information to improve the sensitivity / specificity of variant calling. More specifically, the machine learning engine 220 trains the models 225 using samples from healthy individuals to model the expected noise rate for each position of a sequence read. In some embodiments, variant reads corresponding to different read stratums can be treated differently by different models, each specific to a particular read stratum. The results of each modeling can be combined to generate a combined result. Details of stratifying read stratums and models are further discussed in Figures 8A-11B.
[0065] Furthermore, multiple different models can be stored in the model database 215 or retrieved for post-application training. For example, a first model is trained to model the SNV noise rate, and a second model is trained to model the indel noise rate. Furthermore, the score engine 235 can use the parameters of the model 225 to determine the likelihood of one or more true positives in the sequence read. The score engine 235 can determine a quality score (e.g., on a logarithmic scale) based on the likelihood. For example, the quality score can be a Phred quality score Q=-10 log 10 P, where P is the likelihood of a candidate incorrect variant call (e.g., a false positive).
[0066] In step 320, the score engine 235 scores the variant reads based on the model 225 or the corresponding likelihood or quality score of a true positive. The training and application of the model 225 is described in more detail below.
[0067] In step 325, processing system 200 outputs the analysis results for the variants. In some embodiments, processing system 200 outputs some or all of the determined candidate variants along with the corresponding scores. Downstream systems, e.g., external to processing system 200 or other components of processing system 200, can use the variants and scores for various applications, including, but not limited to, predicting the presence of cancer, disease, or germline mutations. IV. Illustrative Model
[0068] FIG. 4 is a diagram of an example application of the Bayesian hierarchical model 225, according to various embodiments. Mutation A and mutation B are illustrated as examples for illustrative purposes. In most of this disclosure, including FIG. 4 , mutations are represented as SNVs, but in some embodiments, the description in this disclosure is also applicable to indels or other types of mutations. A first variant read from a first sample corresponds to exemplary mutation A, which is a C>T mutation at position 4 of a first reference allele. The first sample has a first AD of 10 and a first total sequencing depth of 1000. A second variant read from a second sample corresponds to exemplary mutation B, which is a T>G mutation at position 3 of a second reference allele. The second sample has a second AD of 1 and a second total depth of 1200. Based solely on the AD (or AF), mutation A may appear to be a true positive, and mutation B may appear to be a false positive, because the AD (or AF) of the former is greater than that of the latter. However, mutations A and B may have different relative levels of noise per allele and / or per allele position. For example, mutation A may be a false positive and mutation B may be a true positive as the relative noise levels at these different positions are accounted for. Thus, the model 225 described herein models this noise for proper identification of true positives.
[0069] The probability mass function (PMF) shown in Figure 4 indicates the probability (or likelihood) that a sample from a subject will have a given AD count at a position. Using sequencing data from a sample of a reference individual (e.g., stored in sequence database 210) as a training dataset, processing system 200 trains a model 225 from which a PMF for the reference sample can be derived. The reference individual can be an individual not known or known to be associated with a mutation at a particular variant location, and is sometimes referred to as a healthy individual, although the healthy individual may be associated with a different mutation at another variant location not associated with the model trained for the particular variant location. The PMF models the expected average AD count per allele per position in normal tissue (e.g., of the reference individual), λ p and r, which models the expected variation (e.g., variability) in this AD count. p In other words, λ p and / or r p represents the baseline level of noise, per allele, per position, in the sequencing data for normal tissue. In some embodiments, the sequencing data of the reference individual can be stratified into different read stratums such that multiple models 225 are trained, each corresponding to a particular read stratum. Each model corresponding to a read stratum uses a different λ p and different r p and
[0070] To further illustrate using the example of Figure 4, a sample from a reference individual may be irepresents a subset of the human population modeled by i, where i is the index of healthy individuals in the training set. For purposes of example, assuming that model 225 has already been trained, the PMF generated by model 225 visually illustrates the likelihood of measured AD for each mutation, thus providing an indication of which are true positives and which are false positives. The exemplary PMF on the left side of FIG. 4 associated with mutation A indicates that there is an approximately 20% probability that a first sample has an AD count of 10 for a mutation at position 4. Additionally, the exemplary PDF on the right side associated with mutation B indicates that there is an approximately 1% probability that a second sample has an AD count of 1 for a mutation at position 3 (Note: the PMFs in FIG. 4 are not precisely to scale). Thus, the noise ratios corresponding to these probabilities in the PMF indicate that mutation A is more likely to occur than mutation B, even though mutation B has lower AD and AF. Thus, in this example, mutation B can be a true positive, and mutation A can be a false positive. Thus, processing system 200 can use model 225 to perform improved variant calling by distinguishing true positives from false positives with a more accurate rate and provide additional numerical confidence regarding these likelihoods.
[0071] 5A illustrates the dependencies between parameters and sub-models of a Bayesian hierarchical model 225 for determining true single nucleotide variants, according to various embodiments. The Bayesian hierarchical model illustrated in FIG. 5A can correspond to a particular read stratum of variant reads. The parameters of the model can be stored in a parameter database 230. In the example illustrated in FIG. 5A,
[0072]
number
[0073] represents the vector of weights assigned to each mixture component.
[0074]
number
[0075] takes values within a simplex in K dimensions and can be learned or updated via posterior sampling during training. For such training, a uniform prior on the simplex can be given. The mixture component to which position p belongs is then sampled using one or more different multinomial distributions for the latent variable z p It can be modeled by
[0076]
number
[0077] latent variable z p , a vector of mixture components
[0078]
number
[0079] Together, α, β, and β allow the model for μ, i.e., a sub-model of the Bayesian hierarchical model 225, to have parameters that "pool" knowledge about noise, i.e., they represent the similarity of noise characteristics across multiple positions. Thus, the positions of sequence reads can be pooled or grouped by the model into latent classes. Advantageously, samples from any of these "pooled" positions can also help train these shared parameters. The advantage of this is that the processing system 200 can determine a model of noise in healthy samples even when there is little to no direct evidence that alternative alleles have previously been observed for a given position (e.g., in healthy tissue samples used to train the model).
[0080] Covariate x pThe predictor (e.g., predictor) encodes known context information about position p, which can include, but is not limited to, information such as trinucleotide context, mappability, segment duplication, or other information associated with the sequence read. Trinucleotide context can be based on a reference allele and can be assigned a numerical (e.g., integer) representation. For example, "AAA" is assigned 1, "ACA" is assigned 2, "AGA" is assigned 3, etc. Mappability represents the level of uniqueness of the alignment of a read to a particular target region of the genome. For example, mappability is calculated as the reciprocal of the number of positions to which a sequence read uniquely maps. Segment duplication corresponds to long nucleic acid sequences (e.g., having a length greater than about 1000 base pairs) that are nearly identical (e.g., greater than 90% matching) and occur in multiple locations within the genome as a result of natural duplication events (e.g., not associated with cancer or disease).
[0081] The expected mean AD frequency of an SNV at position p is determined by the parameter μ p In some embodiments, the parameter μ p is the average AD count per sequencing depth
[0082]
number
[0083] Since SNV is an example of a variant, the parameter μ p For simplicity herein, μ p and y p The term μ refers to a location-specific sub-model of the Bayesian hierarchical model 225. In some embodiments, μ p is the shape parameter
[0084]
number
[0085] and rate parameters
[0086]
number
[0087] The eigenvalues are modeled as gamma-distributed random variables with
[0088]
number
[0089] In some embodiments, μ p Other functions can be used to represent
[0090]
number
[0091] and logarithmic standard deviation
[0092]
number
[0093] and log-normal distribution, Weibull distribution, power law, exponentially-modulated power law, or mixtures of the foregoing.
[0094]
number
[0095] is sometimes expressed as the variability parameter r in the distribution. p This can be an example of.
[0096] The variance of the distribution is the mean variant frequency μ p and the variability parameter r p For example, in the case of a gamma distribution, the variance v p teeth,
[0097]
number
[0098] Lambda λ p can be the average variant count, which is a function of the sequencing depth
[0099]
number
[0100] μ multiplied by p Also, lambda λ p is as follows:
[0101]
number
[0102] By the shape parameter
[0103]
number
[0104] and rate parameters
[0105]
number
[0106] can be related to.
[0107] In the example illustrated in FIG. 5A, the shape and rate parameters are each determined by the covariate xp and the latent variable z p However, in some embodiments, the dependency may be different based on different degrees of information pooling during training. For example, the model may alternatively
[0108]
number
[0109] can be structured so that it depends on the latent variables but not on the covariates. The distribution of AD counts of SNVs at position p in a human population sample i (of healthy individuals) is the random variable
[0110]
number
[0111] It is modeled by the random variable
[0112]
number
[0113] In some embodiments, the distribution is determined by the sequencing depth of the sample at a location.
[0114]
number
[0115] is given, it is a Poisson distribution.
[0116]
number
[0117] In some embodiments, the other function is:
[0118]
number
[0119] Examples include, but are not limited to, the negative binomial, the Conway-Maxwell-Poisson distribution, the Zeta distribution, and the zero-excess Poisson distribution. For example, the random variable
[0120]
number
[0121] can be modeled by a negative binomial distribution.
[0122]
number
[0123] Mean variant frequency μ p , average variant count
[0124]
number
[0125] , and the variability parameter r p are random variables with variant counts.
[0126]
number
[0127] Since it affects the distribution of
[0128] 5B illustrates the dependencies between parameters and sub-models of a Bayesian hierarchical model for determining true insertions or deletions, according to various embodiments. In contrast to the SNV model illustrated in FIG. 5A, the model for indels illustrated in FIG. 5B includes different levels of hierarchy. The covariates x p Lattice variables σ encode known features at position p and can include, for example, distance to homopolymers, distance to RepeatMasker repeats, or other information previously associated with observed sequence reads.
[0129]
number
[0130] is a vector
[0131]
number
[0132] The indel length distribution can be modeled by a Dirichlet distribution based on the parameters of
[0133]
number
[0134] are locations that share the same covariate values
[0135]
number
[0136] Thus, for example, a latent variable can represent information such as that homopolymer indels occur at positions 1, 2, 3, etc. base pairs from the anchor position, and that trinucleotide indels occur at positions 3, 6, 9, etc. from the anchor position.
[0137] The expected mean total indel frequency at position p is given by the distribution μ p In some embodiments, the parameter μ p is the average indel count per sequencing depth
[0138]
number
[0139] Since indels are an example of variants, the parameter μ p In some embodiments, the distribution is based on covariates and has a shape parameter
[0140]
number
[0141] and rate parameters
[0142]
number
[0143] and has a gamma distribution with
[0144]
number
[0145] In some embodiments, μ pOther functions can be used to represent , examples include, but are not limited to, negative binomial, Conway-Maxwell-Poisson distribution, Zeta distribution, and zero-excess Poisson.
[0146]
number
[0147] is sometimes expressed as the variability parameter r in the distribution. p This can be an example of.
[0148] The variance of the distribution is the mean variant frequency μ p and the variability parameter r p For example, in the case of a gamma distribution, the variance v p teeth,
[0149]
number
[0150] It can be determined by
[0151] Lambda λ p can be the average variant count, which is a function of the sequencing depth
[0152]
number
[0153] μ multiplied by p Also, lambda λ p is as follows:
[0154]
number
[0155] By the shape parameter
[0156]
number
[0157] and rate parameters
[0158]
number
[0159] can be related to.
[0160] The observed indel at position p in human population sample i (healthy individuals) is given by the distribution
[0161]
number
[0162] It is modeled by the random variable
[0163]
number
[0164] , which may also be referred to as the variant count or observed variant count. As in the example of FIG. 5A, in some embodiments, the distribution of indel intensities is proportional to the sequencing depth of the sample at position
[0165]
number
[0166] is the Poisson distribution given by
[0167]
number
[0168] In some embodiments,
[0169]
number
[0170] Other functions can be used to represent , examples include, but are not limited to, negative binomial, Conway-Maxwell-Poisson distribution, Zeta distribution, and zero-excessive Poisson. For example, in some examples, the random variable
[0171]
number
[0172] is modeled by a negative binomial distribution.
[0173]
number
[0174] Mean variant frequency μ p , average variant count
[0175]
number
[0176] , and the variability parameter r p are random variables with variant counts.
[0177]
number
[0178] Since it affects the distribution of
[0179] Because indels can be of variable length, an additional length parameter exists in the indel model that is not present in the model for SNV. As a result, the exemplary model illustrated in Figure 5B has an additional hierarchical level (e.g., another submodel) that is also not present in the SNV model discussed above. The observed count of an indel of length l (e.g., an insertion or deletion of up to 100 or more base pairs) at position p in sample i is determined by the random variable
[0180]
number
[0181] The distribution is modeled by: , which represents the indel distribution under noise conditional on the parameters: given indel intensities of the polynomials of the sample,
[0182]
number
[0183] and distribution of indel lengths at positions
[0184]
number
[0185] It can be said that:
[0186]
number
[0187] In some embodiments,
[0188]
number
[0189] Dirichlet polynomial functions or other types of models can be used to represent .
[0190] By designing the model in this manner, the machine learning engine 220 can decouple learning of indel intensity (i.e., noise rate) from learning of indel length distribution. Independently determining inferences for predictions about whether indels will occur in healthy samples and predictions about the length of indels at positions can improve the sensitivity of the model. For example, the length distribution can be more stable compared to the indel intensity at some positions or regions within the genome, or vice versa.
[0191] 6A-6B illustrate diagrams associated with the Bayesian hierarchical model 225, according to various embodiments. The graph shown in FIG. 6A shows the distribution of noise rates μ p , i.e., the likelihood (or strength) of an SNV or indel for a given position characterized by the model. The continuous distribution depicts the mean variant frequency μ of non-cancer or non-disease mutations (e.g., mutations occurring naturally in healthy tissue) based on training data of observed healthy samples from healthy individuals (e.g., retrieved from sequence database 210). p Although not shown in FIG. 6A, in some embodiments, the shape and rate parameters are determined by the covariate x p or the latent variable z p The graph shown in Figure 6B shows the distribution of AD at a given location for a subject's sample, the sequencing depth d at a given location, and other variables. p Draw a given parameter of the sample such as μ p The discrete probability of sampling is given by the expected mean distribution μ p The predicted true mean AD count for the human population is determined based on the predicted true mean AD count for the human population ....
[0192] 7A is a diagram of an exemplary process for determining parameters by fitting a Bayesian hierarchical model 225, according to various embodiments. To train the model, the machine learning engine 220 iteratively samples from the posterior distribution of expected noise rates (e.g., the graph illustrated in FIG. 6B) for each location in the set of locations. The machine learning engine 220 can use Markov Chain Monte Carlo (MCMC) methods for sampling, such as the Metropolis-Hastings (MH) algorithm, a custom MH algorithm, a Gibbs sampling algorithm, sampling based on Hamiltonian mechanics, and random sampling, among other sampling algorithms. During Bayesian inference training, the parameters are determined by fitting all (or some) parameters and latent variables (e.g.,
[0193]
number
[0194] , z p ,
[0195]
number
[0196]
number
[0197] , μ p ) are sampled from the joint posterior distribution to iteratively update
[0198] In some embodiments, the machine learning engine 220 p Model fitting is performed by storing a collection of μ in a parameter database 230. The model is trained or fitted through posterior sampling, as previously described. In some examples, μ pThe samplings of are stored in a matrix data structure with a row for each location in the set of sampled locations and a column for each sampling from the simultaneous posterior sampling (e.g., of all parameters conditional on the observed data). The number of rows R can be greater than 6 million, and the number of columns for the N repetitions of the sample can be in the thousands. In some embodiments, the row and column designations differ from the embodiment illustrated in FIG. 7A, e.g., each row represents a sampling from the posterior sampling and each column represents a sampled location (e.g., a transpose of the example matrix illustrated in FIG. 7A).
[0199] 7B is a diagram illustrating determining the likelihood of a false positive using parameters from a Bayesian hierarchical model 225, according to various embodiments. The machine learning engine 220 can reduce the R-by-N matrix illustrated in FIG. 7A to the R-by-2 matrix illustrated in FIG. 7B. In some examples, the machine learning engine 220 can reduce the posterior sampling μ p The position-by-position variability parameter r p (e.g., shape parameter) and mean parameter λ p (mean rate parameter λ p Determine various noise parameters, such as the variability parameter r p teeth,
[0200]
number
[0201] where λ p and v p are the μ at the position p are the mean and variance of the sampled values of r. Those skilled in the art will understand the p It will be appreciated that other functions for determining λ may also be used. Different noise parameters may be determined for different lead layers. For example, each lead layer may have a function λ p and r p can have different values of
[0202] The machine learning engine 220 can also perform variability re-estimation of the variability parameters in the reduced matrix given the rate parameters. In some embodiments, following Bayesian training and posterior approximation, the machine learning engine 220 estimates the variability parameters based on a position-wise negative binomial maximum likelihood estimator.
[0203]
number
[0204] Variation re-estimation is performed by retraining for . The rate parameters can remain fixed during retraining. In some embodiments, the machine learning engine 220 re-estimates the variability by re-training the original AD counts of the training data (e.g., based on reference samples, stratified by lead stratum).
[0205]
number
[0206] and
[0207]
number
[0208] ) for each position, the dispersion parameter r' p The machine learning engine 220 determines
[0209]
number
[0210] Determine
[0211]
number
[0212] Those skilled in the art will appreciate that methods such as moment estimators, posterior means, or posterior modes
[0213]
number
[0214] It will be appreciated that other functions for determining .times. ...
[0215] During application of the trained model, the processing system 200:
[0216]
number
[0217] and λ p To determine the function parameterized by
[0218]
number
[0219] and the rate parameter λ p 7B , which can be accessed using the reduced matrix. The function can be used to determine a posterior predictive probability mass function (or probability density function) for a new sample from a subject. Based on the predicted probability of a certain AD count at a given position, the processing system 200 can interpret the site-specific noise rate per position of the sequencing read when detecting a true positive from the sample. Returning to the exemplary use case described with reference to FIG. 4 , the PMFs illustrated for mutations A and B can be determined using parameters from the reduced matrix in FIG. 7B . The posterior predictive probability mass function can be used to determine the probability that a sample for mutation A or B will have an AD count at a certain position.
[0220] The Bayesian hierarchical model and the distributions used to model the various parameters within the Bayesian hierarchical model can be trained separately for different read stratums of variant reads. For example, each read stratum can be trained separately for different read stratums of variant reads.
[0221]
number
[0222] ,
[0223]
number
[0224] , μ p It can have its own Bayesian hierarchical model with its own parameters, such as:
[0225] For more detailed information regarding the training and use of Bayesian hierarchical models to model noise levels in sequencing datasets, U.S. Patent Application Publication No. 2018 / 0129998, filed October 5, 2018, entitled "Site-Specific Noise model for Targeted Sequencing," is incorporated herein by reference for all purposes.
[0226] V. Exemplary Lead Demographics 8A includes a diagram illustrating different categories or tiers of sequence reads, according to various embodiments. As contemplated herein, sequence reads can be associated with different tiers representing reads of different quality levels, whereby the quality levels can be based on variant locations relative to overlapping segments of the sequence read. A higher quality tier corresponds to a lower noise level or a lower error rate, and a lower quality tier corresponds to a higher noise level or a higher error rate.
[0227] It is noted that in a sequence amplification process (e.g., massively parallel sequencing), one or more sequences of a sample (e.g., an individual) can be cleaved into different fragments in a pseudorandom manner. In some cases, not all fragments are ligated with a UMI, and therefore some of the fragments are washed away before the ligated fragments are enriched. Thus, the enriched fragments are at least partially random in each sequencing run. The extent of overlap between different fragments can vary. For example, some enriched fragments can have overlapping regions that can be stitched with other enriched fragments. Some enriched fragments also have enriched complementary sequences (e.g., forward and reverse sequences, positive and negative sequences, top and bottom sequences, 5' to 3' sequences and 3' to 5' sequences), thus generating double-stranded reads for all or part of the entire sequence read. As a result, variant reads at different sequence locations can, in some instances, include complementary and / or overlapping sequence reads to confirm the variant. Thus, each variant read can correspond to different read layer qualities. For example, a variant read at a location where both complementary strands of a fragment are enriched will often be of better quality than another variant read at a second location where only a single fragment is enriched. There is an increased likelihood that a variant read at a location that is not within an overlapping or complementary region is due to noise rather than an actual variant present in the subject's sample.
[0228] 8A illustrates four different examples of read tiers. In some embodiments, sequence reads are separated into read tiers based on potential variant locations of interest within the sequence read relative to overlapping and complementary locations within the sequence read. In other words, sequence reads are classified into one of four read tiers based on whether the potential variant location read is contained within or completely embedded in an overlapping region (i.e., a stitched region) and whether the variant location is contained within or completely embedded in a complementary region (i.e., a double-stranded region, a duplex region).
[0229] For example, in FIG. 8A , potential variant locations are shaded. The first exemplary read stratum 810 includes variant reads that fall alongside both double-stranded (also called "duplex" or "complementary") sequence reads and stitched sequence reads. For example, at least two 5' to 3' sequence reads have an overlapping region and can be stitched together. Similarly, at least two 3' to 5' sequence reads have an overlapping region and can also be stitched together. Within the exemplary first read stratum 810, potential variant locations are located within or completely embedded within the overlapping region or stitched region, and thus the sequence reads include the stitched region. Similarly, at least a portion of the 5' to 3' sequence read and a portion of the 3' to 5' sequence read are complementary to each other, and the potential variant location is located within the complementary region (e.g., the potential variant location is completely embedded in both the top sequence read and the bottom sequence read, in the common region of their overlap). Thus, in addition to containing the stitched region, the sequence read also contains a double-stranded region, and the potential variant read belongs to the first read stratum 810, which represents the double-stranded stitched read.
[0230] In FIG. 8A , the second exemplary read stratum 820 includes variant reads located within portions of the sequence reads that are double-stranded but not stitched. In the second read stratum 820, a portion of the 5' to 3' sequence read and a portion of the 3' to 5' sequence read are complementary to each other, and the potential variant locations are located within the complementary regions. Thus, the sequence reads include double-stranded regions. However, the potential variant locations are not included in either the overlapping or stitched regions of the sequence reads. Specifically, this exemplary stratification is despite the fact that the two 5' to 3' sequence reads can be stitched together because the potential variant locations are not included in the overlapping or stitchable regions. Thus, although the sequence reads include double-stranded regions, the sequence reads do not include stitched regions, and the potential variant reads belong to the second read stratum 820, which represents double-stranded but not stitched reads.
[0231] The third exemplary read stratum 830 includes variant reads located or entirely embedded within single-stranded (e.g., non-duplex) stitched reads. In the third read stratum 830, potential variant locations are contained within overlapping regions of two or more sequence reads, and thus the sequence reads include stitched regions. However, because the sequence reads (such as the two illustrated 5' to 3' sequence reads) do not include complementary regions (e.g., the sequence reads are based only on the 5' to 3' strand and are not supported by complementary 3' to 5' strands), the sequence reads are single-stranded. In some cases (not illustrated), one or more complementary sequence reads (e.g., 3' to 5' sequence reads) are discoverable in the exemplary read stratum 3 but do not include potential variant locations. Thus, the potential variant reads belong to the third read stratum 830, which represents single-stranded but stitched reads.
[0232] 8A , the fourth exemplary read stratum 840 includes variant reads located within single-stranded, unstitched reads. Like the third read stratum 830, the illustrated sequence reads do not include complementary regions containing variant locations (or, in some cases (not illustrated), further include complementary regions, but potential variant locations are not located within or completely embedded within the complementary regions), so the fourth read stratum 840 represents single-stranded reads. Thus, the potential variant locations are not contained within the overlapping regions of the two sequence reads, so the fourth read stratum 840 represents unstitched reads.
[0233] In some embodiments, the sequence reads for a sample can be stratified into four read tiers illustrated in FIG. 8A. In some embodiments, there can be an additional, fifth read tier corresponding to the lowest quality of variant reads, such as single-stranded, unstitched sequence reads that contain potential variant locations near one or more ends of the sequence read. For example, if a single-stranded, unstitched sequence read contains a potential variant location within a predetermined threshold number of bases (e.g., within about 7 bases or within about 30 bases) from either end of the sequence read, the sequence read can be classified into the fifth read tier. In some embodiments, each of the four read tiers illustrated in FIG. 8A can be subdivided into two subtiers: a first, low-quality subtier corresponding to sequence reads that contain potential variant locations near either end of one or more sequence reads, and a second, high-quality subtier that contains potential variant locations at a distance greater than a threshold from one or more ends of the sequence read.
[0234] FIG. 8B illustrates experimental results illustrating different qualities of the read stratum of FIG. 8A , according to various embodiments. A high-quality read stratum corresponds to a low error rate and / or a low noise level. In other words, variant reads (e.g., detected SNVs or indels at potential variant locations) within sequence reads stratified into a high-quality read stratum are more likely to be due to actual mutations in the sample, as opposed to random events (e.g., due to noise), than variant reads within sequence reads stratified into a low-quality read stratum. FIG. 8B is a log10 plot of the average error rate of a reference sample (e.g., a healthy individual) for different read stratums t1-t5. Tier 1 (t1) refers to the first read stratum 810 of FIG. 8A , stratum 2 (t2) refers to the second read stratum 820, and stratum 5 (t5) refers to single-stranded, unstitched read stratums and / or potential variant locations located near either end of the sequence read (e.g., within 7 bases from the end). For layer 1, the logarithm base 10 of the mean error rate mu is approximately -6.3 to -7. In other words, the average error rate for healthy individuals is approximately 10 per sequencing depth. -6.3 From 10 -7The average error rate for variant reads is approximately 1 / 1,000,000 to -1 / 10,000,000. Meanwhile, for tier 4, the log base 10 of that average error rate, mu, is approximately -4.7 to -5.5. In some embodiments, FIG. 8B generally illustrates that the average error rate increases across tiers 1 to 4, from a tier 1 average error rate of approximately 1 / 1,000,000 to -1 / 10,000,000, to a tier 2 average error rate closer to approximately 1 / 1,000,000, to an average error rate of approximately <1 / 1,000,000 in tier 3, and still increasing to a tier 4 error rate of approximately 1 / 100,000. Thus, variant reads detected in tier 4 are approximately 100 times more likely to be erroneous alleles than variant reads detected in tier 1. In other words, the fourth read tier is relatively noisier and more error-prone than the first read tier. In other words, the fourth stratum of reads is relatively noisier and more error-prone than the first stratum of reads. Furthermore, as shown for example in stratum 5, there is no mean error rate μ (or no meaningful mean error rate) because the sequence reads there are of such low quality that they are unimportant or discarded from analysis.
[0235] Sequencing reads can additionally or alternatively be classified into different read tiers by other classification methods. For example, if the variants are SNVs, each read tier can be further subdivided into 12 additional subtiers based on the type of nucleotide substitution (e.g., A>C, A>T, G>C, etc.) (see, e.g., Figure 8H, discussed below). Because there are four nucleotides, and each nucleotide is substituted by a different nucleotide within the SNV, there are a total of 12 different types of SNVs.
[0236] Figures 8C through 8G illustrate experimental results for the read tier of Figure 8A when sequence reads are first stratified into read tiers in the manner described in Figure 8A and further stratified into 12 sub-read tiers based on the type of nucleotide substitution. Specifically, Figures 8C through 8G illustrate the mean error rate mu (μ) and size parameters for a statistical model (e.g., negative binomial) for the error distribution of alternative reads at a position, such that the error distribution at each position (e.g., each point) is conditional on the actual read depth seen in a given sample. The model is stratified by different categories of reads (e.g., tiers); thus, Figure 8C illustrates the results of different types of nucleotide substitutions for the first read tier (i.e., double-stranded stitched reads), Figure 8D illustrates the results of different types of nucleotide substitutions for the second read tier (i.e., double-stranded but unstitched reads), Figure 8E illustrates the results of different types of nucleotide substitutions for the third read tier (i.e., single-stranded stitched reads), and Figure 8F illustrates the results of different types of nucleotide substitutions for the fourth read tier. Figure 8G illustrates results from the fifth read tier, which corresponds to the lowest quality reads and a high error rate beyond the illustrated axis. Note that for Figures 8C-8G, the nucleotide bases A, C, G, and T horizontally across the top of the plot refer to alternative bases, while the nucleotide bases A, C, G, and T vertically along the right side of the plot refer to reference bases.
[0237] Referring to Figure 8C, 12 different distributions of representative variant frequencies within different SNVs illustrate that the first read tier can be further divided into 12 subtiers to further refine the noise model. Rows correspond to original nucleotides, and columns correspond to changed nucleotides. For example, the cells in the third row and first column can correspond to a G to A SNV. Experiments illustrate that the subread tier for C to A (i.e., second row, first column), whose μ distribution is centered in the range of -7 to -8 on a base 10 logarithmic scale, is likely to be less noisy than the subread tier for T to C (i.e., fourth row, second column), whose μ distribution spans the range of -5 to -7 on a base 10 logarithmic scale.
[0238] Comparing the differences between Figures 8C through 8G, the distribution of the mean error rate μ on a logarithmic scale illustrates a general trend in which μ shifts toward zero (i.e., μ becomes larger) and eventually exceeds zero (i.e., in Figure 8G) as the lead layer changes from the first lead layer to the fifth lead layer. For example, focusing on the sub-lead layer for T to G (i.e., the fourth row, third column), the distribution of μ on a logarithmic scale shifts from between -6 and -7 in the first lead layer to between -4 and -5 in the fourth lead layer. Thus, Figures 8C through 8G demonstrate that as the lead layer becomes noisier, the mean error rate μ also increases.
[0239]
[0041] Referring now to Figure 8H, experimental results are shown illustrating different mean error rates μ at particular SNV nucleotide substitutions obtained across read strata t1-t4 of Figure 8A, according to various embodiments described herein. Specifically, Figure 8B is a plot of the log10 of the mean error rate μ of a reference sample (e.g., healthy individual) on a logarithmic scale for different SNVs observed across read strata t1-t4. It is noted that the statistical models for error distribution described herein can be stratified by different categories of reads (strata) and / or by different SNVs as further illustrated in Figure 8H.
[0240] VI. Exemplary Data Processing with Stratified Reads 9 is a flowchart illustrating a process for analyzing a sample's DNA data set using a stratified noise model according to various embodiments. The process can be used to process a DNA sequencing data set of a sample, such as an individual's sample containing cfDNA, to generate a quality score that represents the likelihood that the individual has a variant at a potential variant location. The higher the quality score determined by the process, the more likely the variant read is the result of a real mutation rather than noise.
[0241] In step 910, the processing system can access a DNA sequencing dataset generated by DNA sequencing. For example, the DNA sequencing can be a type of massively parallel DNA sequencing, such as next-generation sequencing (NGS). The DNA sequencing dataset includes a plurality of processed sequence reads that include variant locations of interest (e.g., specific genetic locations within a DNA sequence). At least some of the processed sequence reads can be generated from collapsing and stitching raw sequence reads within the DNA sequencing, such as those generated by the process described in FIG. 3. For example, a typical NGS run can generate millions or even billions of sequence reads. Some of the raw sequence reads can include genetic loci that include the variant locations of interest. The raw sequence reads can be processed by collapsing and stitching to generate processed sequence reads. While DNA sequencing is described in this example, it is noted that RNA sequencing can also be implemented for analysis herein.
[0242] The processed sequence reads containing variant locations of interest can be of different base pair lengths and of different extents of overlap and / or complementarity. In step 920, the processing system can stratify the multiple processed sequence reads into different read tiers. The different read tiers can be stratified based on the quality of the sequence reads. For example, the processed sequence reads can be stratified based on whether the variant locations are included in overlapping and / or complementary regions, as discussed in connection with FIG. 8A . Other ways of stratifying the processed sequence reads are also possible. For example, the processed sequence reads can also be stratified based on the type of nucleotide substitution, whether the variant locations are near the ends of the sequence, etc. In some embodiments, the different read tiers include at least four read tiers. In some examples, the four lead layers are: (1) a double-stranded, sewn lead layer; (2) a double-stranded, unsewn lead layer; (3) a single-stranded, sewn lead layer; and (4) a single-stranded, unsewn lead layer.
[0243] In step 930, the processing system can determine a stratified sequencing depth at the variant location for each read stratum. For each read stratum, the stratified sequencing depth can be the sequencing depth of the sequence reads stratified into the read stratum. In other words, the stratified sequencing depth can be the total number of sequence reads stratified into the read stratum. The processing system can also determine an actual variant count for each read stratum. For example, for a read stratum, the majority of sequence reads may not contain an actual variant (whether it is an SNV or an indel) at the variant location. In some cases, only a small number of sequence reads contain an actual variant at the variant location. The stratified variant count can be the total number of actual variant counts for a particular read stratum.
[0244] In step 940, the processing system can determine, for each read stratum, one or more noise parameters conditional on the stratified sequencing depth of the read stratum. The noise parameters can be parameters of a noise model specific to the read stratum. For example, the processing system can include multiple stratified noise models, each specific to a read stratum. The stratified noise models (or some of them) can correspond to the Bayesian hierarchical models described in FIGS. 5A through 7B. In other words, in some embodiments, each read stratum has its own Bayesian hierarchical model. Each noise model for a read stratum can be trained using a different training set of DNA sequencing samples. For example, DNA sequencing datasets of multiple reference individuals, such as healthy individuals, can be collected. The processed sequence reads of the reference individuals' datasets can be stratified by read stratum. The stratified processed sequence reads for each read stratum can be used as a stratified training set of DNA sequencing samples to train the stratified noise model for the read stratum. Various distributions for the stratified noise model, such as the gamma and Poisson distributions discussed in connection with FIGS. 5A and 5B, can be determined based on the stratified training set.
[0245] The probability distribution of stratified variant counts for each read stratum can be modeled by a noise distribution. The probability distribution of stratified variant counts can depend on the type of distribution used and one or more parameters that define the noise distribution. For example, in the case of the Bayesian hierarchical model discussed, the distribution of stratified variant counts can correspond to a posterior distribution conditional on two parameters. The parameters can be stratified mean variant counts conditional on stratified sequencing depth and variability parameters. Each of the parameters can further correspond to one or more prior distributions that affect the parameter. For example, the stratified mean variant counts conditional on stratified sequencing depth can be modeled by a gamma distribution. Because the prior distribution can describe the distribution of the parameter, the prior distribution can also be referred to as a parameter distribution.
[0246] For each lead stratum, the processing system can determine one or more noise parameters conditional on the stratified sequencing depth by inputting the stratified sequencing depth obtained from the subject's dataset into a trained noise model. For example, the trained noise model can access a parameter distribution (e.g., a prior distribution) specific to the lead stratum. The parameter distribution can be formed based on a stratified training set of reference individuals and can describe the distribution of the stratified training set. The trained noise model can use the parameter distribution to determine noise parameters conditional on the stratified sequencing depth corresponding to the lead stratum.
[0247] Although a Bayesian hierarchical model is used as an example of a noise model, in various embodiments, different types of trained machine learning models can be used as the noise model, and depending on the model used, different noise parameters can be used to model the noise distribution.
[0248] In step 950, the processing system can generate an output for a read stratum-specific noise model based on one or more noise parameters conditional on the stratified sequencing depth of the read stratum. Output generation can be repeated for different read stratums. Depending on the embodiment, different types of output can be generated. For example, in some embodiments, each stratified noise model performs no further calculations after the noise parameters are determined. The output of the noise model can be one or more noise parameters conditional on the stratified sequencing depth determined for each stratum. In cases where a negative binomial distribution is used as the noise distribution to model the stratified variant count, the output of the noise model can be a stratified mean variant count conditional on the stratified sequencing depth and the variability parameter. In some embodiments, after determining the noise parameters, each stratified noise model can generate a posterior distribution. In such embodiments, the output of the read stratum-specific noise model can be the likelihood, due to noise, that the variant count of the read stratum for subsequently observed data is greater than or equal to the total variant count observed in the subject's individual DNA dataset. Other suitable outputs are possible.
[0249] In step 960, the processing system can combine the generated noise model outputs to produce a combined result. The combined result can be a representation of the overall processing result of the subject's individual DNA sequencing dataset. The combined result can take any suitable form. In some embodiments, the combined result can include a likelihood that the total variant count for the subsequently observed data is greater than or equal to the total variant count observed in the multiple processed sequence reads due to noise. Stated differently, the likelihood can represent the likelihood that an event is as extreme as or more extreme than the total variant count observed in the multiple processed sequence reads of the subject's individual DNA dataset. In some cases, the likelihood can correspond to a p-value used in the null hypothesis. How the stratified noise model outputs can be combined to produce a combined result can depend on different embodiments. In some embodiments, the moment matching method discussed in detail in FIG. 10 can be used. In some embodiments, the integration method (e.g., the summation method) discussed in detail in FIGS. 11A and 11B can be used.
[0250] In step 970, the processing system may determine a quality score for the combined results. In some examples, the combined results, such as in the form of a likelihood P (e.g., a p-value), may be converted to a Phred-scale quality score, where Q=-10 log 10P. For example, a Phred quality score of 20 indicates a 1 / 100 chance of P = 1 in 1000 chance of an incorrect variant call, and a Phred quality score of 60 indicates a 1 / 1,000,000 chance of P = 1 in 1,000,000 chance of an incorrect variant call. Thus, a higher Phred quality score corresponds to a greater confidence in the detection of an actual mutation. The quality score can be used to distinguish true positives from false positives. In some embodiments, in response to a quality score higher than a predetermined threshold, the processing system can indicate that the individual is statistically likely to have a mutation at the variant location.
[0251] VII. Moment Matching for Combining Stratified Outputs FIG. 10 is a flowchart depicting a process for combining the outputs of stratified noise models for different read stratums using moment matching, according to various embodiments. The process depicted in FIG. 10 may correspond to steps 950 and / or 960 of FIG. 9. In step 1010, the processing system can combine the mean variant count and the variance of the variant counts from each noise model output to produce an overall mean variant count and an overall variability parameter. The output of each noise model can take the form of a noise parameter conditional on the stratified sequencing reads. The processing system can evaluate the overall likelihood (e.g., an overall p-value) of the total variant count across different read stratums given the total observed sequencing reads by first matching the individual moments of each read stratum to generate overall moments. The processing system can use the overall moments to model an overall distribution conditional on the total observed sequencing reads. The individual noise distributions for the read stratums can be negative binomial distributions. Similarly, the overall noise distribution across multiple read stratums can also be negative binomial distributions matching the moments of the individual read stratums.
[0252] Step 1010 can include several substeps. For each read stratum, the processing system can determine a stratified sequencing depth. The first and second moments of the noise distribution for each stratum can be used as noise parameters to define the noise distribution. In step 1012, based on the stratified sequencing depth, the processing system can determine the first moment (e.g., average variant count) of each read stratum. For example, in the case of the Bayesian hierarchical model discussed above, the variant frequency for a particular read stratum can be determined using the shape parameter
[0253]
number
[0254] and rate parameters
[0255]
number
[0256] The eigenvalues can be modeled as gamma-distributed random variables with
[0257]
number
[0258] Each lead stratum can have its own shape and rate parameters determined based on the reference sample dataset, and therefore the variant frequencies of each lead stratum, conditional on stratified sequencing depth, can be different.
[0259] The processing system calculates the first moment, or stratified mean variant count, λ, for each stratum by multiplying the variant frequency by the stratified sequencing depth. p can be determined.
[0260]
number
[0261] In step 1014, the processing system may also determine the second moment, or variance, of each lead stratum. In cases where the Bayesian hierarchical model has gamma-distributed variant frequencies, the variance of each lead stratum is the mean variant count, or λ p and the variability parameter, i.e., r p For example, the variance, i.e., v p teeth,
[0262]
number
[0263] It can be determined by
[0264] In step 1016, the processing system can determine the overall mean variant count (global first moment) and the overall variance (global second moment) by moment matching. In some cases, the processing system can perform moment matching by summing the moments for different read stratums to obtain an overall moment. For example, the overall mean variant count across all read stratums, conditional on total sequencing depth, can be calculated as:
[0265]
number
[0266] Similarly, the overall variance across all lead layers can be determined by summing the variances of each lead layer.
[0267]
number
[0268] The processing system can model the likelihood of the overall observed variant counts conditional on the total sequencing depth by an overall noise distribution, which is the overall mean λ all and the overall variability parameter r all The overall dispersion parameters can be determined by the overall mean and the overall variance.
[0269]
number
[0270] In step 1020, the processing system can determine an overall likelihood using the overall noise distribution modeled by the global first moment and the global second moment. For example, the random variable
[0271]
number
[0272] is modeled by a negative binomial distribution.
[0273]
number
[0274] A random variable representing the likelihood that an event is as extreme as or more extreme than the observed overall variant count, conditional on the total sequencing depth
[0275]
number
[0276] can be the combined result of the processing system. In some cases, the random variables
[0277]
number
[0278] can be used as a p-value to verify or reject the null hypothesis that the variant reads are due to random events (e.g., noise). The processing system also
[0279]
number
[0280] Negative binomial tail probability can be applied to obtain a p-value based on the Phred scale quality score.
[0281] VIII. Integration Methods for Combining Stratified Outputs FIG. 11A is a flowchart depicting a process for combining the output of stratified noise models for different read stratums using an integration approach to combine the likelihood of each read stratum, according to some embodiments. The process depicted in FIG. 11A may correspond to steps 950 and / or 960 of FIG. 9. In the process illustrated in FIG. 11A, the processing system is intended to determine a p-value as the combined result of the overall system. The p-value can then represent the likelihood that the total variant count for the observed data is greater than or equal to the total variant count observed in the actual data due to noise. Stated differently, the p-value can represent that an event occurs as extreme as or more extreme than the observed sequence read of an individual subject under the null hypothesis. In some cases, an event can be more extreme (i.e., less likely) than the observed variant count of an individual subject when the overall variant count of the event is greater than the observed variant count. Similarly, for a read stratum, an event can be more extreme than the stratified observed variant count when the event has a stratified variant count higher than the stratified observed variant count. Because variant reads are often unusual, events with variant counts are more extreme. The likelihood of observing a variant read at a variant location is often significantly lower than the likelihood of observing a non-variant read at a variant location.
[0282] FIG. 11B illustrates counting of more extreme events in a multidimensional space according to various embodiments. For simplicity, two dimensions are shown, but processing systems according to various embodiments can handle higher dimensions using the principles illustrated in FIG. 11B. Each of the two dimensions in FIG. 11B can represent two lead layers. For example, the first lead layer represents a double-stranded lead layer, and the second lead layer represents a single-stranded lead layer.
[0283] Focusing first only on the read stratum for double-stranded variant count (x-axis), in the example, the observed stratified variant count from an individual subject is 2. For the same read stratum, because variant reads at potential variant locations are less likely than non-variant reads at potential variant locations, an event with a stratified variant count of 3 is less likely (more extreme) than the actual observed stratified variant count. Similarly, another event with a stratified variant count of 4 is even less likely than the actual observed stratified variant count. In other words, the combination of less likely (more extreme) events occupies a larger space than the observed variant count, ranging all the way to infinity. Conversely, an event with a stratified variant count of 1 or 0 is more likely than the actual observed stratified variant count of 2.
[0284] Considering both read tiers, there may be different combinations of observed stratified variant counts that can be assumed to be equally or nearly equally likely. During NGS sample preparation, the subject's individual nucleic acid sequence may be cleaved in a partially random manner. As a result, some of the processed sequence reads may not contain complementary sequence reads. Therefore, some of the processed sequence reads may be single-stranded sequence reads. In other words, for the same nucleic acid sequence sample, different NGS runs may produce different combinations of sequence reads in different read tiers. The stratified variant count of the first read tier may be equal to the stratified variant count of the second read tier based on some ratio. In some embodiments, this ratio is modeled as a predetermined value. For example, one double-stranded variant count may be considered to be equal to two single-stranded variant counts, although in some embodiments, numbers other than two may be used.
[0285] Based on the observed stratified variant counts of the different read strata, the coordinates in the graph shown in FIG. 11B can be divided into points representing events that are as extreme as or more extreme than the observed variant counts of the different read strata and points representing events that are less extreme than the observed variant counts. For example, assume that the actual observed variant count of an individual subject is (1,2) (i.e., one double-stranded variant count and two single-stranded variant counts). Coordinates (0,4) and (2,0) are assumed to be equivalent to the combination (1,2) and can represent events that are as extreme as the observed data. All coordinates that fall outside the sampling boundary by (1,2), (0,4), and (2,0) can be classified as more extreme than the observed data. For example, coordinate (3,3) can be considered a more extreme case. All coordinates that fall within the boundary and closer to the origin can be considered less extreme (more likely) than the observed data. For example, coordinates (1,1), (0,2), (1,0), etc. can be classified as less extreme than the observed data.
[0286] The combined results of the processing system can be in the form of a p-value, which represents the likelihood of an event as extreme as or more extreme than the observed data. The processing system can integrate by summing the probabilities corresponding to all coordinates representing an event as extreme as or more extreme than the observed data to determine the p-value. However, because coordinates may include points that can ultimately approach infinity, the processing system can instead calculate the statistical complement of the p-value. In other words, the processing system can sum the probabilities corresponding to all coordinates representing an event less extreme than the observed data to determine the complement of the p-value. The processing system can then determine the p-value by subtracting the complement from 1.0. In some embodiments, because adding floating-point numbers on a computer can be numerically unstable, the processing system can use probabilities on a logarithmic scale for numerical stability.
[0287] Returning to FIG. 11A , a processing unit, according to some embodiments, can determine p-values based on the process illustrated in the figure. In step 1110, the processing unit can determine, within each lead stratum, possible events that are more likely than the lead stratum's observed stratified variant count. These events can define a multidimensional box whose coordinates can be associated with a greater likelihood than the observed data. The likelihood can be on a logarithmic scale. In step 1120, the processing unit can determine, for one of the multiple lead stratums, whether each of the stratified variant count combinations corresponds to a higher or lower likelihood than the observed data. In step 1130, the processing unit can identify possible event combinations associated with a higher likelihood of occurrence than the lead stratum's observed stratified variant count. The processing unit can repeat steps 1120 and 1130 for each lead stratum. In step 1140, the processing unit can sum the probabilities of the identified combinations to determine a statistical complement. In step 1150, the processing unit may subtract the statistical complement from 1.0 to determine an overall p-value, which may be the combined result corresponding to step 960 in FIG.
[0288] Other ways of determining the overall p-value are also possible. For example, tail probability techniques can be used. In some embodiments, the integration method can be replaced by one or more machine learning models. For example, a random forest regression model can be trained to determine a Phred-scale quality score or p-value from a set of training sample data. The integration process described in FIG. 11A can be used to generate multiple training set samples. The training set samples can be used to train a machine learning model, and the model can then be used to determine a quality score.
[0289] IX. Experimental Results 12A and 12B illustrate plots of observed quality scores against default quality scores in an experiment performed using a process according to some embodiments. In FIG. 12A, a calibrated set of simulated data for a sample from a particular individual #1355 is used and analyzed using the stratum-level noise model process described in FIGS. 9 and 11A. The individual has a variant symbolized as "chr6" found in the sample. The data simulates random events of variant reads. Each simulated event can have several variant reads that can be stratified into one or more read stratums. The x-axis represents the actual quality score values of the simulated events calculated using the simulated data. The y-axis represents the observed quality score values determined by the process described in FIG. 11A. The results illustrate that, barring some discretization, the observed quality scores mostly fall on the diagonal. This indicates that the process described in FIG. 11A successfully determines the likelihood of possible events that are as extreme as or more extreme than the observed data.
[0290] In Figure 12B, the actual dataset for individual #1355 is used and analyzed using the stratum-level noise model process described in Figures 9 and 11A. The dataset can include data for various potential variant locations. Each point corresponds to a potential variant location, and a quality score is determined based on the process described in Figure 11A. In some cases, many of the locations in the actual dataset may not have variant counts found in any read stratum, resulting in several zero-quality score points at the beginning. Therefore, the Phred-scale quality score is equal to zero for all of those locations. Because individuals largely lack mutations at the majority of sequence locations, the majority of the remaining points fall mostly along the diagonal. Stated differently, the variant counts at those locations along the diagonal may be largely due to noise. Outliers with observed quality scores significantly higher than the default quality score (e.g., points at a default Phred score of approximately 55 in Figure 12B) can indicate possible non-noise events that may be flagged for further evaluation.
[0291] 13A and 13B illustrate experimental results comparing quality score results using read tiers with quality scores using a noise model that does not separate read tiers for sequence reads or variant reads. The y-axis in FIGS. 13A and 13B represents quality scores determined using a method for stratifying sequence reads into different read tiers according to various embodiments. The x-axis in FIGS. 13A and 13B represents quality scores determined using a similar noise model, but this noise model does not distinguish sequence reads by read tier. FIG. 13A illustrates experimental results determining quality scores using the integration method described in FIG. 11A. The results illustrate that for sequence reads containing duplex reads (e.g., dark points marked "true"), data points are shifted upward compared to sequence reads containing only simplex reads (e.g., bright points marked "false"). This indicates that the read layer noise model improves the overall quality score for sequence reads that include double-stranded reads, because double-stranded reads often contain more evidence than single-stranded reads. Figure 13B illustrates experimental results for determining quality scores using the moment matching method described in Figure 10. Similar to Figure 13A, in this case, the moment matching method also improves the quality score for sequence reads that include double-stranded reads.
[0292] X. Variant Identification FIG. 14 is a flowchart depicting a process for identifying potential variant locations for an individual, according to various embodiments. In step 1410, the system may receive a DNA sample from the individual. In step 1420, the system may perform DNA sequencing to generate processed sequence reads. In step 1430, the system may assign the processed sequence reads with different variant locations. In step 1440, for each variant location, the system may stratify the processed sequence reads assigned to the variant location into multiple read strata. In step 1450, the system may determine quality scores of likelihood at different variant locations. Each quality score may be determined based on the process described above using a noise model to stratify read strata. In step 1460, the system may identify variant locations with quality scores above a predetermined threshold. These variant locations may be flagged for further investigation of potential mutations or potential diagnoses.
[0293] In step 1470, the system can generate a diagnosis of disease based on the identified variant locations. In some embodiments, variants or mutations that are indicative of some cancers and / or can serve as biomarkers for some therapeutic agents include ACVR1B, AKT3, AMER1, APC, ARID1A, ARID1B, ARID2, ASXL1, ASXL2, ATM, ATR, BAP1, BCL2, BCL6, BCORL1, BCR, BLM, BRAF, BRCA1, BTG1, CASP8, CBL, CCND3, CCNE1, CD74, CDC73, CDK12, CDKN2A, CHD2, CJD2, CREBBP, CSF1R, CTCF, CTNNB1, DICER1, DNAJB1, DNMT1, DNMT3A, DNMT3B, DOT1L, EED, EGFR, EIF1AX, EP300, EPHA3, EP HA5, EPHB1, ERBB2, ERBB4, ERCC2, ERCC3, ERCC4, ESR1, FAM46C, FANCA, FANCC, FANCD2, FANCE, FAT1, FBXW7, FGFR3, FL CN, FLT1, FOXO1, FUBP1, FYN, GATA3, GPR124, GRIN2A, GRM3, H3F3A, HIST1H1C, IDH1, IDH2, IKZF1, IL7R, INPP4B, IRF4 , IRS1, IRS2, JAK2, KAT6A, KDM6A, KEAP1, KIF5B, KIT, KLF4, KLH6, KMT2C, KRAS, LMAP1, LRP1B, LZTR1, MAP3K1, MCL1, M GA, MSH2, MSH6, MST1R, MTOR, MYD88, NPM1, NRAS, NTRK1, NTRK2, NUP93, NUTM1, PAX3, PAX8, PBRM1, PGR, PHOX2B, PIK3C A, POLE, PTCH1, PTEN, PTPN11, PTPRT, RAD21, RAF1, RANBP2, RB1, REL, RFWD2, RHOA, RPTOR, RUNX1, RUNX1T1, SDHA, SHQ1, SLIT2, SMAD4, SMARCA4, SMARCD1, SNCAIP, SOCS1, SPEN, SPTA1, SUZ12, TET1, TET2, TGFBR, and TNFRSF14.In some embodiments, the cancer immunotherapy may target OX40, LAG3, and / or ICOS.
[0294] In step 1480, treatment for the disease may be provided. Prior to providing the treatment, a companion diagnostic operation may also be performed. The companion diagnostic operation may identify one or more criteria, including variants or mutations, using the processes described herein. Providing the treatment may take the form of having or recommending to a medical professional administer a particular dose of a medication to the patient.
[0295] For example, the systems and methods described herein can be used to detect variants or mutations that are biomarkers for cancer therapies, such as some immunotherapies and targeted therapeutics. Such therapies can include, for example, immunoglobulins, proteins, peptides, small molecules, nanoparticles, or nucleic acids. In some embodiments, the therapeutic comprises an antibody, or a functional fragment thereof. In some embodiments, the antibody is selected from the group consisting of Rituxan® (rituximab), Herceptin® (trastuzumab), Erbitux® (cetuximab), Vectibix® (panitumumab), Arzerra® (ofatumumab), Benlysta® (belimumab), Yervoy® (ipilimumab), Perjeta® (pertuzumab), Tremelimumab®, Opdivo® (nivolumab), Dacetuzu®, and others. mab®, Urelumab®, Tecentriq® (atezolizumab, MPDL3280A), Lambrolizumab®, Blinatumomab®, CT-011, Keytruda® (pembrolizumab, MK-3475), BMS-936559, MED14736, MSB0010718C, Imfinzi® (durvalumab), Bavencio® (avelumab), and margetuximab (MGAH22).
[0296] In some embodiments, immunotherapy and targeted therapy include PD-1 inhibition, PD-L1 inhibition, or CTL-4 inhibition. PD-1 inhibition targets programmed death receptors on T cells and other immune cells. Examples of PD-1 inhibitor immunotherapies include pembrolizumab, Keytruda, nivolumab, Opdivo, cemiplimab, and Libtayo®. PD-L1 inhibition targets programmed death receptor ligands expressed by tumor cells and regulatory immune cells. Examples of PD-L1 inhibitor immunotherapies include atezolizumab, Tecentriq, avelumab, Bavencio, durvalumab, and Imfinzi. CTL-4 inhibition targets T cell activation. Examples of CTL-4 inhibitor immunotherapies include ipilimumab and Yervoy.
[0297] For non-small cell lung cancer indications, variants or mutations that can be biomarkers for immunotherapy treatment include EGFR exon 19 deletions & EGFR exon 21 L858R alterations (e.g., for treatments such as Gilotrif® (afatinib), Iressa® (gefitinib), Tagrisso® (osimertinib), or Tarceva® (erlotinib)), EGFR exon 20 T790M alterations (e.g., may be treated with Tagrisso® (osimertinib)), ALK rearrangements (e.g., may be treated with Alecensa® (alectinib), Xalkori® (crizotinib), or Zykadia® (ceritinib)); BRAF These may include V600E (which may be treated, for example, with Tafinlar® (dabrafenib) in combination with Mekinist® (trametinib)), single nucleotide variants (SNVs) and indels that result in MET exon 14 skipping (which may be treated, for example, with Tabrecta™ (capmatinib)).
[0298] For melanoma indications, variants or mutations that can be biomarkers for immunotherapy treatment include BRAF V600E (which may be treated, for example, with Tafinlar® (dabrafenib) or Zelboraf® (vemurafenib)), BRAF V600E or V600K (which may be treated, for example, with Mekinist® (trametinib) or Cotellic® (cobimetinib) in combination with Zelboraf® (vemurafenib)).
[0299] For breast cancer indications, variants or mutations that may be biomarkers for immunotherapy treatment include ERBB2 (HER2) amplification (e.g., which may be treated with Herceptin® (trastuzumab), Kadcyla® (ado-trastuzumab emtansine), or Perjeta® (pertuzumab)), PIK3CA alterations (e.g., which may be treated with Piqray® (alpelisib)).
[0300] Variants or mutations that may be biomarkers for immunotherapy treatment in colorectal cancer indications include KRAS wild-type (lack of mutations in codons 12 and 13) (e.g., may be treated with Erbitux® (cetuximab)), KRAS wild-type (lack of mutations in exons 2, 3, and 4), and NRAS wild-type (lack of mutations in exons 2, 3, and 4) (e.g., may be treated with Vectibix® (panitumumab)).
[0301] For ovarian cancer indications, variants or mutations that may be biomarkers for immunotherapy treatment may include BRCA1 / 2 alterations (e.g., which may be treated with Lynparza® (olaparib) or Rubraca® (rucaparib)).
[0302] For prostate cancer indications, variants or mutations that may be biomarkers for immunotherapy treatment can include homologous recombination repair (HRR) gene (BRCA1, BRCA2, ATM, BARD1, BRIP1, CDK12, CHEK1, CHEK2, FANCL, PALB2, RAD51B, RAD51C, RAD51D, and RAD54L) alterations (e.g., which may be treated with Lynparza® (olaparib)).
[0303] For solid tumor cancer indications, variants or mutations that may be biomarkers for immunotherapy treatment may include a tumor mutational burden (TMB) greater than or equal to 10 mutations per megabase (e.g., as may be treated with Keytruda® (pembrolizumab)).
[0304] XI. Computing Machine Architecture 15 is a block diagram illustrating components of an exemplary computing machine capable of reading instructions from a computer-readable medium and executing them within a processor (or controller). The computers described herein may include a single computing machine as illustrated in FIG. 15, a virtual machine, a distributed computing system including multiple nodes of the computing machine as illustrated in FIG. 15, or any other suitable arrangement of computing devices.
[0305] 15 illustrates a diagrammatic representation of a computing machine in the exemplary form of a computer system 1500 within which instructions 1524 (e.g., software, program code, or machine code) may be stored on a computer-readable medium to cause the machine to perform any one or more of the processes discussed herein. In some embodiments, the computing machine may operate as a standalone device or may be connected (e.g., networked) to other machines. In a networked deployment, the machine may operate as a server machine or a client machine in a server-client network environment, or as a peer machine in a peer-to-peer (or distributed) network environment.
[0306] The structure of the computing machine illustrated in Figure 15 may correspond to any software, hardware, or combined components (e.g., processing units as illustrated in Figure 2 or described herein), including, but not limited to, any engine, module, computing server, machine used to perform one or more processes described herein. While Figure 15 illustrates various hardware and software elements, each of the components described herein may include additional or fewer elements.
[0307] By way of example, the computing machine may be a personal computer (PC), a tablet PC, a set-top box (STB), a personal digital assistant (PDA), a cellular phone, a smartphone, a web appliance, a network router, an Internet of Things (IoT) device, a switch or bridge, or any machine capable of executing instructions 1524 that specify actions to be taken by the machine. Additionally, although only a single machine is illustrated, the terms "machine" and "computer" may also be used to include any collection of machines that individually or jointly execute instructions 1524 to perform any one or more of the methodologies discussed herein.
[0308] An exemplary computer system 1500 includes one or more processors 1502, such as a CPU (Central Processing Unit), GPU (Graphics Processing Unit), TPU (Tensor Processing Unit), DSP (Digital Signal Processor), System on a Chip (SOC), controller, state machine, Application Specific Integrated Circuit (ASIC), Field Programmable Gate Array (FPGA), or any combination thereof. Portions of the computing system 1500 may also include memory 1504 that stores computer code including instructions 1524 that, when executed directly or indirectly by the processor 1502, may cause the processor 1502 to perform certain actions. Instructions may be any instructions, commands, or orders that may be stored in different forms, such as machine-readable instructions, programming instructions including source code, and other communication signals and orders. Instructions may be used in a general sense and is not limited to machine-readable code.
[0309] One or more methods described herein improve the operating speed of the processor 1502 and reduce the space required in the memory 1504. For example, the machine learning methods described herein reduce the computational complexity of the processor 1502 by applying one or more novel techniques that simplify the steps in training, reaching convergence, and generating results for the processor 1502. The algorithms described herein also reduce the size of models and datasets to reduce the storage space requirements of the memory 1504.
[0310] Performance of some of the operations may be distributed across processors that reside within a single machine or are spread across several machines. In some exemplary embodiments, one or more processors or processor-implemented modules may be located within a single geographic location (e.g., within a home environment, an office environment, or a server farm). In other exemplary embodiments, one or more processors or processor-implemented modules may be distributed across several geographic locations. While the specification or claims may refer to several processes to be performed by a processor, this should be interpreted to include the cooperative operation of multiple distributed processors.
[0311] Computer system 1500 may include a main memory 1504 and a static memory 1506, which are configured to communicate with each other via a bus 1508. Computer system 1500 may further include a graphics display unit 1510 (e.g., a plasma display panel (PDP), a liquid crystal display (LCD), a projector, or a cathode ray tube (CRT)). Graphics display unit 1510, controlled by processor 1502, displays a graphical user interface (GUI) to display one or more results and data generated by the processes described herein. The computer system 1500 may also include an alphanumeric input device 1512 (e.g., a keyboard), a cursor control device 1514 (e.g., a mouse, trackball, joystick, motion sensor, or other pointing instrument), a storage unit 1516 (a hard drive, solid-state drive, hybrid drive, memory disk, etc.), a signal generating device 1518 (e.g., a speaker), and a network interface device 1520, which are also configured to communicate via the bus 1508.
[0312] Storage unit 1516 includes computer-readable medium 1522 on which are stored instructions 1524 embodying any one or more of the methodologies or functions described herein. The instructions 1524 reside, completely or at least partially, within main memory 1504 or within processor 1502 (e.g., within the processor's cache memory) during their execution by computer system 1500, main memory 1504, and processor 1502, which also constitute computer-readable media. The instructions 1524 may be transmitted or received over network 1526 via network interface device 1520. While computer-readable medium 1522 is illustrated as being a single medium in the exemplary embodiment, the term “computer-readable medium” should be used to include a single medium or multiple media (e.g., a centralized or distributed database, or associated caches and servers) capable of storing instructions (e.g., instructions 1524). Computer-readable media may include any medium capable of storing instructions (e.g., instructions 1524) for execution by a processor (e.g., processor 1502) and causing the processor to perform any one or more of the methodologies disclosed herein. Computer-readable media may include, but are not limited to, data repositories in the form of solid-state memory, optical media, and magnetic media. Computer-readable media does not include transitory media such as propagated signals or carrier waves.
[0313] XII. Additional Considerations Beneficially, various embodiments described herein improve the accuracy and efficiency of existing techniques in the field of sequencing, such as PCR and massively parallel DNA sequencing (e.g., NGS). The embodiments provide a solution to the challenging problem of identifying errors introduced by the sequencing and amplification processes. Massively parallel DNA sequencing may begin with one or more DNA samples, which are randomly sheared and amplified, typically using PCR. The parallel nature of massively parallel DNA sequencing results in copies of the nucleotide sequence of each allele. The extent of copying and sequencing at each allele site may vary. For example, some sequences are overlapping and / or double-stranded, while others are not. Both the PCR amplification process and the sequencing process have a significant error rate. Sequence errors can serve to obscure the nucleotide sequence of the true allele. The embodiments may be used to determine one or more alleles analyzed by a massively parallel DNA sequencing instrument. By considering the read-layer intrinsic noise model, massively parallel DNA sequencing workflows exhibit sufficient fidelity to generate proper sequencing by more accurately distinguishing true alleles from erroneous sequences.
[0314] Conventionally, to reduce the error rate in determining the correct sequence, it is necessary to increase the sequencing depth of the sample. This means that fewer samples can be analyzed in a sequencing batch because more resources are dedicated to the sample. Embodiments improve sequencing accuracy without increasing the sequencing depth of specific allele sites, thereby enabling more allele sites or patient samples to be sequenced simultaneously in the case of massively parallel DNA sequencing. The described embodiments can reduce the required sequencing depth while increasing the accuracy of massively parallel DNA sequencing used to read the nucleotide sequences generated within the amplification.
[0315] The foregoing description of embodiments of the invention has been presented for purposes of illustration. It is not intended to be exhaustive or to limit the invention to the precise form disclosed. Those skilled in the art will appreciate that many modifications and variations are possible in light of the above disclosure.
[0316] Some portions of this specification describe embodiments of the invention in terms of algorithms and symbolic representations of operations on information. These algorithmic descriptions and representations are commonly used by those skilled in the data processing arts to effectively convey the substance of their work to others skilled in the art. While these operations are described functionally, computationally, or logically, they are understood to be implemented by computer programs or equivalent electrical circuits, microcode, or the like. Moreover, it has proven convenient at times to refer to these arrangements of operations as modules, without loss of generality. The described operations and their associated modules may be embodied in software, firmware, hardware, or any combination thereof.
[0317] Any of the steps, operations, or processes described herein may be performed or implemented with one or more hardware or software modules, alone or in combination with other devices. In some embodiments, software modules are implemented with a computer program product that includes a non-transitory computer-readable medium containing computer program code, the computer program code being executable by a computer processor to perform any or all of the described steps, operations, or processes.
[0318] Embodiments of the present invention also relate to products produced by the computing processes described herein. Such products may include information resulting from the computing processes, which information is stored on a non-transitory, tangible computer-readable storage medium, and may include any embodiment of a computer program product or other data combination described herein.
[0319] One or more processes described herein may be described with one or more steps, but the use of the term "step" does not imply a particular order. For example, the present disclosure may describe a process that includes multiple steps in sequence, but the steps within the process do not require that they be performed in the particular order claimed or described in the present disclosure. Some steps may be performed before other steps, even though the other steps are claimed or described first in the present disclosure.
[0320] Finally, the language used herein has been selected primarily for ease of reading and educational purposes, and may not be selected to elaborate or limit the subject matter of the present invention. Accordingly, the scope of the present invention is intended not to be limited by this detailed description, but rather by any claims issued in an application based thereon. Accordingly, the disclosure of embodiments of the present invention is intended to be illustrative, but not limiting, of the scope of the invention, which is set forth in the following claims.
Claims
1. 1. A computer-implemented method for processing a DNA sequencing dataset of a sample, comprising: accessing a DNA sequencing dataset generated by DNA sequencing, the DNA sequencing dataset including a plurality of processed sequence reads including variant locations; stratifying the plurality of processed sequencing reads into a plurality of read strata; For each read stratum, determining a stratified sequencing depth at the variant location; determining, for each lead stratum, one or more noise parameters conditional on the stratified sequencing depth of the lead stratum, the one or more noise parameters corresponding to a noise model specific to the lead stratum, training the noise model comprising: stratifying a training DNA dataset of a plurality of reference healthy individuals; selecting stratified sequencing reads for the read stratum as a stratified training set; setting initial values of the one or more noise parameters that model a noise distribution representative of the noise model; iteratively adjusting values of the one or more noise parameters based on the noise distribution of the stratified training set from the plurality of reference healthy individuals; determining, for each lead stratum, generating an output of the noise model specific to the lead stratum based on the one or more noise parameters conditional on the stratified sequencing depth of the lead stratum; combining the generated noise model outputs to produce a combined result representing the likelihood that the sample is associated with a total variant count; Including, The computer-implemented method, wherein the plurality of lead layers includes one or more of: (1) a double-stranded, sewn lead layer; (2) a double-stranded, unsewn lead layer; (3) a single-stranded, sewn lead layer; and (4) a single-stranded, unsewn lead layer.
2. 2. The computer-implemented method of claim 1, wherein the mutation at the variant location is one of a single base variant, an insertion, and a deletion.
3. determining a quality score of the combined result, the quality score being a Phred scale score. The computer-implemented method of claim 1 , further comprising:
4. In response to the quality score being greater than a predetermined threshold, indicating that the sample has a high likelihood of carrying a mutation at the variant location. The computer-implemented method of claim 3 , further comprising:
5. Determining, for a lead stratum, the one or more noise parameters conditional on the stratified sequencing depth of the lead stratum, comprises: accessing a parameter distribution specific to the lead stratum, the parameter distribution describing a distribution of a set of DNA sequencing samples associated with the lead stratum, and the noise parameter being determined from the parameter distribution; 2. The computer-implemented method of claim 1, comprising:
6. 6. The computer-implemented method of claim 5, wherein for each lead stratum, the set of DNA sequencing samples associated with the lead stratum includes sequence reads stratified into the lead stratum and correspond to one or more healthy individuals.
7. The computer-implemented method of claim 5 , wherein for each lead stratum, the noise model specific to that lead stratum is a Bayesian hierarchical model and the parameter distribution is based on a gamma distribution.
8. 2. The computer-implemented method of claim 1, wherein a first noise parameter corresponding to a noise model specific to a first lead layer has a different value than a corresponding second noise parameter corresponding to a noise model specific to a second lead layer.
9. 2. The computer-implemented method of claim 1, wherein for each lead stratum, the determined one or more noise parameters comprise a mean of the noise distribution conditional on the stratified sequencing depth of the lead stratum.
10. 10. The computer-implemented method of claim 9, wherein each noise distribution is a negative binomial distribution conditional on the stratified sequencing depth of each read stratum.
11. The computer-implemented method of claim 10 , wherein for each lead stratum, the determined one or more noise parameters further include a variability parameter.
12. 2. The computer-implemented method of claim 1, wherein the generated output of each noise model is the one or more noise parameters conditional on the stratified sequencing depth determined for the lead stratum.
13. 2. The computer-implemented method of claim 1, wherein the generated output of each noise model comprises a likelihood that a stratified variant count for the lead stratum exceeds a threshold.
14. 2. The computer-implemented method of claim 1, wherein combining the generated noise model outputs comprises combining an average variant count and a variance from each noise model output to produce an overall average variant count and an overall variability parameter that represents an overall noise distribution for the combined result.
15. The global noise distribution is modeled based on a negative binomial distribution, and determining the global mean variant count and the global variability parameters comprises: determining the average variant count for each lead stratum based on the stratified sequencing depth of the lead stratum; determining the variance for each lead layer; summing the average variant counts for each lead stratum to determine the overall average variant count; combining the variances for each lead layer to determine an overall variance; determining the overall variability parameter based on the overall average variant count and the overall variance; 15. The computer-implemented method of claim 14, comprising:
16. Combining the generated noise model outputs to produce the combined result comprises: determining an observed stratified variant count for each lead stratum; determining, within each lead stratum, a possible event that is more likely than the observed stratified variant count for each lead stratum; identifying the possible event combinations associated with a likelihood of occurrence higher than the observed stratified variant count for each lead stratum; summing the probabilities of the identified combinations to determine a statistical complement; determining a likelihood value by subtracting said statistical complement from 1.0; 2. The computer-implemented method of claim 1, comprising:
17. 20. The computer-implemented method of claim 16, wherein a first identified combination comprising one double-stranded read is equal to a second identified combination comprising two single-stranded reads.
18. 20. The computer-implemented method of claim 16, wherein the determined likelihood values are equal to or greater than the likelihood of occurrence of the observed stratified variant counts for each lead stratum.
19. 20. The computer-implemented method of claim 16, further comprising training a machine learning model to determine the likelihood value.
20. receiving a bodily fluid sample from an individual; performing said DNA sequencing on the cfDNA of said bodily fluid sample; generating raw sequence reads based on the results of the DNA sequencing; collapsing and stitching the raw sequence reads to generate the plurality of processed sequence reads; The computer-implemented method of claim 1 , further comprising:
21. 21. The computer-implemented method of claim 20, wherein the bodily fluid sample is one of a sample of blood, whole blood, plasma, serum, urine, cerebrospinal fluid, feces, saliva, tears, tissue biopsy, pleural fluid, pericardial fluid, or peritoneal fluid of the individual.
22. 21. The computer-implemented method of claim 20, wherein the plurality of processed sequencing reads are sequenced from a tumor biopsy.
23. 21. The computer-implemented method of claim 20, wherein the plurality of processed sequencing reads are sequenced from an isolate of cells from blood, the isolate of cells including at least buffy coat white blood cells or CD4+ cells.
24. The computer-implemented method of claim 1 , wherein the DNA sequencing comprises a massively parallel DNA sequencing operation.
25. 2. The computer-implemented method of claim 1, wherein the DNA sequencing dataset is a cfDNA sequencing dataset of a bodily fluid sample of an individual.
26. providing a diagnosis that the subject has a variant based on the combined result. The computer-implemented method of claim 1 , further comprising:
27. The variants include ACVR1B, AKT3, AMER1, APC, ARID1A, ARID1B, ARID2, ASXL1, ASXL2, ATM, ATR, and BAP1. BCL2, BCL6, BCORL1, BCR, BLM, BRAF, BRCA1, BTG1, CASP8, CBL, CCND3, CCNE1, CD74, CDC73, CDK12, CDKN2A, CHD2, CJD2, C REBBP, CSF1R, CTCF, CTNNB1, DICER1, DNAJB1, DNMT1, DNMT3A, DNMT3B, DOT1L, EED, EGFR, EIF1AX, EP300, EPHA3, EPHA5, E PHB1, ERBB2, ERBB4, ERCC2, ERCC3, ERCC4, ESR1, FAM46C, FANCA, FANCC, FANCD2, FANCE, FAT1, FBXW7, FGFR3, FLCN, FLT1, FOXO1, FUBP1, FYN, GATA3, GPR124, GRIN2A, GRM3, H3F3A, HIST1H1C, IDH1, IDH2, IKZF1, IL7R, INPP4B, IRF4, IRS1, IRS2, JAK2, KAT6A, KDM6A, KEAP1, KIF5B, KIT, KLF4, KLH6, KMT2C, KRAS, LMAP1, LRP1B, LZTR1, MAP3K1, MCL1, MGA, MSH2, MSH6, M ST1R, MTOR, MYD88, NPM1, NRAS, NTRK1, NTRK2, NUP93, NUTM1, PAX3, PAX8, PBRM1, PGR, PHOX2B, PIK3CA, POLE, PTCH1, PTEN , PTPN11, PTPRT, RAD21, RAF1, RANBP2, RB1, REL, RFWD2, RHOA, RPTOR, RUNX1, RUNX1T1, SDHA, SHQ1, SLIT2, SMAD4, SMARCA4, SMARCD1, SNCAIP, SOCS1, SPEN, SPTA1, SUZ12, TET1, TET2, TGFBR, and TNFRSF14.
28. providing instructions for administering a treatment to said subjects identified as having said variant.
27. The computer-implemented method of claim 26, further comprising:
29. 30. The computer-implemented method of claim 28, wherein the treatment comprises administering an agent selected from the group consisting of Rituxan, Herceptin, Erbitux, Vectibix, Arzerra, Benlysta, Yervoy, Perjeta, Tremelimumab, Opdivo, Dacetuzumab, Urelumab, Tecentriq, Lambrolizumab, Blinatumomab, CT-011, Keytruda, BMS-936559, MED14736, MSB0010718C, Imfinzi, Bavencio, and Margetuximab.
30. 2. The computer-implemented method of claim 1, wherein the likelihood represents that a total variant count for subsequently observed data is greater than or equal to a total variant count observed within the plurality of processed sequencing reads due to noise.
31. When executed by one or more processors, the one or more processors are accessing a DNA sequencing dataset generated by DNA sequencing, the DNA sequencing dataset including a plurality of processed sequence reads including variant locations; stratifying the plurality of processed sequencing reads into a plurality of read strata; determining a stratified sequencing depth at said variant locations for each read stratum; determining, for each lead stratum, one or more noise parameters conditional on the stratified sequencing depth of the lead stratum, the one or more noise parameters corresponding to a noise model specific to the lead stratum, training the noise model comprising: stratifying a training DNA dataset of a plurality of reference healthy individuals; selecting stratified sequencing reads for the read stratum as a stratified training set; setting initial values of the one or more noise parameters that model a noise distribution representative of the noise model; iteratively adjusting values of the one or more noise parameters based on the noise distribution of the stratified training set from the plurality of reference healthy individuals; determining a value for the number of times ... generating, for each lead stratum, an output of the noise model specific to the lead stratum based on the one or more noise parameters conditional on the stratified sequencing depth of the lead stratum; then combining the generated noise model outputs to produce a combined result that represents the likelihood that a total variant count for the observed data is greater than or equal to a total variant count observed in the plurality of processed sequencing reads due to noise; comprising instructions to perform steps including: The plurality of lead layers includes one or more of: (1) a double-stranded, sewn lead layer; (2) a double-stranded, unsewn lead layer; (3) a single-stranded, sewn lead layer; and (4) a single-stranded, unsewn lead layer.
32. 32. The non-transitory computer-readable medium of claim 31 , wherein combining the generated noise model outputs comprises combining an average variant count and a variance from each noise model output to produce an overall average variant count and the overall variability parameter representing an overall noise distribution for the combined result.
33. The global noise distribution is modeled based on a negative binomial distribution, and determining the global mean variant count and the global variability parameters comprises: determining the average variant count for each lead stratum based on the stratified sequencing depth of the lead stratum; determining the variance for each lead layer; summing the average variant counts for each lead stratum to determine the overall average variant count; combining the variances for each lead layer to determine an overall variance; determining the overall variability parameter based on the overall average variant count and the overall variance; 33. The non-transitory computer readable medium of claim 32, comprising:
34. Combining the generated noise model outputs to produce the combined result comprises: determining an observed stratified variant count for each lead stratum; determining, within each lead stratum, a possible event that is more likely than the observed stratified variant count for each lead stratum; identifying the possible event combinations associated with a likelihood of occurrence higher than the observed stratified variant count for each lead stratum; summing the probabilities of the identified combinations to determine a statistical complement; determining a likelihood value by subtracting said statistical complement from 1.0; 32. The non-transitory computer readable medium of claim 31, comprising:
35. The steps include: providing a diagnosis that the subject has a variant based on the combined result.
32. The non-transitory computer readable medium of claim 31, further comprising:
36. The variants include ACVR1B, AKT3, AMER1, APC, ARID1A, ARID1B, ARID2, ASXL1, ASXL2, ATM, ATR, and BAP1. BCL2, BCL6, BCORL1, BCR, BLM, BRAF, BRCA1, BTG1, CASP8, CBL, CCND3, CCNE1, CD74, CDC73, CDK12, CDKN2A, CHD2, CJD2, CR EBBP, CSF1R, CTCF, CTNNB1, DICER1, DNAJB1, DNMT1, DNMT3A, DNMT3B, DOT1L, EED, EGFR, EIF1AX, EP300, EPHA3, EPHA5, EPH B1, ERBB2, ERBB4, ERCC2, ERCC3, ERCC4, ESR1, FAM46C, FANCA, FANCC, FANCD2, FANCE, FAT1, FBXW7, FGFR3, FLCN, FLT1, FO XO1, FUBP1, FYN, GATA3, GPR124, GRIN2A, GRM3, H3F3A, HIST1H1C, IDH1, IDH2, IKZF1, IL7R, INPP4B, IRF4, IRS1, IRS2, JAK 2, KAT6A, KDM6A, KEAP1, KIF5B, KIT, KLF4, KLH6, KMT2C, KRAS, LMAP1, LRP1B, LZTR1, MAP3K1, MCL1, MGA, MSH2, MSH6, MST1 R, MTOR, MYD88, NPM1, NRAS, NTRK1, NTRK2, NUP93, NUTM1, PAX3, PAX8, PBRM1, PGR, PHOX2B, PIK3CA, POLE, PTCH1, PTEN, PTP 36. The non-transitory computer readable medium of claim 35, wherein the nucleic acid sequence is selected from the group consisting of N11, PTPRT, RAD21, RAF1, RANBP2, RB1, REL, RFWD2, RHOA, RPTOR, RUNX1, RUNX1T1, SDHA, SHQ1, SLIT2, SMAD4, SMARCA4, SMARCD1, SNCAIP, SOCS1, SPEN, SPTA1, SUZ12, TET1, TET2, TGFBR, and TNFRSF14.
37. The steps include: providing instructions for administering a treatment to said subjects identified as having said variant.
36. The non-transitory computer readable medium of claim 35, further comprising:
38. 38. The non-transitory computer readable medium of claim 37, wherein the treatment comprises administering an agent selected from the group consisting of Rituxan, Herceptin, Erbitux, Vectibix, Arzerra, Benlysta, Yervoy, Perjeta, Tremelimumab, Opdivo, Dacetuzumab, Urelumab, Tecentriq, Lambrolizumab, Blinatumomab, CT-011, Keytruda, BMS-936559, MED14736, MSB0010718C, Imfinzi, Bavencio, and Margetuximab.
39. 32. The non-transitory computer-readable medium of claim 31 , wherein the likelihood represents that a total variant count for subsequently observed data is greater than or equal to a total variant count observed in the plurality of processed sequencing reads due to noise.
40. 1. A system comprising a computer processor and a memory storing computer program instructions, the computer program instructions, when executed by the computer processor, causing the computer processor to: accessing a DNA sequencing dataset generated by DNA sequencing, the DNA sequencing dataset including a plurality of processed sequence reads including variant locations; stratifying the plurality of processed sequencing reads into a plurality of read strata; determining a stratified sequencing depth at said variant locations for each read stratum; determining, for each lead stratum, one or more noise parameters conditional on the stratified sequencing depth of the lead stratum, the one or more noise parameters corresponding to a noise model specific to the lead stratum, training the noise model comprising: stratifying a training DNA dataset of a plurality of reference healthy individuals; selecting stratified sequencing reads for the read stratum as a stratified training set; Steps and iteratively adjusting values of the one or more noise parameters based on the noise distribution of the stratified training set from the plurality of reference healthy individuals; determining a value for the number of times ... generating, for each lead stratum, an output of the noise model specific to the lead stratum based on the one or more noise parameters conditional on the stratified sequencing depth of the lead stratum; then combining the generated noise model outputs to produce a combined result that represents the likelihood that a total variant count for the observed data is greater than or equal to a total variant count observed in the plurality of processed sequencing reads due to noise; performing a step including The system, wherein the multiple lead layers include one or more of: (1) a double-stranded, sewn lead layer; (2) a double-stranded, non-sewn lead layer; (3) a single-stranded, sewn lead layer; and (4) a single-stranded, non-sewn lead layer.
41. 41. The system of claim 40, wherein combining the generated noise model outputs comprises combining an average variant count and a variance from each noise model output to produce an overall average variant count and an overall variability parameter that represents an overall noise distribution for the combined result.
42. The global noise distribution is modeled based on a negative binomial distribution, and determining the global mean variant count and the global variability parameters comprises: determining the average variant count for each lead stratum based on the stratified sequencing depth of the lead stratum; determining the variance for each lead layer; summing the average variant counts for each lead stratum to determine the overall average variant count; combining the variances for each lead layer to determine an overall variance; determining the overall variability parameter based on the overall average variant count and the overall variance; 42. The system of claim 41, comprising:
43. Combining the generated noise model outputs to produce the combined result comprises: determining an observed stratified variant count for each lead stratum; determining, within each lead stratum, a possible event that is more likely than the observed stratified variant count for each lead stratum; identifying the possible event combinations associated with a likelihood of occurrence higher than the observed stratified variant count for each lead stratum; summing the probabilities of the identified combinations to determine a statistical complement; determining a likelihood value by subtracting said statistical complement from 1.0; 41. The system of claim 40, comprising:
44. The steps include: providing a diagnosis that the subject has a variant based on the combined result.
41. The system of claim 40, further comprising:
45. The variants include ACVR1B, AKT3, AMER1, APC, ARID1A, ARID1B, ARID2, ASXL1, ASXL2, ATM, ATR, and BAP1. BCL2, BCL6, BCORL1, BCR, BLM, BRAF, BRCA1, BTG1, CASP8, CBL, CCND3, CCNE1, CD74, CDC73, CDK12, CDKN2A, CHD2, CJD2, C REBBP, CSF1R, CTCF, CTNNB1, DICER1, DNAJB1, DNMT1, DNMT3A, DNMT3B, DOT1L, EED, EGFR, EIF1AX, EP300, EPHA3, EPHA5, EPHB1, ERBB2, ERBB4, ERCC2, ERCC3, ERCC4, ESR1, FAM46C, FANCA, FANCC, FANCD2, FANCE, FAT1, FBXW7, FGFR3, FLCN, FLT 1, FOXO1, FUBP1, FYN, GATA3, GPR124, GRIN2A, GRM3, H3F3A, HIST1H1C, IDH1, IDH2, IKZF1, IL7R, INPP4B, IRF4, IRS1, IR S2, JAK2, KAT6A, KDM6A, KEAP1, KIF5B, KIT, KLF4, KLH6, KMT2C, KRAS, LMAP1, LRP1B, LZTR1, MAP3K1, MCL1, MGA, MSH2, MS H6, MST1R, MTOR, MYD88, NPM1, NRAS, NTRK1, NTRK2, NUP93, NUTM1, PAX3, PAX8, PBRM1, PGR, PHOX2B, PIK3CA, POLE, PTCH1 , PTEN, PTPN11, PTPRT, RAD21, RAF1, RANBP2, RB1, REL, RFWD2, RHOA, RPTOR, RUNX1, RUNX1T1, SDHA, SHQ1, SLIT2, SMAD4, SMARCA4, SMARCD1, SNCAIP, SOCS1, SPEN, SPTA1, SUZ12, TET1, TET2, TGFBR, and TNFRSF14.
46. The steps include: providing instructions for administering a treatment to said subjects identified as having said variant.
45. The system of claim 44, further comprising:
47. 47. The system of claim 46, wherein the treatment comprises administering an agent selected from the group consisting of Rituxan, Herceptin, Erbitux, Vectibix, Arzerra, Benlysta, Yervoy, Perjeta, Tremelimumab, Opdivo, Dacetuzumab, Urelumab, Tecentriq, Lambrolizumab, Blinatumomab, CT-011, Keytruda, BMS-936559, MED14736, MSB0010718C, Imfinzi, Bavencio, and Margetuximab.
48. 41. The system of claim 40, wherein the likelihood represents that a total variant count for subsequently observed data is greater than or equal to a total variant count observed in the plurality of processed sequencing reads due to noise.
Citation Information
Patent Citations
How to adjust molecular label counts
JP2019522268A
Improved variant detection
JP2021516962A
Identifying false positive variants using a significance model
US20190073445A1
Site-specific noise model for targeted sequencing
US20190108311A1
Methods and materials for assessing and treating cancer
US20190256924A1
