CN114746947A - 用于分析dna数据的读数层特定噪声模型 - Google Patents

用于分析dna数据的读数层特定噪声模型 Download PDF

Info

Publication number
CN114746947A
CN114746947A CN202080077722.0A CN202080077722A CN114746947A CN 114746947 A CN114746947 A CN 114746947A CN 202080077722 A CN202080077722 A CN 202080077722A CN 114746947 A CN114746947 A CN 114746947A
Authority
CN
China
Prior art keywords
variant
noise
read
layer
hierarchical
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202080077722.0A
Other languages
English (en)
Inventor
E·哈贝尔
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Greer Co ltd
Original Assignee
Greer Co ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Greer Co ltd filed Critical Greer Co ltd
Publication of CN114746947A publication Critical patent/CN114746947A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16BBIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B20/00ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • CCHEMISTRY; METALLURGY
    • C12BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12QMEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6869Methods for sequencing
    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16BBIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/20Sequence assembly
    • GPHYSICS
    • G16INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16BBIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B40/00ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
    • G16B40/30Unsupervised data analysis

Landscapes

  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Biophysics (AREA)
  • Biotechnology (AREA)
  • General Health & Medical Sciences (AREA)
  • Organic Chemistry (AREA)
  • Medical Informatics (AREA)
  • Analytical Chemistry (AREA)
  • Theoretical Computer Science (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Zoology (AREA)
  • Wood Science & Technology (AREA)
  • Molecular Biology (AREA)
  • Genetics & Genomics (AREA)
  • Data Mining & Analysis (AREA)
  • Immunology (AREA)
  • General Engineering & Computer Science (AREA)
  • Biochemistry (AREA)
  • Microbiology (AREA)
  • Databases & Information Systems (AREA)
  • Bioethics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Artificial Intelligence (AREA)
  • Epidemiology (AREA)
  • Evolutionary Computation (AREA)
  • Public Health (AREA)
  • Software Systems (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
  • Apparatus Associated With Microorganisms And Enzymes (AREA)

Abstract

用于处理核酸数据集的噪声模型可以将已处理的序列读数分层到不同的读数层。每个读数层可以是基于潜在的变体位置是否在所述序列读数的重叠区域和/或互补区域来定义的。处理系统可以针对每个读数层确定在所述变体位置处的分层测序深度。所述处理系统可以针对每个读数层确定以所述读数层的所述分层测序深度为条件的一个或多个噪声参数。所述噪声参数可以与噪声分布相关联。所述处理系统可以基于以所述分层测序深度为条件的所述噪声参数生成每个噪声模型的输出。所述处理系统可以组合每个分层噪声模型的输出以生成组合结果,所述组合结果可以表示事件与观察到的所述数据一样或者比观察到的所述数据更极端的可能性。

Description

用于分析DNA数据的读数层特定噪声模型
技术领域
本公开大体上涉及用于确定核酸测序数据集的质量分数的噪声模型。
背景技术
计算技术可以被用于DNA测序数据,以标识可能对应于各种类型的癌症或其他疾病的DNA中的突变或变体。因此,癌症诊断或预测可以通过分析生物样本来执行,诸如从个体、动物、植物等中抽取的组织活检或血液。从血液样本中检测出源于肿瘤细胞的DNA是困难的,这是因为循环肿瘤DNA(ctDNA)相对于从血液中提取的无细胞DNA(cfDNA)中的其他分子而言以低水平存在。现有方法无法从信号噪声中标识真阳性(例如,指示受试者中的癌症)降低了已知系统和未来系统区分真阳性与由噪声源引起的假阳性的能力,这可能会导致变体调用或其他类型的分析的不可靠结果。
发明内容
本文公开了用于训练和应用位点特定噪声模型的系统和方法,这些位点特定噪声模型被分类为多个读数层。噪声模型可以确定靶向测序中真阳性的可能性。真阳性可以包括单核苷酸变体、插入或碱基对的缺失。具体地,模型可以使用贝叶斯推理来确定核酸序列的每个位置的噪声的比率或水平,例如指示某些突变的预期可能性。每个模型可以特定于读数层。读数层可以是基于潜在的变体位置是否在已处理的测序读数的重叠区域和/或互补区域来确定的。特定于读数层的每个模型可以是考虑了协变量(例如三核苷酸上下文、可映射性或片段重复)和特定于该读数层的各种类型的参数(例如混合物组分或序列读数的深度)的层次模型。模型可以从健康受试者的序列读数训练得到,这些序列读数也通过读数层分层。不同噪声模型的输出可以被组合以生成总体质量分数。与不会通过读数层区分序列读数的单个模型相比,包含各种读数层模型的总体管线可以以较高灵敏度标识出真阳性,并且过滤掉假阳性。
通过示例,在各种实施例中,一种用于处理样本(例如个体)的DNA测序数据集的方法可以包括:访问由DNA测序生成的测序数据集,该DNA测序数据集包括多个已处理的序列读数,该多个已处理的序列读数包括变体位置。该方法还可以包括:将多个已处理的序列读数分层到多个读数层。该方法还可以包括:针对每个读数层,确定在变体位置处的分层测序深度。该方法还可以包括:针对每个读数层,确定以读数层的分层测序深度为条件的一个或多个噪声参数,该一个或多个噪声参数对应特定于读数层的噪声模型。该方法还可以包括:针对每个读数层,基于以读数层的分层测序深度为条件的一个或多个噪声参数,生成特定于读数层的噪声模型的输出。该方法还可以包括:组合所生成的噪声模型输出以产生组合结果。该组合结果可以表示后续观察的数据的总变体计数大于或等于多个已处理序列读数中观察到的总变体计数可归因于噪声的可能性。
在一个或多个实施例中,多个读数层包括以下一项或多项:(1)双链的缝合读数层,(2)双链的未缝合读数层,(3)单链的缝合读数层,以及(4)单链的未缝合读数层。
在一个或多个实施例中,在变体位置处的突变是以下之一:单核苷酸变体、插入和缺失。
在一个或多个实施例中,该方法还可以包括:确定组合结果的质量分数,该质量分数是Phred等级分数。
在一个或多个实施例中,该方法还可以包括:响应于质量分数高于预定阈值,指示样本可能在变体位置处有突变。
在一个或多个实施例中,针对每个读数层,确定以读数层的分层测序深度为条件的一个或多个噪声参数可以包括:访问特定于读数层的参数分布,该参数分布描述了与读数层相关联的一组DNA测序样本的分布。噪声参数是从参数分布中确定的。
在一个或多个实施例中,针对每个读数层,与读数层相关联的一组DNA测序样本包括分层到读数层的序列读数,并且对应于一个或多个健康个体。
在一个或多个实施例中,针对每个读数层,特定于读数层的噪声模型是贝叶斯层次模型,并且参数分布基于伽马分布。
在一个或多个实施例中,与特定于第一读数层的噪声模型相对应的第一噪声参数具有与特定于第二读数层的噪声模型相对应的对应第二噪声参数不同的值。
在一个或多个实施例中,针对每个读数层,所确定的一个或多个噪声参数包括以读数层的分层测序深度为条件的噪声分布的均值。
在一个或多个实施例中,每个噪声分布均是以每个读数层的分层测序深度为条件的负二项分布。
在一个或多个实施例中,针对每个读数层,所确定的一个或多个噪声参数还包括分散参数。
在一个或多个实施例,每个噪声模型的输出是针对读数层确定的以分层测序深度为条件的一个或多个噪声参数。
在一个或多个实施例,每个噪声模型的所生成的输出是针对读数层确定的以分层测序深度为条件的一个或多个噪声参数。
在一个或多个实施例中,每个噪声模型的所生成的输出包括读数层的分层变体计数超过阈值的可能性。
在一个或多个实施例中,组合所生成的噪声模型输出包括组合来自每个噪声模型输出的平均变体计数和方差以产生表示组合结果的总体噪声分布的总体平均变体计数和总体分散参数。
在一个或多个实施例中,总体噪声分布是基于负二项分布而被建模的。确定总体平均变体计数和总体分散参数可以包括:基于读数层的分层测序深度,确定针对每个读数层的平均变体计数。确定步骤还可以包括:确定每个读数层的方差。确定步骤还可以包括:对每个读数层的平均变体计数进行求和以确定总体平均变体计数。确定步骤还可以包括:组合每个读数层的方差以确定总体方差。确定步骤还可以包括:基于总体平均变体计数和总体方差,确定总体分散参数。
在一个或多个实施例中,组合每个噪声模型的输出以生成组合结果可以包括:确定每个读数层的观察到的分层变体计数。组合步骤还可以包括:在每个读数层中确定比每个读数层的观察到的分层变体计数更有可能的可能事件。组合步骤还可以包括:标识与比每个读数层的观察到的分层变体计数更高的发生可能性相关联的可能事件的组合。组合步骤还可以包括:对所标识的组合的概率进行求和以确定统计补数。组合步骤还可以包括:通过从1.0中减去统计补数来确定可能性值。
在一个或多个实施例中,包括一个双链读数的第一标识组合等同于包括两个单链读数的第二标识组合。
在一个或多个实施例中,所确定的可能性值等于或大于每个读数层的观察到的分层变体计数的发生可能性。
在一个或多个实施例中,该方法还可以包括:训练机器学习模型以确定可能性值。
在一个或多个实施例中,该方法还可以包括:接收个体的体液样本。该方法还可以包括:对体液样本的cfDNA执行DNA测序。该方法还可以包括:基于DNA测序的结果,生成原始序列读数。该方法还可以包括:折叠并且缝合原始序列读数以生成多个已处理的序列读数。
在一个或多个实施例中,体液样本是以下之一的样本:个体的血液、全血、血浆、血清、尿液、脑脊液、粪便、唾液、眼泪、组织活检、胸膜液、心包液或腹膜液。
在一个或多个实施例中,多个已处理序列读数是从肿瘤活检测序的。
在一个或多个实施例中,多个已处理的序列读数是从来自血液的细胞分离物测序的,该细胞分离物至少包括血沉棕黄层白血细胞或CD4+细胞。
在一个或多个实施例中,DNA测序是一种大规模并行DNA测序。
在各种实施例中,一种非瞬态计算机可读介质包括指令,在由一个或多个处理器执行时,该指令使一个或多个处理器执行上面描述的和本文公开的任何步骤。
进一步地,在各种实施例中,提供了一种系统,该系统具有计算机处理器和存储器,该存储器存储计算机程序指令,通过计算机处理器执行指令使处理器执行上面描述的和本文公开的任何步骤。
根据本发明的实施例具体地在针对方法和计算机程序产品的所附权利要求中公开,其中在一个权利要求类别(例如方法)中提及的任何特征也可以在另一权利要求类别(例如计算机程序产品、系统、存储介质)中要求保护。所附权利要求中的相关性或引用仅仅是出于形式上的原因而选择的。然而,由对任何先前权利要求(具体地多种相关性)的有意引用而产生的任何主题也可以要求保护,从而不管所附权利要求中选择的相关性,公开了权利要求的任何组合及其特征,并且可以要求保护。可以要求保护的主题不仅包括所附权利要求中陈述的特征的组合,还包括权利要求中的特征的任何其他组合,其中权利要求中提及的每个特征可以与权利要求中的任何其他特征或其他特征的组合组合。此外,本文描述或描绘的任何实施例和特征可以在单独的权利要求中和/或与本文描述或描绘的任何实施例或特征或者所附权利要求的任何特征的任何组合中要求保护。
附图说明
图1是根据本公开的各种实施例的制备用于测序的核酸样本的方法的流程图。
图2是根据本公开的各种实施例的用于处理序列读数的处理系统的框图。
图3是根据本公开的各种实施例的用于确定序列读数的变体的方法的流程图。
图4是根据本公开的各种实施例的应用贝叶斯层次模型的图。
图5A示出了根据本公开的各种实施例的用于确定真实单核苷酸变体的贝叶斯层次模型的参数与子模型之间的相关性。
图5B示出了根据本公开的各种实施例的用于确定真实插入或缺失的贝叶斯层次模型的参数与子模型之间的相关性。
图6A图示了根据本公开的各种实施例的与贝叶斯层次模型相关联的噪声比率的分布的图。
图6B图示了根据本公开的各种实施例的与贝叶斯层次模型相关联的变体给定参数的分布的图。
图7A是根据本公开的各种实施例的通过拟合贝叶斯层次模型来确定参数的图。
图7B是根据本公开的各种实施例的使用来自贝叶斯层次模型的参数来确定假阳性的可能性的图。
图8A包括图示了根据本公开的各种实施例的序列读数的不同读数层。
图8B示出了根据本公开的各种实施例的图示图8A的读数层的不同质量的实验结果。
图8C示出了根据本公开的各种实施例的通过基于核苷酸替换的类型将序列读数分层到读数层并且进一步分层到子读数层而获得的第一读数层的实验结果。
图8D示出了根据本公开的各种实施例的第二读数层的实验结果该第二读数层通过基于核苷酸替换的类型将序列读数分层到读数层并且进一步分层到子读数层而获得的。
图8E示出了根据本公开的各种实施例的第三读数层的实验结果,该第三读数层通过基于核苷酸替换的类型将序列读数分层到读数层并且进一步分层到子读数层而获得的。
图8F示出了根据本公开的各种实施例的第四读数层的实验结果,该第四读数层通过基于核苷酸替换的类型将序列读数分层到读数层并且进一步分层到子读数层而获得的。
图8G示出了根据本公开的各种实施例的第五读数层的实验结果,该第五读数层通过基于核苷酸替换的类型将序列读数分层到读数层并且进一步分层到子读数层而获得的。
图8H示出了根据本公开的各种实施例的跨基于交替等位基因的类型的四个读数层的平均误差率。
图9是描绘根据本公开的各种实施例的使用分层噪声模型分析DNA测序样本的过程的流程图。
图10是描绘根据本公开的各种实施例的使用矩匹配来组合不同读数层的分层噪声模型的输出的过程的流程图。
图11A是描绘根据本公开的各种实施例的使用积分来组合不同读数层的分层噪声模型的输出的过程的流程图。
图11B图示了根据本公开的各种实施例的多维空间中的更极端事件的计数。
图12A图示了根据本公开的各种实施例的观察到的质量分数与默认质量分数的示例图表。
图12B图示了根据本公开的各种实施例的观察到的质量分数与默认质量分数的另一示例图表。
图13A图示了根据本公开的各种实施例的使用读数层的质量分数的实验结果。
图13B图示了根据本公开的各种实施例的使用不区分读数层的噪声模型的质量分数的实验结果。
图14是根据本公开的各种实施例的描绘标识个体的潜在突变位置的示例过程的流程图。
图15是根据本公开的各种实施例的示例计算设备的框图。
附图仅出于图示的目的描绘本发明的实施例。本领域技术人员将从以下讨论容易地认识到,在不脱离本文描述的本发明的原理的情况下,本文图示的结构和方法的替代实施例可以被采用。
具体实施方式
I.定义
术语“个体”是指人类个体。术语“健康个体”是指假定为没有癌症或疾病的个体。术语“受试者”是指正在接受癌症或疾病测试的个体。
术语“序列读数”是指来自从个体获得的样本中的核苷酸序列读数。序列读数可以通过本领域已知的各种方法获得。
术语“读数片段”或“读数”是指任何核苷酸序列,包括从个体获得的序列读数和/或源自从个体获得的样本的初始序列读数的核苷酸序列。例如,读数片段可以指比对的序列读数、折叠的序列读数或缝合的读数。此外,读数片段可以指单个核苷酸碱基,诸如单个核苷酸变体。
术语“单核苷酸变体”或“SNV”是指在核苷酸序列(例如来自个体的序列读数)的位置(例如位点)处将一个核苷酸替换为不同核苷酸。从第一核碱基X到第二核碱基Y的替换可以被表示为“X>Y”。例如,胞嘧啶到胸腺嘧啶SNV可以被表示为“C>T”。
术语“插入/缺失”是指在序列读数中具有长度和位置(也可以被称为锚定位置)的一个或多个碱基对的任何插入或缺失。插入对应于正长度,而缺失对应于负长度。
术语“变体”是指一个或多个SNV或插入/缺失。变体位置是指DNA测序中可能包含SNV或插入/缺失的感兴趣位置。
术语“真阳性”是指指示真实生物学的突变,例如个体中存在潜在的癌症、疾病或种系突变。真阳性不是由健康个体中自然发生的突变(例如频发突变)或其他伪影来源(诸如核酸样本测定制备期间的过程错误)引起的。
术语“假阳性”是指错误地确定为真阳性的突变。通常,在处理与更大的平均噪声率或更大的噪声率不确定性相关联的序列读数时,更可能发生假阳性。
术语“无细胞DNA”或“cfDNA”是指在个体体内(例如血液)循环并且源自一个或多个健康细胞和/或一个或多个癌细胞的核酸碎片。
术语“循环肿瘤DNA”或“ctDNA”是指源自肿瘤细胞或其他类型的癌细胞的核酸碎片,其可以由于诸如凋亡或垂死细胞坏死等生物学过程而释放到个体的血液中或由活的肿瘤细胞主动释放。
术语“交替等位基因”或“ALT”是指相对于参考等位基因(例如对应于已知基因)具有一个或多个突变的等位基因。
术语“测序深度”或“深度”是指来自从个体获得的样本的读数片段的总数。
术语“交替深度”或“AD”是指样本中支持ALT(例如包括ALT的突变)的多个读数片段。
术语“交替频率”或“AF”是指给定ALT的频率。AF可以通过将样本的对应AD除以给定ALT的样本深度来确定。
II.示例测定协议
图1是根据各种实施例的用于制备用于测序的核酸样本的方法100的流程图。方法100包括但不限于以下步骤。例如,方法100的任何步骤可以包括用于质量控制的定量子步骤或本领域技术人员已知的其他实验室测定程序。方法100可以对应于一种类型的大规模并行DNA测序,例如下一代测序(NGS)。
在步骤110中,从受试者提取核酸样本(DNA或RNA)。受试者可以是个体。样本可以是人类基因组的任何子集或整个基因组。样本可以从已知患有或怀疑患有癌症的受试者中提取。样本可以包括血液、血浆、血清、尿液、粪便、唾液、其他类型的体液或其任何组合。在一些实施例中,用于抽取血液样本的方法(例如注射器或手指点刺)可以比用于获得组织活检的程序(这可能需要手术)侵入性更小。提取的样本可以包括cfDNA和/或ctDNA。针对健康的个体,人体可以自然地清除cfDNA和其他细胞碎片。如果受试者患有癌症或疾病,则提取的样本中的ctDNA可以以可检测的水平存在以进行诊断。
在步骤120中,测序文库被制备。在文库制备期间,核酸样本被随机切割为数千或数百万个碎片。唯一分子标识符(UMI)通过接合体连接(adapter ligation)添加到核酸碎片(例如DNA碎片)。UMI是在接合体连接期间添加到DNA碎片末端的短核酸序列(例如4至10个碱基对)。在一些实施例中,UMI是简并碱基对,其用作可以被用于标识源自具体DNA碎片的序列读数的独特标签。在接合体连接后的PCR扩增期间,UMI与附接的DNA碎片一起被复制,这提供了一种在下游分析中标识来自同一原始碎片的序列读数的方式。
在步骤130中,靶向DNA序列从文库中富集。在富集期间,杂交探针(在本文中也称为“探针”)被用于靶向和下拉核酸碎片,这些核酸碎片可以提供关于癌症(或疾病)的存在或不存在、癌症状况或癌症分类(例如癌症类型或起源组织)的信息。针对给定的工作流程,探针可以被设计为与DNA或RNA的靶标(互补)链重组(或杂交)。靶标链可以是“正”链(例如转录为mRNA并且随后翻译为蛋白质的链)或互补的“负”链。探针的长度范围可以是10、100或1000个碱基对。在一些实施例中,探针是基于基因板设计的,以分析怀疑对应于某些癌症或其他类型疾病的(例如人类或另一生物体的)基因组的特定突变或靶标区域。而且,探针可以覆盖靶标区域的重叠部分。通过使用靶向基因板而不是对基因组的所有表达基因进行测序,也称为“全外显子组测序”,方法100可以被用于增加靶标区域的测序深度,其中深度指样本内的给定靶标序列已被测序的次数的计数。增加测序深度减少所需的核酸样本输入量。在杂交步骤之后,杂交的核酸碎片被捕获,并且也可以使用PCR进行扩增。
在步骤140中,从富集的DNA序列生成序列读数。测序数据可以通过本领域已知方式从富集的DNA序列中获取。例如,方法100可以包括下一代测序(NGS)技术,包括合成技术(依诺米那)、焦磷酸测序(454生命科学)、离子半导体技术(离子激流测序)、单分子实时测序(太平洋生物科学)、连接测序(SOLiD测序)、纳米孔测序(牛津纳米孔技术)或双端测序。在一些实施例中,大规模并行测序使用具有可逆染料终止剂的合成测序执行。
在一些实施例中,序列读数可以使用本领域已知的方法与参考基因组比对,以确定比对位置信息。比对位置信息可以指示参考基因组中对应于给定序列读数的开始核苷酸碱基和结束核苷酸碱基的区域的开始位置和结束位置。比对位置信息还可以包括序列读数长度,其可以从开始位置和结束位置确定。参考基因组中的区域可以与基因或基因碎片相关联。
在各种实施例中,序列读数由表示为R1和R2的读数对组成。例如,第一读数R1可以从核酸碎片的第一末端测序,而第二读数R2可以从核酸碎片的第二末端测序。因此,第一读数R1和第二读数R2的核苷酸碱基对可以与参考基因组的核苷酸碱基一致地比对(例如以相反的定向)。从读数对R1和R2获取的比对位置信息可以包括参考基因组中对应于第一读数(例如R1)的末端的开始位置和参考基因组中对应于第二读数(例如R2)的末端的结束位置。换句话说,参考基因组中的开始位置和结束位置表示核酸碎片对应于的参考基因组内的可能位置。具有SAM(序列比对图)格式或BAM(二进制)格式的输出文件可以被生成并且输出以供进一步分析,诸如变体调用,如下面关于图2描述的。
III.示例处理设备
图2是根据各种实施例的用于处理序列读数的处理系统200的框图。处理系统200包括序列处理器205、序列数据库210、模型数据库215、机器学习引擎220、模型225(例如对应于不同读数层的贝叶斯层次模型)、参数数据库230、评分引擎235和变体调用器240。图3是根据各种实施例的用于确定序列读数的变体的方法300的流程图。在一些实施例中,处理系统200执行方法300,以基于输入测序数据执行变体调用(例如,针对SNV和/或插入/缺失)。进一步地,处理系统200可以从与使用上述方法100制备的核酸样本相关联的输出文件获得输入测序数据。方法300包括但不限于以下步骤,这些步骤是相对于处理系统200的组件进行描述的。在一些实施例中,方法300的一个或多个步骤可以由用于例如使用变体调用格式(VCF)(诸如HaplotypeCaller、VarScan、Strelka或SomaticSniper)生成变体调用的不同过程的步骤取代。
在步骤300中,序列处理器205折叠输入测序数据的序列读数。在一些实施例中,折叠序列读数包括使用UMI以及可选地来自输出文件的测序数据(例如来自图1所示的方法100)的比对位置信息以将多个序列读数折叠为共有序列,以确定核酸碎片或其部分的最可能的序列。因为UMI通过富集和PCR与连接的核酸碎片一起复制,序列处理器205可以确定某些序列读数源自核酸样本中的相同分子。在一些实施例中,具有相同或类似的比对位置信息(例如阈值偏移内的开始和结束位置)并且包括共同UMI的序列读数被折叠,并且序列处理器205生成折叠读数(在本文中也称为共有读数)来表示核酸碎片。如果对应的折叠读数对具有共同UMI,则序列处理器205将共有读数指定为“双工的(duplex)”,这指示起始核酸分子的正链和负链都被捕获;否则,折叠的读数被指定为“非双工的”。在一些实施例中,序列处理器205可以执行其他类型的序列读数的纠错,作为折叠序列读数的备选或补充。
在步骤305中,序列处理器205基于两个或多个序列读数之间的重叠核苷酸序列的部分缝合折叠的读数。在一些实施例中,序列处理器205比较第一读数和第二读数之间的核苷酸序列,以确定第一读数和第二读数的核苷酸碱基对是否重叠。这两个序列读数也可以与参考基因组进行比较。在示例用例中,响应于确定第一读数和第二读数之间的(例如给定数量的核苷酸碱基的)重叠大于阈值长度(例如核苷酸碱基的阈值数量),序列处理器205将第一读数和第二读数指定为“缝合的”;否则,折叠的读数被指定为“未缝合的”。在一些实施例中,如果重叠大于阈值长度并且如果重叠不是滑动重叠,则缝合第一读数和第二读数。例如,滑动重叠可以包括均聚物运行(例如单个重复的核苷酸碱基)、二核苷酸运行(例如二核苷酸碱基序列)或三核苷酸运行(例如三核苷酸碱基序列),其中均聚物运行、二核苷酸运行或三核苷酸运行至少具有碱基对的阈值长度。
在步骤310中,序列处理器205将读数组装为路径。在一些实施例中,序列处理器205组装读数以生成针对靶标区域(例如基因)的有向图,例如德布莱英图。有向图的单向边缘表示靶标区域中的k个核苷酸碱基(在本文中也称为“k-mers”)的序列,并且边缘由顶点(或节点)连接。序列处理器205将折叠的读数与有向图进行比对,使得任何折叠的读数可以由边缘和对应顶点的子集按顺序表示。
在一些实施例中,序列处理器205确定描述有向图的参数集,并且处理有向图。附加地,该参数集可以包括从折叠读数到由有向图中的节点或边缘表示的k-mer的成功比对的k-mer的计数。序列处理器205例如在序列数据库210中存储有向图和对应的参数集,这可以被取回以更新图形或生成新图形。例如,序列处理器205可以基于参数集生成有向图的压缩版本(例如或修改现有图形)。在示例用例中,为了过滤掉具有较低重要性水平的有向图的数据,序列处理器205移除(例如“修整”或“修剪”)计数小于阈值的节点或边缘,并且维持计数大于或等于阈值的节点或边缘。
在步骤315中,变体调用器240从由序列处理器205组装的路径生成候选变体读数。变体可以对应于SNV或插入/缺失。在一些实施例中,变体调用器240可以通过将有向图(其可能已经通过在步骤310中修剪边缘或节点而被压缩)与基因组的靶标区域的参考序列进行比较来生成候选变体读数。变体调用器240可以将有向图的边缘与参考序列比对,并且记录错配边缘的基因组位置和与边缘相邻的错配核苷酸碱基作为候选变体的位置。附加地,变体检出器240可以基于靶标区域的测序深度生成候选变体读数。具体地,变体调用器240可以更有信心标识具有更大测序深度的靶标区域中的变体,例如因为更大数量的序列读数有助于解决(例如,使用冗余)错配或序列之间的其他碱基对变化。
在一些实施例中,变体读数可以基于变体读数的质量分类为不同的读数层。与折叠序列的重叠和/或互补位置相比,变体读数的质量可以对应于潜在变体位置的位置。在大规模变形测序中的样本制备(例如文库制备过程)中,受试者个体的核酸样本可以在并行测序被执行之前随机切割。核酸序列的相同副本可以被不同且随机地切割。因此,一些富集碎片可以具有可以与其他富集碎片缝合的重叠区域,而其他富集碎片则没有。一些富集碎片还可以具有同样富集的互补序列,从而在序列处理中生成双链碎片。因此,不同序列位置的变体读数可以对应于不同的质量。例如,在碎片的两个互补链都富集的位置处的变体读数通常比在仅从单链碎片获得支持的另一位置处的另一变体读数具有更好的质量。变体读数的读数层的细节将在图8A至8B中进一步讨论。
在一些实施例中,变体调用器240使用模型225生成候选变体读数,以确定来自受试者的序列读数的预期噪声率。模型225中的每个模型225可以是贝叶斯层次模型。贝叶斯层次模型可以是许多可能的模型架构中的一个架构,其可以被用于生成候选变体,并且它们彼此相关,因为它们都对特定于位置的噪声信息进行建模,以提高变体调用的灵敏度/特异性。更具体地,机器学习引擎220使用来自健康个体的样本来训练模型225,以对序列读数的每个位置的预期噪声率进行建模。在一些实施例中,对应于不同读数层的变体读数可以由不同模型不同地处理,每个模型特定于特定读数层。每个模型的结果可以被组合以生成组合结果。对读数层和模型进行分层的细节在下面将与图8A至11B相关联地进一步讨论。
进一步地,多个不同的模型可以被存储在模型数据库215中或取回用于应用后期训练。例如,第一模型被训练以对SNV噪声率进行建模,并且第二模型被训练以对插入/缺失噪声率进行建模。进一步地,评分引擎235可以使用模型225的参数来确定序列读数中的一个或多个真阳性的可能性。评分引擎235可以基于可能性确定质量分数(例如在对数尺度上)。例如,质量分数是Phred质量分数Q=-10·log10P,其中P是不正确的候选变体调用的可能性(例如假阳性)。
在步骤320中,评分引擎235基于模型225或真阳性的对应可能性或质量分数对变体读数进行评分。模型225的训练和应用在下面更详细地描述。
在步骤325,处理系统200输出相对于变体的分析结果。在一些实施例中,处理系统200输出一些或所有确定的候选变体连同对应的分数。下游系统(例如在处理系统200或处理系统200的其他组件外部的系统)可以将变体和分数用于各种应用,包括但不限于预测癌症、疾病或种系突变的存在。
IV.示例模型
图4是根据各种实施例的贝叶斯层次模型225的应用的图。出于解释的目的,突变A和突变B被示出为示例。在包括在图4中的本公开的大部分内容中,突变被表示为SNV,尽管在一些实施例中,本公开中的描述也适用于插入/缺失或其他类型的突变。第一样本的第一变体读数对应于示例突变A,其是第一参考等位基因的位置4处的C>T突变。第一样本具有为10的第一AD和为1000的第一总测序深度。第二样本的第二变体读数对应于示例突变B,其是第二参考等位基因的位置3处的T>G突变。第二样本具有为1的第二AD和为1200第二总深度。仅基于AD(或AF),突变A可能看起来是真阳性,而突变B可能看起来是假阳性,因为前者的AD(或AF)大于后者。然而,突变A和B可以具有每个等位基因和/或等位基因的每个位置的不同相对水平的噪声率。例如,一旦考虑了这些不同位置的相对噪声水平,突变A可能是假阳性,并且突变B可能是真阳性。本文描述的模型225对该噪声进行建模,以便相应地适当标识真阳性。
图4所图示的概率质量函数(PMF)指示来自受试者的样本在某个位置具有给定AD计数的概率(或可能性)。使用来自参考个体样本的测序数据(例如存储在序列数据库210中)作为训练数据集,处理系统200训练模型225,参考样本的PMF可以从中导出。参考个体可以是未知的或被发现与具体变体位置的突变相关联的个体,并且有时可以被称为健康个体,尽管健康个体可以与另一变体位置的不同突变相关联,其与针对具体变体位置训练的模型无关。PMF基于λp和rp,λp对(例如参考个体的)正常组织中的每个位置的每个等位基因的预期平均AD计数进行建模,rp对该AD计数中的预期变化(例如分散)进行建模。换句话说,λp和/或rp表示正常组织的测序数据中基于每个位置每个等位基因的基线噪声水平。在一些实施例中,参考个体的测序数据可以被分层到不同的读数层中,使得分别对应于特定读数层的多个模型225被训练。对应于读数层的每个模型可以具有不同的λp和不同的rp
使用图4的示例以进一步说明,来自参考个体的样本表示由yi建模的人口子集,其中i是训练集中的健康个体的索引。出于示例,假设模型225已经被训练,由模型225产生的PMF在视觉上说明了每个突变的测量AD的可能性,因此提供了哪些是真阳性以及哪些是假阳性的指示。与突变A相关联的图4左侧的示例PMF指示针对位置4的突变,第一样本的AD计数为10的概率约为20%。附加地,与突变B相关联的右侧的示例PDF指示针对位置3的突变,第二样本的AD计数为1的概率约为1%(注意:图4的PMF不完全按比例)。因此,与PMF的这些概率相对应的噪声率指示,尽管突变B具有较低的AD和AF,但突变A比突变B更可能发生。因此,在该示例中,突变B可以是真阳性,并且突变A可以是假阳性。因此,处理系统200可以通过使用模型225以更准确的速率区分真阳性和假阳性来执行改进的变体调用,并且进一步提供关于这些可能性的数值置信度。
图5A示出了根据各种实施例的用于确定真正的单核苷酸变体的贝叶斯层次模型225的参数和子模型之间的相关性。图5A所示的贝叶斯层次模型可以对应于变体读数的具体读数层。模型的参数可以被存储在参数数据库230中。在图5A所示的示例中,
Figure BDA0003632395070000163
表示分配给每个混合分量的权重向量。向量
Figure BDA0003632395070000164
采用K维单纯形内的值,并且可以在训练期间经由后验采样来学习或更新。针对这种训练,可以在所述单纯形上给予统一的先验。位置p所属的混合分量可以通过潜变量zp使用一个或多个不同的多项分布建模:
Figure BDA0003632395070000161
潜变量zp、混合分量向量
Figure BDA0003632395070000162
α和β共同允许μ的模型(即,贝叶斯层次模型225的子模型)具有参数,其“集用(pool)”关于噪声的知识,即,它们表示多个位置的噪声特点的相似性。因此,序列读数的位置可以由模型集用或者分组为潜在类别。同样有利地,任何这些“集用的”位置的样本都可以帮助训练这些共享参数。这样做的益处是处理系统200可以确定健康样本中的噪声模型,即使很少或没有先前观察到给定位置的交替等位基因(例如在用于训练模型的健康组织样本中)直接证据表明。
协变量xp(例如预测器)编码关于位置p的已知上下文信息,其可以包括但不限于诸如三核苷酸上下文、可映射性、片段重复或与序列读数相关联的其他信息等信息。三核苷酸上下文可以基于参考等位基因,并且可以被分配数字(例如整数)表示。例如,“AAA”被分配1,“ACA”被分配2,“AGA”被分配3等。可映射性表示读数与基因组的特定靶标区域的比对的唯一性水平。例如,可映射性被计算为序列读数将唯一映射的(多个)位置数量的倒数。片段重复对应于长核酸序列(例如长度大于约1000个碱基对),其几乎相同(例如,大于90%匹配)并且由于自然复制事件(例如,与癌症或疾病不关联)而出现在基因组中的多个位置。
位置p处的SNV的预期平均AD频率由参数μp建模。在一些实施例中,参数μp对应于每个测序深度
Figure BDA0003632395070000171
的平均AD计数。因为SNV是变体的示例,参数μp也可以被称为平均变体频率。为了在该描述中清晰起见,术语μp和yp指的是贝叶斯层次模型225的位置特定子模型。在一些实施例中,μp被建模为具有形状参数
Figure BDA0003632395070000172
和速率参数
Figure BDA0003632395070000173
的伽马分布随机变量:
Figure BDA0003632395070000174
在一些实施例中,其他函数可以被用于表示μp,其示例包括但不限于:具有对数均值
Figure BDA0003632395070000175
和对数标准偏差
Figure BDA0003632395070000176
的对数正态分布、威布尔分布、幂律、指数调制幂律,或前述的混合。形状参数
Figure BDA0003632395070000177
有时可以是分布中的分散参数rp的示例。
分布的方差可以由平均变体频率μp和分散参数rp确定。例如,在伽马分布的情况下,方差vp可以通过以下确定:
Figure BDA0003632395070000181
Lambda λp可以是平均变体计数,其可以通过μp乘以测序深度
Figure BDA0003632395070000182
确定。而且,lambda λp可以通过以下与形状参数
Figure BDA0003632395070000183
和速率参数
Figure BDA0003632395070000184
相关:
Figure BDA0003632395070000185
在图5A所示的示例中,形状和速率参数分别取决于协变量xp和潜变量zp,尽管在一些实施例中,相关性可以基于训练期间信息集用的不同程度而不同。例如,模型可以交替地构建,使得
Figure BDA0003632395070000186
取决于潜变量而不是协变量。(健康个体的)人口样本i中的位置p的SNV的AD计数分布由随机变量
Figure BDA0003632395070000187
建模。随机变量
Figure BDA0003632395070000188
也可以被称为变体计数或观察到的变体计数。在一些实施例中,给定样本在位置处的测序深度dip,该分布是泊松分布:
Figure BDA0003632395070000189
在一些实施例中,其他函数可以被用于表示
Figure BDA00036323950700001810
其示例包括但不限于:负二项式、康威-麦克斯韦-泊松分布、正态分布和零膨胀泊松。例如,随机变量
Figure BDA00036323950700001811
可以用负二项分布建模:
Figure BDA00036323950700001812
平均变体频率μp、平均变体计数
Figure BDA00036323950700001813
和分散参数rp可以被称为噪声参数,因为这些参数影响变体计数
Figure BDA00036323950700001814
的随机变量分布。
图5B示出了根据各种实施例的用于确定真实插入或缺失的贝叶斯层次模型的参数和子模型之间的相关性。与图5A所示的SNV模型相反,图5B所示的插入/缺失的模型包括层次的不同水平。协变量xp编码位置p处的已知特征,并且可以包括例如到均聚物的距离、到RepeatMasker重复的距离或与先前观察到的序列读数相关联的其他信息。潜变量
Figure BDA00036323950700001815
可以通过基于向量
Figure BDA00036323950700001816
的参数的狄利克雷分布建模,该向量表示某个位置的插入/缺失长度分布并且可以基于协变量。在一些实施例中,
Figure BDA0003632395070000191
还在共享(多个)相同协变量值的位置
Figure BDA0003632395070000192
之间共享。因此,例如潜变量可以表示以下信息,诸如均聚物插入/缺失出现在距锚定位置的位置1、2、3等碱基对处,而三核苷酸插入/缺失出现在距锚定位置的位置3、6、9等处。
位置p处的预期平均总插入/缺失频率由分布μp建模。在一些实施例中,参数μp对应于每个测序深度
Figure BDA0003632395070000193
的平均插入/缺失计数。因为插入/缺失是变体的示例,参数μp也可以被称为平均变体频率。在一些实施例中,分布基于协变量,并且具有具有形状参数
Figure BDA0003632395070000194
和速率参数
Figure BDA0003632395070000195
的伽马分布:
Figure BDA0003632395070000196
在一些实施例中,其他函数可以被用于表示μp,其示例包括但不限于:负二项式、康威-麦克斯韦-泊松分布、正态分布和零膨胀泊松。形状参数
Figure BDA0003632395070000197
有时可以是分布中的分散参数rp的示例。
分布的方差可以由平均变体频率μp和分散参数rp确定。例如,在伽马分布的情况下,方差vp可以通过以下确定:
Figure BDA0003632395070000198
Lambda λp可以是平均变体计数,其可以通过μp乘以测序深度
Figure BDA0003632395070000199
确定。而且,lambda λp可以通过以下与形状参数
Figure BDA00036323950700001910
和速率参数
Figure BDA00036323950700001911
相关:
Figure BDA00036323950700001912
(健康个体的)人口样本i中的位置p的观察到的插入/缺失由分布
Figure BDA00036323950700001913
建模。随机变量
Figure BDA00036323950700001914
也可以被称为变体计数或观察到的变体计数。与图5A中的示例类似,在一些实施例中,给定样本在位置的测序深度
Figure BDA00036323950700001915
插入/缺失强度的分布是泊松分布:
Figure BDA00036323950700001916
在一些实施例中,其他函数可以被用于表示
Figure BDA00036323950700001917
其示例包括但不限于:负二项式、康威-麦克斯韦-泊松分布、正态分布和零膨胀泊松。例如,在一些示例中,可以通过负二项分布对随机变量
Figure BDA0003632395070000201
建模:
Figure BDA0003632395070000202
平均变体频率μp、平均变体计数
Figure BDA0003632395070000203
和分散参数rp可以被称为噪声参数,因为这些参数影响变体计数
Figure BDA0003632395070000204
的随机变量的分布。
因为插入/缺失可以有不同的长度,所以插入/缺失模型中存在附加长度参数,而SNV的模型中不存在该参数。结果,图5B所示的示例模型具有附加的层次水平(例如,另一子模型),这在上面讨论的SNV模型中同样不存在。在样本i中的位置p处观察到的长度为l的插入/缺失计数(例如最多100个或多个碱基对的插入或缺失)由随机变量
Figure BDA0003632395070000205
建模,该变量表示以参数为条件的噪声下的插入/缺失分布。给定样本的插入/缺失强度
Figure BDA0003632395070000206
以及插入/缺失长度
Figure BDA0003632395070000207
在该位置的分布,该插入/缺失分布可以是多项式:
Figure BDA0003632395070000208
在一些实施例中,狄利克雷多项式函数或其他类型的模型可以被用于表示
Figure BDA0003632395070000209
通过以这种方式构架模型,机器学习引擎220可以将插入/缺失强度(即,噪声率)的学习与插入/缺失长度分布的学习分离。独立地确定对健康样本中是否会出现插入/缺失的预期和对某个位置的插入/缺失长度的预期的推断可以提高模型的灵敏度。例如,相对于基因组中的多个位置或区域的插入/缺失强度,长度分布可能更稳定,反之亦然。
图6A至6B图示了根据各种实施例的与贝叶斯层次模型225相关联的图。图6A所示的图形描绘了噪声率的分布μp,即,由模型表征的给定位置的SNV或插入/缺失的可能性(或强度)。连续分布表示基于来自健康个体的观察到的健康样本的训练数据(例如从序列数据库210取回)的非癌症或非疾病突变(例如健康组织中自然发生的突变)的平均变体频率μp。虽然在图6A中未示出,在一些实施例中,形状和速率参数可以基于其他变量(诸如协变量xp或潜变量zp)。图6B所示的图形描绘了在给定样本的给定参数(诸如给定位置处的测序深度dp)的情况下受试者的样本的给定位置处的AD分布。μp的绘图的离散概率是根据基于预期平均分布μp的人口的预测真实平均AD计数确定的。
图7A是根据各种实施例的用于通过拟合贝叶斯层次模型225来确定参数的示例过程的图。为了训练模型,机器学习引擎220从一组位置中的每个位置的预期噪声率的后验分布(例如图6B所示的图形)迭代地采样。机器学习引擎220可以使用马尔可夫链蒙特卡洛(MCMC)方法进行采样,例如Metropolis-Hastings(MH)算法、自定义MH算法、吉布斯采样算法、基于哈密顿力学的采样、随机采样以及其他采样算法。在贝叶斯推理训练期间,参数从联合后验分布中抽取以迭代更新模型的所有(或一些)参数和潜变量(例如
Figure BDA0003632395070000211
zp
Figure BDA0003632395070000212
Figure BDA0003632395070000213
μp等)。
在一些实施例中,机器学习引擎220通过将μp的绘图存储在参数数据库230中来执行模型拟合。如先前描述的,模型通过后验采样训练或拟合。在一些示例中,μp的绘图被存储在矩阵数据结构中,该矩阵数据结构具有在采样的位置集合中的每个位置的行和在来自(例如以观察数据为条件的所有参数的)联合后验的每个绘图的列。行数R可以大于600万,并且N次样本迭代的列数可以是数千。在一些实施例中,行和列指定不同于图7A所示的实施例,例如每行表示来自后验样本的绘图,并且每列表示采样位置(例如图7A所示的矩阵示例的转置)。
图7B是根据各种实施例的使用来自贝叶斯层次模型225的参数来确定假阳性的可能性的图。机器学习引擎220可以将图7A所示的R行×N列矩阵减少为图7B所图示的R行×2列矩阵。在一些示例中,机器学习引擎220确定各种噪声参数,诸如跨后验样本μp的每个位置的分散参数rp(例如,形状参数)和平均参数λp,其也可以被称为平均速率参数λp。分散参数rp可以被确定为
Figure BDA0003632395070000214
其中λp和vp分别为该位置处的μp的采样值的均值和方差。本领域技术人员将了解,也可以使用用于确定rp的其他函数,诸如最大似然估计。不同的噪声参数可以针对不同的读数层而被确定。例如,每个读数层可以有不同的λp和rp值。
在给定速率参数的情况下,机器学习引擎220还可以对约化矩阵中的分散参数执行分散重新估计。在一些实施例中,在贝叶斯训练和后验近似之后,机器学习引擎220通过基于每个位置的负二项式最大似然估计量对分散参数
Figure BDA0003632395070000221
重新训练来执行分散重新估计。在重新训练期间,速率参数可以保持固定。在一些实施例中,机器学习引擎220为训练数据的原始AD计数(例如基于参考样本的
Figure BDA0003632395070000222
Figure BDA0003632395070000223
通过读数层分层)确定每个位置处的分散参数r′p。机器学习引擎220确定
Figure BDA0003632395070000224
并且将
Figure BDA0003632395070000225
存储在约化矩阵中。本领域技术人员将了解,用于确定
Figure BDA0003632395070000226
的其他函数也可以被使用,诸如矩估计量、后验均值或后验模式的方法。
在训练模型的应用期间,处理系统200可以访问分散(例如形状)参数
Figure BDA0003632395070000227
和速率参数λp以确定由
Figure BDA0003632395070000228
和λp参数化的函数。该函数可以被用于确定受试者的新样本的后验预测概率质量函数(或概率密度函数)。基于给定位置的某个AD计数的预测概率,处理系统200可以在从样本中检测真阳性时考虑序列读数的每个位置的位点特异性噪声率。返回参照相对于图4描述的示例用例,突变A和B所示的PMF可以使用来自图7B的约化矩阵的参数来确定。后验预测概率质量函数可以被用于确定突变A或B的样本在某个位置具有AD计数的概率。
贝叶斯层次模型和用于对贝叶斯层次模型中的各种参数进行建模的分布可以针对变体读数的不同读数层进行单独训练。例如,每个读数层都可以有自己的贝叶斯层次模型,该模型具有自己的参数,诸如
Figure BDA0003632395070000229
μp等。
关于对测序数据集的噪声水平建模的贝叶斯层次模型的训练和使用的更多信息,于2018年10月5日提交的名称为“Site-Specific Noise model for TargetedSequencing”的美国专利申请号16/153,593出于所有目的通过引用并入本文。
V.示例读数层
图8A包括图示了根据各种实施例的序列读数的不同类别或读数层的图。如本文设想的,序列读数可以与表示读数的不同质量水平的不同读数层相关联,由此质量水平可以基于相对于序列读数的重叠片段的变体位置。更高质量的读数层对应于更低的噪声水平或更低的误差率,并且较低质量的读数层对应于较高的噪声水平或较高的误差率。
要注意的是,在序列扩增过程(例如大规模并行测序)中,样本(例如个体)的一个或多个序列可以以伪随机方式切割为不同的碎片。在一些情况下,并非所有碎片都与UMI连接,因此一些碎片在连接的碎片被富集之前洗掉。因此,富集碎片在每次测序运行中至少部分是随机的。不同碎片之间的重叠程度可能会有所不同。例如,一些富集碎片可以具有可以与其他富集碎片缝合的重叠区域。一些富集碎片还可以具有富集的互补序列(例如正向和反向序列、正负序列、顶部和底部序列、5'到3'和3'到5'序列),从而生成整个序列读数的全部或部分的双链读数。因此,在一些示例中,不同序列位置的变体读数可以包括互补和/或重叠的序列读数以确认变体。因此,每个变体读数可以对应于不同的读数层质量。例如,在碎片的两个互补链富集的位置处的变体读数通常比在仅单个碎片富集的第二位置处的另一变体读数具有更好的质量。未被包括在重叠区域或互补区域内的位置处的变体读数可归因于噪声而不是受试者的样本中存在的实际变体的可能性增加。
图8A图示了读数层的四个不同示例。在一些实施例中,基于相对于序列读数内的重叠和互补位置的序列读数内的感兴趣的潜在变体位置,序列读数被分离为读数层。换言之,基于潜在变体位置读数是否被包括或以其他方式完全嵌入重叠区域(即,缝合区域)内以及变体位置是否被包括或以其他方式完全嵌入互补区域(即,双链区域、双工区域)内,序列读数被分类为四个读数层中的一个读数层。
通过示例,在图8A中,潜在变体位置被阴影化。第一示例读数层810包括同时具有双链(也称为“双工”或(互补的”)和缝合序列读数的变体读数。例如,至少两个5'到3'序列读数具有重叠区域,并且可以被缝合在一起。类似地,至少两个3'到5'序列读数具有重叠区域,并且也可以被缝合在一起。在示例第一读数层810中,潜在变体位置位于或者完全嵌入重叠或其他缝合区域内,因此序列读数包括缝合区域。同样地,5'到3'序列读数的至少一部分和3'到5'序列读数的一部分彼此互补,并且潜在变体位置位于互补区域内(例如潜在变体位置在公共重叠区域处完全嵌入顶部和底部序列读数内)。因此,除了包括缝合区域外,序列读数包括双链区域,并且潜在变体读数属于表示双链和缝合读数的第一读数层810。
在图8A中,第二示例读数层820包括位于双链但未缝合的序列读数部分内的变体读数。在第二读数层820中,5'到3'序列读数的一部分和3'到5'序列读数的一部分彼此互补,并且潜在变体位置位于互补区域内。因此,序列读数包括双链区域。然而,潜在变体位置未被包括在序列读数的任何重叠区域或其他缝合区域内。具体地,尽管两个5'到3'序列读数可以被缝合在一起这一事实,但这个示例分层是因为潜在变体位置未被包括在重叠区域或其他可缝合区域内。因此,虽然序列读数包括双链区域,但序列读数不包括缝合区域,并且潜在变体读数属于表示双链但未缝合读数的第二读数层820。
第三示例读数层830包括位于或以其他方式完全嵌入单链(例如非双工)缝合读数内的变体读数。在第三读数层830中,潜在变体位置被包括在两个或多个序列读数的重叠区域内,因此序列读数包括缝合区域。然而,序列读数是单链的,因为序列读数(诸如两个所图示的5'到3'序列读数)不包括互补区域(例如序列读数仅基于5'到3'链并且不受互补的3'到5'链的支持)。在一些情况下(未图示),一个或多个互补序列读数(例如3'到5'序列读数)可以在示例读数层3处找到,但不包括潜在变体位置。因此,潜在变体读数属于表示单链但缝合读数的第三读数层830。
进一步如图8A所示,第四示例读数层840包括位于单链和未缝合读数处的变体读数。与第三读数层830一样,第四读数层840表示单链读数,因为所图示的序列读数不包括包含变体位置的互补区域(或在一些情况下(未图示),还包括互补区域,但潜在变体位置未位于或完全嵌入互补区域内)。因此,第四读数层840表示未缝合读数,因为潜在变体位置未被包括在两个序列读数的重叠区域内。
在一些实施例中,样本的序列读数可以被分层到图8A所图示的四个读数层中。在一些实施例中,可以存在对应于最低质量的变体读数的附加的第五读数层,诸如单链和未缝合的序列读数,其包括靠近一个或多个序列读数末端的潜在变体位置。例如,如果单链和未缝合的序列读数包括位于距序列读数任一端的预定阈值碱基数量内(例如在约7个碱基内或约30个碱基内)的潜在变体位置,则序列读数可以被分类为第五读数层。在一些实施例中,图8A所示的四个读数层中的每个读数层可以被细分为两个子层:第一较低质量的子层,对应于序列读数,包括靠近一个或多个序列读数的任一端的潜在变体位置;以及第二较高质量的子层,包括距一个或多个序列读数的末端超过阈值距离的潜在变体位置。
图8B示出了根据各种实施例的图示图8A的读数层的不同质量的实验结果。更高质量的读数层对应于更低的误差率和/或更低的噪声水平。换言之,与被分层到高质量读数层的序列读数中的变体读数(例如在潜在变体位置检测到的SNV或插入/缺失)相比,分层到较低质量的读数层的序列读数中的变体读数更可能可归因于样本的实际突变,而不是随机事件(例如,由于噪声)。图8B是不同读数层t1至t5的参考样本(例如健康个体)的平均误差率的log 10的图表。层1(t1)指的是图8A中的第一读数层810,层2(t2)是指第二读数层820等。层5(t5)是指单链和未缝合读数层和/或位于序列读数任一端附近(例如与末端相距7个碱基内)的潜在变体位置。针对层1,其平均误差率mu的以10为底的对数约为-6.3到-7。换句话说,针对健康个体,每个测序深度的平均误差率约为10-6.3到10-7个变体读数。另一方面,针对层4,其平均误差率mu的以10为底的对数约为-4.7到-5.5。在一些方面中,图8B通常示出平均误差率在层1至4之间增加,从约1/1,000,000-1/10,000,000的层1平均误差率到接近约1/1,000,000的层2平均误差率,在层3中平均误差率再次增加到约为<1/1,000,000,并且仍会再次增加到约为1/100,000的层4误差率。因此,在层4检测到的变体读数是误差等位基因的可能性是在层1检测到的变体读数的大约100倍。换言之,第四读数层比第一读数层噪声相对更大且更容易出错。换句话说,第四读数层比第一读数层噪声相对更大且更容易出错。进一步地,如示例层5所示,没有平均误差率mu(或没有有意义的平均误差率),因为其序列读数质量是如此低质量的,以至于它们不显著或以其他方式从分析中丢弃。
附加地或备选地,序列读数可以通过其他分类方法分类为不同的读数层。例如,如果变体是SNV,则每个读数层可以基于核苷酸替换的类型进一步细分为12个附加子层(例如A>C、A>T、G>C等)(参见例如下面讨论的图8H)。因为在SNV中有四个核苷酸并且每个核苷酸都被不同的核苷酸替换,所以总共有12种不同类型的SNV。
图8C到8G示出了当序列读数通过图8A中描述的方式首先被分层到读数层并且基于核苷酸替代的类型进一步分层到十二个子读数层时图8A的读数层的实验结果。具体地,图8C至8G图示了针对某个位置的交替读数的误差分布的统计模型(例如,负二项式)的平均误差率mu(μ)和大小参数,使得每个位置(例如,每个点)的误差分布以在给定样本中看到的实际读数深度为条件。该模型由不同类别的读数(例如,层)分层,使得图8C图示了第一读数层(即,双链和缝合读数)的不同类型的核苷酸替换的结果。图8D图示了第二读数层(即,双链但未缝合的读数)的不同类型的核苷酸替换的结果。图8E图示了第三读数层(即,单链和缝合读数)的不同类型的核苷酸替换的结果,并且图8F图示了第四读数层的不同类型的核苷酸替换的结果。图8G图示了来自第五读数层的结果,对应于超过所示轴的最低质量读数和高误差率。要注意的是,针对图8C至8G,在图表顶部水平的核苷酸碱基A、C、G和T是指交替碱基,而竖直沿着图表右侧的核苷酸碱基A、C、G和T是指参考碱基。
参照图8C,不同SNV中的典型变体频率的十二个不同分布示出了第一读数层可以被进一步划分为十二个子层,以进一步改进噪声模型。行对应于原始核苷酸,而列对应于改变的核苷酸。例如,第三行第一列的单元格可以对应于G到A的SNV。实验示出,C到A的子读数层(即,第二行第一列)(其μ分布集中在以10为底的对数标度的-7到-8范围内)可能比其μ分布散布在以10为底的对数标度的-5到-7范围内的T到C的子读数层(即,第四行第二列)的噪声更小。
比较图8C到图8G之间的差异,对数标度中的平均误差率μ的分布示出随着读数层从第一读数层改变到第五读数层而向零移位(即,μ变得更大)并且最终超过零(即,在图8G处)的总体趋势。例如,关注T到G的子读数层(即,第四行第三列),对数标度的μ的分布从第一读数层中的-6和-7之间移位到第四读数层中的-4和-5之间。因此,图8C至图8G表明,随着读数层变得更加嘈杂,平均误差率μ变得更高。
现在转到图8H,示出了根据本文描述的各种实施例的实验结果,图示了在跨图8A的读数层t1至t4取得的具体SNV核苷酸替换处的不同平均误差率μ。具体地,图8B是在读数层t1至t4上观察到的不同SNV的对数标度的参考样本(例如健康个体)的平均误差率μ的log10的图表。要注意的是,本文描述的误差分布的统计模型可以由不同类别的读数(层)和/或进一步由不同的SNV分层,如图8H所示。
VI.使用分层读数的示例数据处理
图9是描绘了根据一些实施例的使用分层噪声模型分析样本的DNA数据集的过程的流程图。该过程可以被用于处理样本的DNA测序数据集,诸如包括cfDNA的个体样本,以生成质量分数,表示个体在潜在变体位置具有变体的可能性。由该过程确定的质量分数越高,变体读数就越有可能是来自实际突变而不是噪声的结果。
在步骤910中,处理系统可以访问由DNA测序生成的DNA测序数据集。例如,DNA测序可以是一种大规模并行DNA测序,诸如下一代测序(NGS)。DNA测序数据集包括多个已处理的序列读数,其包括感兴趣的变体位置(例如DNA序列中的具体基因位置)。至少一些已处理的序列读数可以通过折叠并且缝合诸如通过图3中描述的过程生成的DNA测序中的原始序列读数而生成。例如,NGS的典型运行可能会生成数百万甚或数十亿的序列读数。一些原始序列读数可以被包括在包括感兴趣的变体位置的基因座中。反过来,原始序列读数可以通过折叠和缝合进行处理,以生成已处理的序列读数。要注意的是,虽然在本示例中描述了DNA测序,但RNA测序也可以被实施用于本文的分析。
包括感兴趣的变体位置的已处理序列读数可以具有不同的碱基对长度和不同程度的重叠和/或互补。在步骤920中,处理系统可以将多个已处理的序列读数分层到不同的读数层。不同的读数层可以基于序列读数的质量进行分层。例如,已处理的序列读数可以基于变体位置是否被包括在重叠区域中和/或被包括在互补区域中来分层,如结合图8A讨论的。对已处理的序列读数进行分层的其他方法也是可能的。例如,已处理的序列读数还可以基于核苷酸替换的类型、变体位置是否接近序列末端等进行分层。在一些实施例中,不同的读数层包括至少四个读数层。在一些示例中,四个读数层是(1)双链的缝合读数层,(2)双链的未缝合读数层,(3)单链的缝合读数层,以及(4)单链的未缝合读数层。
在步骤930中,处理系统可以为每个读数层确定变体位置处的分层测序深度。针对每个读数层,分层测序深度可以是被分层到读数层的序列读数的测序深度。换句话说,分层测序深度可以是被分层到读数层的序列读数的总数。处理系统还可以确定每个读数层的实际变体计数。例如,针对读数层,大多数序列读数可能在变体位置不包含实际变体(无论是SNV还是插入/缺失)。在一些情况下,只有少数序列读数在变体位置包括实际变体。分层变体计数可以是特定读数层的实际变体计数的总数。
在步骤940中,处理系统可以为每个读数层确定以读数层的分层测序深度为条件的一个或多个噪声参数。噪声参数可以是特定于读数层的噪声模型的参数。例如,处理系统可以包括多个分层噪声模型,分别特定于读数层。分层噪声模型(或其中一些)可以对应于图5A到图7B中描述的贝叶斯层次模型。换句话说,在一些实施例中,每个读数层具有它自己的贝叶斯层次模型。读数层的每个噪声模型都可以使用DNA测序样本的不同训练集进行训练。通过示例,多个参考个体(诸如健康个体)的DNA测序数据集可以被收集。参考个体的数据集的已处理序列读数可以按读数层进行分层。每个读数层的分层已处理序列读数可以被用作DNA测序样本的分层训练集,以训练读数层的分层噪声模型。分层噪声模型的各种分布(例如结合图5A和5B讨论的伽马分布和泊松分布)可以基于分层训练集确定。
每个读数层的分层变体计数的概率分布可以通过噪声分布建模。分层变体计数的概率分布可以取决于使用的分布类型和定义噪声分布的一个或多个参数。例如,在讨论的贝叶斯层次模型的情况下,分层变体计数的分布可以对应于以两个参数为条件的后验分布。参数可以是以分层测序深度和分散参数为条件的分层平均变体计数。参数中的每个参数还可以对应于影响参数的一个或多个先验分布。例如,以分层测序深度为条件的分层平均变体计数可以通过伽马分布建模。因为先验分布可以描述参数的分布,所以先验分布也可以被称为参数分布。
针对每个读数层,处理系统可以通过将从受试者的数据集获得的分层测序深度输入到训练的噪声模型确定以分层测序深度为条件的一个或多个噪声参数。例如,训练后的噪声模型可以访问特定于读数层的参数分布(例如先验分布)。参数分布可以基于参考个体的分层训练集形成,并且可以描述分层训练集的分布。训练后的噪声模型可以使用参数分布来确定以对应于读数层的分层测序深度为条件的噪声参数。
虽然贝叶斯层次模型被用作噪声模型的示例,但在各种实施例中,不同类型的训练机器学习模型可以被用作噪声模型。而且,取决于所使用的模型,不同的噪声参数可以被用于模拟噪声分布。
在步骤950中,处理系统可以基于以读数层的分层测序深度为条件的一个或多个噪声参数生成特定于读数层的噪声模型的输出。输出的生成可以针对不同的读数层重复。取决于实施例,不同类型的输出可以被生成。例如,在一些实施例中,每个分层噪声模型在噪声参数被确定之后不执行进一步的计算。噪声模型的输出可以是以针对每层确定的分层测序深度为条件的一个或多个噪声参数。在负二项分布被用作噪声分布以对分层变体计数建模的情况下,噪声模型的输出可以是以分层测序深度和分散参数为条件的分层平均变体计数。在一些实施例中,在确定噪声参数之后,每个分层噪声模型可以生成后验分布。在这种实施例中,特定于读数层的噪声模型的输出可以是随后观察到的数据的读数层的变体计数大于或等于在个体受试者的DNA数据集中观察到的总变体计数可归因于噪声的可能性。其他合适的输出也是可能的。
在步骤960中,处理系统可以组合生成的噪声模型输出以产生组合结果。组合结果可以是受试者个体的DNA测序数据集的总体处理结果的表示。组合结果可以采用任何合适的形式。在一些实施例中,组合结果可以包括随后观察到的数据的总变体计数大于或等于在多个已处理的序列读数中观察到的总变体计数可归因于噪声的可能性。换句话说,可能性可以表示事件与在受试者个体的DNA数据集的多个已处理序列读数中观察到的总变体计数一样或更极端的可能性。在一些情况下,可能性可以对应于在零假设中使用的p值。分层噪声模型的输出可以如何被组合以生成组合结果可以取决于不同的实施例。在一些实施例中,将在图10中详细讨论的矩匹配方法可以被使用。在一些实施例中,将在图11A和11B中详细讨论的积分方法(例如求和方法)可以被使用。
在步骤970中,处理系统可以确定组合结果的质量分数。在一些示例中,诸如以可能性P(例如p值)的形式的组合结果可以被转换为Phred等级的质量分数,其中Q=-10·log10P。例如,Phred质量分数为20指示P=1/100的机会出现不正确的变体调用,并且Phred质量分数为60指示P=1/1,000,000的机会出现不正确的变体调用。因此,较高的Phred质量分数对应于检测到实际突变的更大置信度。质量分数可以被用于区分真阳性和假阳性。在一些实施例中,响应于质量分数高于预定阈值,处理系统可以指示个体在统计上可能在变体位置具有突变。
VII.用以组合分层输出的矩匹配
图10是描绘了根据一些实施例的用于使用矩匹配来组合用于不同读数层的分层噪声模型的输出的过程的流程图。图10中描绘的过程可以对应于图9中的步骤950和/或960。在步骤1010中,处理系统可以组合来自每个噪声模型输出的平均变体计数和变体计数的方差,以产生总体平均变体计数和总体分散参数。每个噪声模型的输出可以采用以分层测序读数为条件的噪声参数的形式。处理系统可以通过首先匹配每个读数层的各个矩以生成总矩来评估在给定观察到的总体测序读数的情况下不同读数层的总变体计数(例如总体p值)的总体可能性。处理系统可以使用总体矩来建模以观察到的总体测序读数为条件的总体分布。读数层的单个噪声分布可以是负二项分布。同样地,跨多个读数层的总体噪声分布也可以是与各个读数层的矩相匹配的负二项式分布。
步骤1010可以包括若干子步骤。针对每个读数层,处理系统可以确定分层测序深度。每层的噪声分布的第一矩和第二矩可以被用作定义噪声分布的噪声参数。在步骤1012中,基于分层测序深度,处理系统可以确定每个读数层的第一时刻(例如平均变体计数)。例如,在上面讨论的贝叶斯层次模型的情况下,特定读数层的变体频率可以被建模为具有形状参数
Figure BDA0003632395070000311
和速率参数
Figure BDA0003632395070000312
的伽马分布随机变量:
Figure BDA0003632395070000313
每个读数层可以具有其自己的形状参数和速率参数,这些参数是基于参考样本数据集确定的。因此,以分层测序深度为条件的每个读数层的变体频率可以不同。
处理系统可以通过将变体频率乘以分层测序深度来确定每层的第一时刻,即,分层平均变体计数λp
Figure BDA0003632395070000321
在步骤1014中,处理系统还可以确定每个读数层的第二矩,即,方差。在贝叶斯层次模型具有伽马分布变体频率的情况下,每个读数层的方差可以由平均变体计数λp和分散参数rp确定。例如,方差vp可以通过以下确定:
Figure BDA0003632395070000322
在步骤1016中,处理系统可以通过矩匹配确定总体平均变体计数(总体第一矩)和总体方差(总体第二矩)。在一些情况下,处理系统可以通过对不同读数层的矩求和来执行矩匹配以获得总体矩。例如,以总测序深度为条件的所有读数层的总体平均变体计数可以通过以下确定:
λall=∑λp
同样地,所有读数层的总体方差可以通过对每个读数层的方差求和来确定:
vall=∑vp
处理系统可以通过总体噪声分布对以总测序深度为条件的观察到的总变体计数的可能性进行建模。总体噪声分布可以是负二项分布,其由总体均值λall和总体分散参数rall参数化。总体分散参数可以由总体均值和总体方差确定:
Figure BDA0003632395070000323
在步骤1020中,处理系统可以使用由总体第一矩和总体第二矩建模的总体噪声分布来确定总体似然性。例如,随机变量
Figure BDA0003632395070000324
通过负二项分布建模为:
Figure BDA0003632395070000331
表示事件与以总测序深度为条件的观察到的总变体计数一样或更极端的可能性的随机变量
Figure BDA0003632395070000332
可以是处理系统的组合结果。在一些情况下,随机变量
Figure BDA0003632395070000333
可以被用作p值来验证或拒绝变体读数是由随机事件(例如噪声)引起的零假设。处理系统还可以应用负二项式尾概率来获得基于随机变量
Figure BDA0003632395070000334
的p值,并且可以确定Phred等级质量分数。
VII.组合分层输出的积分方法
图11A是描绘了根据一些实施例的用于使用积分方法组合不同读数层的分层噪声模型的输出以组合每个读数层的可能性的过程的流程图。图11A中描绘的处理可以对应于图9中的步骤950和/或960。在图11A所示的过程中,处理系统打算确定p值作为整个系统的组合结果。p值可以表示随后观察到的数据的总变体计数大于或等于在实际数据中观察到的总变体计数可归因于噪声的可能性。换句话说,p值可以代表与在零假设下观察到的受试者个体的序列读数一样或比其更极端的事件。在一些情况下,当事件的总变体计数大于观察到的变体计数时,事件可能比观察到的受试者个体的变体计数更极端(即,不太可能)。同样地,针对读数层,当事件的分层变体计数高于分层观察到的变体计数时,事件可能比分层观察到的变体计数更极端。具有较高变体计数的事件更为极端,因为变体读数通常是不寻常的。在变体位置观察到变体读数的机会通常明显低于在变体位置观察到非变体读数的机会。
图11B图示了根据各种实施例的多维空间中更极端事件的计数。为了简单起见,两个维度被示出,但是根据各种实施例的处理系统可以使用图11B所图示的原理来处理更高维度。图11B中的两个维度可以分别表示两个读数层。例如,第一读数层表示双链读数层,并且第二读数层表示单链读数层。
首先仅关注双链变体计数(x轴)的读数层,在示例中,观察到的来自受试者个体的分层变体计数为2。针对相同的读数层,分层变体计数为3的事件比实际观察到的分层变体计数更不可能(更极端),因为与潜在变体位置的非变体读数相比,潜在变体位置的变体读数不太可能。同样地,分层变体计数为4的另一事件比实际观察到的分层变体计数更不可能。换句话说,较不可能(更极端)事件的组合占据了大于观察到的变体计数的空间,并且一直跨越到无穷大。相反,具有分层变体计数1或0的事件比实际观察到的分层变体计数为2的可能性更大。
现在考虑两个读数层,可能存在可以被假设为可能相同或大致相同的观察到的分层变体计数的不同组合。在NGS样本制备中,受试者个体的核酸序列可以以部分随机方式切割。因此,一些已处理的序列读数可能不包括互补序列读数。因此,一些已处理的序列读数可以是单链序列读数。换句话说,针对相同的核酸序列样本,不同的NGS运行会在不同的读数层中产生序列读数的不同组合。基于一些比率,第一读数层的分层变体计数可以等同于第二读数层的分层变体计数。在一些实施例中,该比率被建模为预定值。例如,一个双链变体计数可以被认为等同于两个单链变体计数,尽管在一些实施例中,也可以使用除2之外的数字。
基于观察到的不同读数层的分层变体计数,图11B所示的图形中的坐标可以被划分为表示与不同读数层的观察到的变体计数一样或更极端的事件的点和表示比观察到的变体计数更不极端的事件的点。例如,假设受试者个体的实际观察到的变体计数为(1,2)(即,一个双链变体计数和两个单链变体计数)。坐标(0,4)和(2,0)可以被假设为与(1,2)的组合等效,并且表示与观察到的数据一样极端的事件。超出由(1,2)、(0,4)、(2,0)绘制的边界的所有坐标都可以被分类为比观察到的数据更极端。例如,坐标(3,3)可以被认为是更极端的情况。与观察到的数据相比,在边界内并且更靠近原点的所有坐标都可以被视为不太极端(更有可能)的情况。例如,坐标(1,1)、(0,2)、(1,0)等可以被分类为比观察到的数据不太极端。
处理系统的组合结果可以采用p值的形式,p值表示与观察到的数据一样或更极端的事件的可能性。处理系统可以通过将对应于所有坐标的概率求和来进行积分,这些坐标表示与观察到的数据一样或更极端的事件以确定p值。然而,因为坐标可以包括一直到无穷远的点,所以处理系统也可以计算p值的统计补数。换句话说,处理系统可以将对应于表示比观察到的数据较不极端的事件所有坐标的概率求和,以确定p值的补数。然后处理系统可以通过从1.0中减去补数来确定p值。在一些实施例中,处理系统可以使用对数标度的概率来实现数值稳定性,因为在计算机上添加浮点数可能在数值上是不稳定的。
返回参照11A,根据一些实施例,处理单元可以基于附图所示的过程来确定p值。在步骤1110中,处理单元可以在每个读数层中确定比读数层的观察到的分层变体计数更有可能的可能事件。这些事件可以定义多维框,其中坐标可以与比观察到的数据更大的可能性相关联。可能性可以在对数标度上。在步骤1120中,处理单元可以针对多个读数层中的一个读数层确定分层变体计数的组合中的每个组合是否对应于比观察到的数据更高或更低的可能性。在步骤1130中,处理单元可以标识与比观察到的读数层的分层变体计数更高的发生可能性相关联的可能事件的组合。处理单元可以针对读数层中的每个读数层重复步骤1120和1130。在步骤1140中,处理单元可以对所标识的组合的概率求和以确定统计补数。在步骤1150中,处理单元可以从1.0中减去统计补数以确定总体p值。p值可以是对应于图9中的步骤960的组合结果。
其他方式也可以用于确定总体p值。例如,尾概率技术可以被使用。在一些实施例中,积分方法可以由一种或多种机器学习模型代替。例如,随机森林回归模型可以被训练以确定Phred等级质量分数或来自训练样本数据集的p值。图11A中描述的积分过程可以被用于生成多个训练集样本。训练集样本可以被用于训练机器学习模型,使得模型可以被用于确定质量分数。
IX.实验结果
图12A和12B图示了在使用根据一些实施例的过程执行的实验中观察到的质量分数与默认质量分数的图表。在图12A中,特定个体#1355的样本的校准的模拟数据集被使用,并且利用图9和11A中描述的层级噪声模型过程进行分析。该个体在样本中发现了表示为“chr6”的变体。数据模拟变体读数的随机事件。每个模拟事件都可以具有某些变体读数,这些读数可以被分层到一个或多个读数层中。x轴表示通过使用模拟数据计算的模拟事件的实际质量分数的值。y轴表示由图11A中描述的过程确定的观察到的质量分数的值。结果示出,除了一些离散化之外,观察到的质量分数大部分落在对角线上。这指示图11A中描述的过程成功地确定了与观察到的数据一样或更极端的可能事件的可能性。
在图12B中,个体#1355的实际数据集被使用,并且利用图9和11A中描述的层级噪声模型过程分析。数据集可以包括各种潜在变体位置的数据。每个点对应于潜在的变体位置,并且质量分数是基于由图11A中描述的过程确定的过程来确定的。一开始存在一些零质量分数点,因为在一些情况下,实际数据集中的许多位置可能在任何读数层中都没有发现变体计数。因此,所有这些位置的Phred等级质量分数均等于零。其余的大部分点主要沿着对角线落下,因为个体在大多数序列位置基本上没有突变。换句话说,沿着对角线落下的那些位置的变体计数很大程度上可以归因于噪声。观察到的质量分数显著高于默认质量分数的异常值(例如图12B中的默认Phred分数约为55的点)可以指示可能存在可以被标注用于进一步评估的非噪声事件。
图13A和13B图示了将使用读数层的质量分数的结果与使用不区分任何序列读数或变体读数的读数层的噪声模型的质量分数进行比较的实验结果。图13A和13B中的y轴表示根据各种实施例的使用将序列读数分层到不同读数层中的方法确定的质量分数。图13A和13B中的x轴表示使用类似噪声模型确定的质量分数,但噪声模型不按读数层区分序列读数。图13A图示了使用图11A中描述的积分方法以确定质量分数的实验结果。结果示出,针对包括双工读数的序列读数(例如标记为“真”的较暗点),与仅包括单工读数的序列读数(例如标记为“假”的较亮点)相比,数据点向上移位。这指示读数层噪声模型提高了包括双链读数的序列读数的总体质量分数,因为双链读数通常包括比单链读数更多的证据。图13B图示了使用图10中描述的矩匹配方法以确定质量分数的实验结果。与图13A一样,在这种情况下,矩匹配方法还提高了包括双链读数的序列读数的质量分数。
X.变体标识
图14是根据各种实施例的描绘标识个体的潜在突变位置的过程的流程图。在步骤1410中,系统可以接收个体的DNA样本。在步骤1420中,系统可以执行DNA测序以生成已处理的序列读数。在步骤1430中,系统可以通过不同的变体位置分配已处理的序列读数。在步骤1440中,针对每个变体位置,系统可以将分配给变体位置的已处理序列读数分层到多个读数层中。在步骤1450,系统可以确定不同变体位置处的可能性的质量分数。每个质量分数可以使用对读数层进行分层的噪声模型基于上述过程来确定。在步骤1460中,系统可以标识质量分数高于预定阈值的变体位置。这些变体位置可以被标注,以进一步调查潜在的突变或潜在的诊断。
在步骤1470中,系统可以基于所标识的变体位置来生成疾病的诊断。在一些实施例中,可以指示某些癌症和/或用作某些疗法的生物标志物的变体或突变可以包括:ACVR1B、AKT3、AMER1、APC、ARID1A、ARID1B、ARID2、ASXL1、ASXL2、ATM、ATR、BAP1BCL2、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、EPHA5、EPHB1、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、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和TNFRSF14。在一些实施例中,癌症免疫疗法可以靶向OX40、LAG3和/或ICOS。
在步骤1480中,可以提供疾病的治疗。在提供治疗之前,伴随诊断操作还可以被执行。伴随诊断操作可以使用本文描述的过程来标识一个或多个准则,包括变体或突变。提供治疗可以采取引起或推荐医师向患者施用具体剂量的药物的形式。
例如,本文描述的系统和方法可以被用于检测作为癌症治疗(诸如某些免疫疗法和靶向疗法)的生物标志物的变体或突变。这种疗法可以包括例如免疫球蛋白、蛋白质、肽、小分子、纳米颗粒或核酸。在一些实施例中,疗法包括抗体或其功能碎片。在一些实施例中,抗体可以包括:
Figure BDA0003632395070000381
(利妥昔单抗)、
Figure BDA0003632395070000382
(曲妥珠单抗)、
Figure BDA0003632395070000383
(西妥昔单抗)、
Figure BDA0003632395070000384
(帕尼单抗)、
Figure BDA0003632395070000385
(奥法木单抗)、
Figure BDA0003632395070000386
(贝利尤单抗)、
Figure BDA0003632395070000387
(伊匹单抗)、
Figure BDA0003632395070000388
(帕妥珠单抗)、
Figure BDA0003632395070000389
(纳武单抗)、
Figure BDA00036323950700003810
(阿特珠单抗、MPDL3280A)、
Figure BDA00036323950700003811
CT-011、
Figure BDA00036323950700003812
(派姆单抗、MK-3475)、BMS-936559、MED14736、MSB0010718C、
Figure BDA00036323950700003813
(度伐鲁单抗)、
Figure BDA00036323950700003814
(阿维鲁单抗)和玛格妥昔单抗(MGAH22)。
在一些实施例中,免疫疗法和靶向疗法包括PD-1抑制、PD-L1抑制或CTL-4抑制。PD-1抑制靶向T细胞和其他免疫细胞上的程序性死亡受体。PD-1抑制免疫疗法的示例包括派姆单抗;可瑞达;纳武单抗;欧狄沃;西米普利单抗;利塔约。PD-L1抑制靶向由肿瘤和调节免疫细胞表达的程序性死亡受体配体。PD-L1抑制免疫疗法的示例包括阿特珠单抗;特善奇;阿维鲁单抗;巴文西亚;度伐鲁单抗;英飞凡。CTL-4抑制靶向T细胞活化。CTL-4抑制免疫疗法的示例包括伊匹单抗;易普利姆玛。
针对非小细胞肺癌适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括EGFR外显子19缺失和EGFR外显子21L858R改变(例如针对诸如
Figure BDA0003632395070000391
(阿法替尼)、
Figure BDA0003632395070000392
(吉非替尼)、
Figure BDA0003632395070000393
(奥希替尼)或
Figure BDA0003632395070000394
(厄罗替尼)等疗法);EGFR外显子20T790M改变(例如可以用
Figure BDA0003632395070000395
(奥希替尼)治疗);ALK重排(例如可以用
Figure BDA0003632395070000396
(艾乐替尼)、
Figure BDA0003632395070000397
(克唑替尼)或
Figure BDA0003632395070000398
(色瑞替尼)治疗);BRAFV600E(例如可以用
Figure BDA0003632395070000399
(达拉非尼)与
Figure BDA00036323950700003910
(曲美替尼)组合治疗);导致MET外显子14跳跃的单核苷酸变体(SNV)和插入/缺失(例如可以用TabrectaTM(卡马替尼)治疗)。
针对黑色素瘤适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括BRAF V600E(例如可以用
Figure BDA00036323950700003911
(达拉非尼)或
Figure BDA00036323950700003912
(维莫非尼)治疗);BRAFV600E或V600K(例如可以用
Figure BDA00036323950700003913
(曲美替尼)或
Figure BDA00036323950700003914
(考比替尼)与
Figure BDA00036323950700003915
(维莫非尼)组合治疗)。
针对乳腺癌适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括ERBB2(HER2)扩增(例如,其可以用
Figure BDA00036323950700003916
(曲妥珠单抗)、
Figure BDA00036323950700003917
(ado曲妥珠单抗)或
Figure BDA00036323950700003918
(帕妥珠单抗)治疗);PIK3CA改变(例如可以用
Figure BDA00036323950700003919
(阿培利司)治疗)。
针对结直肠癌适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括KRAS野生型(密码子12和13中不存在突变)(例如可以用
Figure BDA00036323950700003920
(西妥昔单抗)治疗);KRAS野生型(外显子2、3和4中不存在突变)和NRAS野生型(外显子2、3和4中不存在突变)(例如可以用
Figure BDA00036323950700003921
(帕尼单抗)治疗)。
针对卵巢癌适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括BRCA1/2改变(例如可以用
Figure BDA0003632395070000401
(奥拉帕尼)或
Figure BDA0003632395070000402
(卢卡帕尼)治疗)。
针对前列腺癌适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括同源重组修复(HRR)基因(BRCA1、BRCA2、ATM、BARD1、BRIP1、CDK12、CHEK1、CHEK2、FANCL、PALB2、RAD51B、RAD51C、RAD51D和RAD54L)改变(例如可以用
Figure BDA0003632395070000403
(奥拉帕尼)治疗)。
针对实体瘤癌症适应症,可以作为免疫疗法的生物标志物的变体或突变可以包括大于或等于每兆碱基10个突变的肿瘤突变负荷(TMB)(例如可以用
Figure BDA0003632395070000404
(派姆单抗)治疗)。
XI.计算机器架构
图15是图示了能够从计算机可读介质读取指令并且在处理器(或控制器)中执行它们的示例计算机器的组件的框图。本文描述的计算机可以包括图15所示的单个计算机器、虚拟机、包括图15所示的计算机器的多个节点的分布式计算系统或者计算设备的任何其他合适的布置。
通过示例,图15示出了计算机系统1500的示例形式的计算机器的图解表示,其中可以被存储在计算机可读介质中以使机器执行本文讨论的任何一个或多个过程的指令1524(例如软件、程序代码或机器代码)可以被执行。在一些实施例中,计算机器作为独立设备操作或者可以被耦合(例如联网)至其他机器。在联网部署中,机器可以在服务器-客户端网络环境中以服务器机器或者客户端机器的资格操作,或者作为对等(或者分布式)网络环境的对等机器操作。
图15中描述的计算机器的结构可以对应于任何软件、硬件或组合组件(例如图2所示的那些或本文描述的处理单元),包括但不限于被用于执行本文描述的一个或多个过程的任何引擎、模块、计算服务器、机器。尽管图15示出了各种硬件和软件元件,本文描述的组件中的每个组件可以包括附加的或更少的元件。
通过示例,计算机器可以是个人计算机(PC)、平板PC、机顶盒(STB)、个人数字助理(PDA)、蜂窝电话、智能手机、web设备、网络路由器、物联网(IoT)设备、交换机或桥接器或者能够执行指令1524的任何机器,该指令1524指定要由该机器采取的动作。进一步地,虽然仅单个机器被图示,但术语“机器”和“计算机”也可以被认为包括单独或联合执行指令1524以执行本文讨论的任何一种或多种方法的机器的任何集合。
示例计算机系统1500包括一个或多个处理器1502,诸如CPU(中央处理单元)、GPU(图形处理单元)、TPU(张量处理单元)、DSP(数字信号处理器)、片上系统(SOC)、控制器、状态设备、专用集成电路(ASIC)、现场可编程门阵列(FPGA)或这些的任何组合。计算系统1500的部分还可以包括存储器1504,该存储器1504存储包括指令1524的计算机代码,当指令被处理器1502直接或间接执行时,指令1524可以使处理器1502执行某些动作。指令可以是任何指导、命令或者可以以不同形式存储的顺序,诸如设备可读指令、包括源代码的编程指令以及其他通信信号和顺序。指令可以在一般意义上使用,并且不被限于机器可读代码。
本文描述的一种或多种方法提高了处理器1502的操作速度,并且减少了存储器1504所需的空间。例如,本文描述的机器学习方法通过应用一种或多种新颖技术降低了处理器1502的计算复杂度,这简化了处理器1502的训练、达到收敛并且生成结果的步骤。本文描述的算法还减小了模型和数据集的大小,以减少对存储器1504的存储空间要求。
某些操作的性能可以被分布在多于一个处理器之间,不仅驻留在单个机器内,而且部署在多个机器上。在一些示例实施例中,一个或多个处理器或处理器实施的模块可以位于单个地理位置(例如在家庭环境、办公室环境或服务器农场内)。在其他示例实施例中,一个或多个处理器或处理器实施的模块可以被分布在多个地理位置。尽管在说明书或权利要求书中可能引用由处理器执行的一些过程,但这应该被解释为包括多个分布式处理器的联合操作。
计算机系统1500可以包括主存储器1504和静态存储器1506,它们被配置为经由总线1508彼此通信。计算机系统1500还可以包括图形显示单元1510(例如等离子体显示面板(PDP)、液晶显示器(LCD)、投影仪或阴极射线管(CRT))。由处理器1502控制的图形显示单元1510显示图形用户界面(GUI),以显示由本文描述的过程生成的一个或多个结果和数据。计算机系统1500还可以包括字母数字输入设备1512(例如键盘)、光标控制设备1514(例如鼠标、轨迹球、操纵杆、运动传感器或其他指向仪器)、存储单元1516(硬盘驱动器、固态驱动器、混合驱动器、存储器磁盘等)、信号生成设备1518(例如扬声器)和网络接口设备1520,它们也被配置为经由总线1508进行通信。
存储单元1516包括计算机可读介质1522,在该计算机可读介质1522上存储有实施本文描述的任何一种或多种方法或功能的指令1524。指令1524还可以在由计算机系统1500执行期间完全或至少部分地驻留在主存储器1504内或处理器1502内(例如在处理器的缓存存储器内),主存储器1504和处理器1502也构成计算机可读介质。指令1524可以经由网络接口设备1520通过网络1526发送或接收。虽然计算机可读介质1522在示例实施例中被示出为单个介质,但是术语“计算机可读介质”应被理解为包括能够存储指令(例如指令1524)的单个介质或者多个介质(例如集中式或分布式数据库或者关联的缓存和服务器)。计算机可读介质可以包括能够存储指令(例如指令1524)以由处理器(例如处理器1502)执行并且使处理器执行本文公开的任何一种或多种方法的任何介质。计算机可读介质可以包括但不限于固态存储器、光学介质和磁性介质形式的数据储存库。计算机可读介质不包括诸如传播信号或载波等瞬态介质。
XII.附加考虑
有益地,本文描述的各种实施例提高了测序领域中现有技术的准确性和效率,诸如PCR和大规模并行DNA测序(例如NGS)。实施例提供了对标识由测序和扩增过程引入的误差的挑战的解决方案。大规模并行DNA测序可以从一个或多个DNA样本开始,这些样本被随机切割并且通常使用PCR扩增。大规模并行DNA测序的并行性质导致每个等位基因的核苷酸序列的复制。每个等位基因位点处的复制和测序程度可能会有所不同。例如,一些序列是重叠的和/或双链的,而其他序列则不是。PCR扩增过程和测序过程以及测序过程都具有非平凡的误差率。序列误差可能会掩盖真正等位基因的核苷酸序列。这些实施例可以被用于确定由大规模并行DNA测序仪器分析的一个或多个等位基因。通过考虑读数层特定噪声模型,大规模并行DNA测序工作流程表现出足够的保真度,以通过更准确地区分真实等位基因和错误序列来生成正确的序列确定。
通常,为了降低确定正确序列的误差率,需要增加样本的测序深度。这意味着可以在一批测序中分析更少的样本,因为更多的资源专用于样本。这些实施例在不增加特定等位基因位点的测序深度的情况下提高了测序的准确性,从而允许在大规模并行DNA测序的情况下同时对更多的等位基因位点或患者样本进行测序。所描述的实施例可以减少所需的测序深度,同时提高被用于读取扩增中生成的核苷酸序列的大规模并行DNA测序的准确性。
本发明的实施例的前述描述已经出于图示的目的而被呈现;它并不旨在穷举或将本发明限制为所公开的精确形式。相关领域的技术人员可以了解,鉴于以上公开内容,许多修改和变化是可能的。
本描述的一些部分根据算法和信息的操作的符号表示来描述本发明的实施例。这些算法描述和表示通常被数据处理领域的技术人员用于将他们工作的实质有效地传达给本领域的其他技术人员。虽然在功能上、计算上或逻辑上进行了描述,但这些操作被理解为由计算机程序或等效电路、微代码等来实施。此外,在不失一般性的情况下,有时将这些操作布置称为模块也被证明是方便的。所描述的操作及其关联模块可以在软件、固件、硬件或其任何组合中实施。
本文描述的任何步骤、操作或过程可以单独或与其他设备组合使用一个或多个硬件或软件模块来执行或实施。在一些实施例中,软件模块是用计算机程序产品实施的,该计算机程序产品包括包含计算机程序代码的非瞬态计算机可读介质,该计算机程序代码可以由计算机处理器执行以执行所描述的任何或所有步骤、操作或过程。
本发明的实施例还可以涉及由本文描述的计算过程产生的产品。这种产品可以包括由计算过程产生的信息,其中该信息被存储在非瞬态、有形的计算机可读存储介质上,并且可以包括计算机程序产品的任何实施例或本文描述的其他数据组合。
虽然本文描述的一个或多个过程可以用一个或多个步骤来描述,但术语“步骤”的使用并不意味着特定的顺序。例如,虽然本公开可以描述包括顺序的多个步骤的过程,但该过程中的步骤不需要按照本公开中要求保护或描述的具体顺序来执行。一些步骤可能会在其他步骤之前执行,尽管其他步骤在本公开中首先要求保护或描述。
最后,说明书中使用的语言主要是为了可读性和指导目的而选择的,并且可能不是为了描写或约束发明主题而选择的。因此,本发明的范围意在不受该详细描述的限制,而是受在基于此的申请上发布的任何权利要求限制。因此,本发明的实施例的公开内容旨在说明而非限制在以下权利要求中陈述的本发明的范围。

Claims (49)

1.一种用于处理样本的DNA测序数据集的计算机实现的方法,所述计算机实现的方法包括:
访问由DNA测序生成的所述DNA测序数据集,所述DNA测序数据集包括多个已处理的序列读数,所述多个已处理的序列读数包括变体位置;
将所述多个已处理的序列读数分层到多个读数层;
针对每个读数层,确定在所述变体位置处的分层测序深度;
针对每个读数层,确定以所述读数层的所述分层测序深度为条件的一个或多个噪声参数,所述一个或多个噪声参数对应特定于所述读数层的噪声模型,其中训练所述噪声模型包括:
对多个参考健康个体的训练DNA数据集分层,
为所述读数层选择分层序列读数作为分层训练集,
初始化所述一个或多个噪声参数,所述一个或多个噪声参数对表示所述噪声模型的噪声分布建模,以及
基于来自所述多个参考健康个体的所述分层训练集的所述噪声分布,迭代地调整所述一个或多个噪声参数的值;
针对每个读数层,基于以所述读数层的所述分层测序深度为条件的所述一个或多个噪声参数,生成特定于所述读数层的所述噪声模型的输出;以及
组合所生成的所述噪声模型输出以产生组合结果,所述组合结果代表所述样本与总变体计数相关联的可能性。
2.根据权利要求1所述的计算机实现的方法,其中所述多个读数层包括以下一项或多项:(1)双链的缝合读数层,(2)双链的未缝合读数层,(3)单链的缝合读数层,以及(4)单链的未缝合读数层。
3.根据权利要求1所述的计算机实现的方法,其中所述变体位置处的突变是以下中的一者:单核苷酸变体、插入和缺失。
4.根据权利要求1所述的计算机实现的方法,还包括:
确定所述组合结果的质量分数,所述质量分数是Phred等级分数。
5.根据权利要求4所述的计算机实现的方法,还包括:
响应于所述质量分数高于预定阈值,指示所述样本可能在所述变体位置处具有突变。
6.根据权利要求1所述的计算机实现的方法,其中针对读数层确定以所述读数层的所述分层测序深度为条件的所述一个或多个噪声参数包括:
访问特定于所述读数层的参数分布,所述参数分布描述了与所述读数层相关联的一组DNA测序样本的分布,其中所述噪声参数是从所述参数分布中而被确定的。
7.根据权利要求6所述的计算机实现的方法,其中针对每个读数层,与所述读数层相关联的所述一组DNA测序样本包括被分层到所述读数层的序列读数,并且对应于一个或多个健康个体。
8.根据权利要求6所述的计算机实现的方法,其中针对每个读数层,特定于所述读数层的所述噪声模型是贝叶斯层次模型,并且所述参数分布基于伽马分布。
9.根据权利要求1所述的计算机实现的方法,其中与特定于第一读数层的噪声模型相对应的第一噪声参数具有与特定于第二读数层的噪声模型相对应的对应第二噪声参数不同的值。
10.根据权利要求1所述的计算机实现的方法,其中针对每个读数层,所确定的所述一个或多个噪声参数包括以所述读数层的所述分层测序深度为条件的所述噪声分布的平均值。
11.根据权利要求10所述的计算机实现的方法,其中每个噪声分布是以每个读数层的所述分层测序深度为条件的负二项分布。
12.根据权利要求11所述的计算机实现的方法,其中针对每个读数层,所确定的所述一个或多个噪声参数还包括分散参数。
13.根据权利要求1所述的计算机实现的方法,其中所生成的每个噪声模型的所述输出是针对所述读数层所确定的以所述分层测序深度为条件的所述一个或多个噪声参数。
14.根据权利要求1所述的计算机实现的方法,其中所生成的每个噪声模型的所述输出包括所述读数层的分层变体计数超过阈值的可能性。
15.根据权利要求1所述的计算机实现的方法,其中组合所生成的所述噪声模型输出包括组合来自每个噪声模型输出的平均变体计数和方差,以产生代表所述组合结果的总体噪声分布的总体平均变体计数和所述总体分散参数。
16.根据权利要求15所述的计算机实现的方法,其中所述总体噪声分布是基于负二项分布而被建模的,并且其中确定所述总体平均变体计数和所述总体分散参数包括:
基于所述读数层的所述分层测序深度,确定针对每个读数层的所述平均变体计数;
确定针对每个读数层的所述方差;
针对每个读数层的所述平均变体计数求和以确定所述总平均变体计数;
组合针对每个读数层的所述方差以确定总体方差;以及
基于所述总体平均变体计数和所述总体方差,确定所述总体分散参数。
17.根据权利要求1所述的计算机实现的方法,其中组合所生成的所述噪声模型输出以产生所述组合结果包括:
确定每个读数层的观察到的分层变体计数;
在每个读数层中确定比每个读数层的观察到的所述分层变体计数更有可能的可能事件;
标识与比每个读数层的观察到的所述分层变体计数更高的发生可能性相关联的所述可能事件的组合;
对所标识的所述组合的概率求和以确定统计补数;以及
通过从1.0中减去所述统计补数来确定可能性值。
18.根据权利要求17所述的计算机实现的方法,其中包括一个双链读数的第一标识组合等同于包括两个单链读数的第二标识组合。
19.根据权利要求17所述的计算机实现的方法,其中所确定的所述可能性值等于或大于每个读数层的观察到的所述分层变体计数的发生可能性。
20.根据权利要求17所述的计算机实现的方法,还包括训练机器学习模型以确定所述可能性值。
21.根据权利要求1所述的计算机实现的方法,还包括:
接收个体的体液样本;
对所述体液样本的cfDNA执行所述DNA测序;
基于所述DNA测序的结果,生成原始序列读数;以及
折叠并且缝合所述原始序列读数以生成所述多个已处理的序列读数。
22.根据权利要求21所述的计算机实现的方法,其中所述体液样本是以下之一的样本:所述个体的血液、全血、血浆、血清、尿液、脑脊液、粪便、唾液、眼泪、组织活检、胸膜液、心包液或腹膜液。
23.根据权利要求21所述的计算机实现的方法,其中所述多个已处理的序列读数是从肿瘤活检测序的。
24.根据权利要求21所述的计算机实现的方法,其中所述多个已处理的序列读数是从来自血液的细胞分离物测序的,所述细胞分离物至少包括血沉棕黄层白细胞或CD4+细胞。
25.根据权利要求1所述的计算机实现的方法,其中所述DNA测序包括大规模并行DNA测序操作。
26.根据权利要求1所述的计算机实现的方法,其中所述DNA测序数据集是个体的体液样本的cfDNA测序数据集。
27.根据权利要求1所述的计算机实现的方法,还包括:
基于所述组合结果,提供对具有变体的受试者的诊断。
28.根据权利要求27所述的计算机实现的方法,其中所述变体选自由以下组成的所述组: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、EPHA5、EPHB1、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、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和TNFRSF14。
29.根据权利要求27所述的计算机实现的方法,还包括:
向被标识为具有所述变体的所述受试者提供施用治疗的指导。
30.根据权利要求29所述的计算机实现的方法,其中所述治疗包括施用选自由以下组成的所述组的药物:利妥昔单抗、曲妥珠单抗、西妥昔单抗、帕尼单抗、奥法木单抗、贝利尤单抗、伊匹单抗、帕妥珠单抗、曲美木单抗、纳武单抗、达西组单抗、乌瑞芦单抗、阿特珠单抗、派姆单抗、博纳吐单抗、CT-011、帕博利珠单抗、BMS-936559、MED14736、MSB0010718C、度伐鲁单抗、阿维鲁单抗和玛格妥昔单抗。
31.根据权利要求1所述的计算机实现的方法,其中所述可能性表示后续观察到的数据的总变体计数大于或等于所述多个已处理序列读数中观察到的总变体计数可归因于噪声。
32.一种非瞬态计算机可读介质,包括指令,所述指令在由一个或多个处理器执行时,使所述一个或多个处理器执行步骤,包括:
访问由DNA测序生成的所述DNA测序数据集,所述DNA测序数据集包括多个已处理的序列读数,所述多个已处理的序列读数包括变体位置;
将所述多个已处理的序列读数分层到多个读数层;
针对每个读数层,确定在所述变体位置处的分层测序深度;
针对每个读数层,确定以所述读数层的所述分层测序深度为条件的一个或多个噪声参数,所述一个或多个噪声参数对应特定于所述读数层的噪声模型,其中所述噪声模型的训练包括:
对多个参考健康个体的训练DNA数据集分层,
为所述读数层选择分层序列读数作为分层训练集,
初始化所述一个或多个噪声参数,所述一个或多个噪声参数对表示所述噪声模型的噪声分布建模,以及
基于来自所述多个参考健康个体的所述分层训练集的所述噪声分布,迭代地调整所述一个或多个噪声参数的值;
针对每个读数层,基于以所述读数层的所述分层测序深度为条件的所述一个或多个噪声参数,生成特定于所述读数层的所述噪声模型的输出;以及
组合所生成的所述噪声模型输出以产生组合结果,所述组合结果代表后续观察到的数据的总变体计数大于或等于所述多个已处理序列读数中观察到的总变体计数可归因于噪声的可能性。
33.根据权利要求32所述的非瞬态计算机可读介质,其中组合所生成的所述噪声模型输出包括组合来自每个噪声模型输出的平均变体计数和方差,以产生代表所述组合结果的总体噪声分布的总体平均变体计数和所述总体分散参数。
34.根据权利要求33所述的非瞬态计算机可读介质,其中所述总体噪声分布是基于负二项分布而被建模的,并且其中确定所述总体平均变体计数和所述总体分散参数包括:
基于所述读数层的所述分层测序深度,确定每个读数层的所述平均变体计数;
确定每个读数层的所述方差;
针对每个读数层的所述平均变体计数求和以确定所述总体平均变体计数;
组合每个读数层的所述方差以确定总体方差;以及
基于所述总体平均变体计数和所述总体方差,确定所述总体分散参数。
35.根据权利要求32所述的非瞬态计算机可读介质,其中组合所生成的所述噪声模型输出以产生所述组合结果包括:
确定每个读数层的观察到的分层变体计数;
在每个读数层中确定比每个读数层的观察到的所述分层变体计数更有可能的可能事件;
标识与比每个读数层的观察到的所述分层变体计数更高的发生可能性相关联的所述可能事件的组合;
对所标识的所述组合的概率求和以确定统计补数;以及
通过从1.0中减去所述统计补数来确定可能性值。
36.根据权利要求32所述的非瞬态计算机可读介质,其中所述步骤还包括:
基于所述组合结果,提供具有变体的受试者的诊断。
37.根据权利要求36所述的非瞬态计算机可读介质,其中所述变体选自由以下组成的所述组:ACVR1B、AKT3、AMER1、APC、ARID1A、ARID1B、ARID2、ASXL1、ASXL2、ATM、ATR、BAP1BCL2、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、EPHA5、EPHB1、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、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和TNFRSF14。
38.根据权利要求36所述的非瞬态计算机可读介质,其中所述步骤还包括:
向被标识为具有所述变体的所述受试者提供施用治疗的指导。
39.根据权利要求38所述的非瞬态计算机可读介质,其中所述治疗包括施用选自由以下组成的所述组的药物:利妥昔单抗、曲妥珠单抗、西妥昔单抗、帕尼单抗、奥法木单抗、贝利尤单抗、伊匹单抗、帕妥珠单抗、曲美木单抗、纳武单抗、达西组单抗、乌瑞芦单抗、阿特珠单抗、派姆单抗、博纳吐单抗、CT-011、帕博利珠单抗、BMS-936559、MED14736、MSB0010718C、度伐鲁单抗、阿维鲁单抗和玛格妥昔单抗。
40.根据权利要求32所述的非瞬态计算机可读介质,其中所述可能性表示后续观察到的数据的总变体计数大于或等于所述多个已处理序列读数中观察到的总变体计数可归因于噪声。
41.一种系统,包括计算机处理器和存储器,所述存储器存储计算机程序指令,所述计算机程序指令在由所述计算机处理器执行时使所述处理器执行包括以下步骤的步骤:
访问由DNA测序生成的所述DNA测序数据集,所述DNA测序数据集包括多个已处理的序列读数,所述多个已处理的序列读数包括变体位置;
将所述多个已处理的序列读数分层到多个读数层;
针对每个读数层,确定在所述变体位置处的分层测序深度;
针对每个读数层,确定以所述读数层的所述分层测序深度为条件的一个或多个噪声参数,所述一个或多个噪声参数对应特定于所述读数层的噪声模型,其中训练所述噪声模型包括:
对多个参考健康个体的训练DNA数据集分层,
为所述读数层选择分层序列读数作为分层训练集,
初始化所述一个或多个噪声参数,所述一个或多个噪声参数对表示所述噪声模型的噪声分布建模,以及
基于来自所述多个参考健康个体的所述分层训练集的所述噪声分布,迭代地调整所述一个或多个噪声参数的值;
针对每个读数层,基于以所述读数层的所述分层测序深度为条件的所述一个或多个噪声参数,生成特定于所述读数层的所述噪声模型的输出;以及
组合所生成的所述噪声模型输出以产生组合结果,所述组合结果代表后续观察到的数据的总变体计数大于或等于所述多个已处理序列读数中观察到的总变体计数可归因于噪声的可能性。
42.根据权利要求41所述的系统,其中组合所生成的所述噪声模型输出包括组合来自每个噪声模型输出的平均变体计数和方差,以产生代表所述组合结果的总体噪声分布的总体平均变体计数和所述总体分散参数。
43.根据权利要求42所述的系统,其中所述总体噪声分布是基于负二项分布而被建模,并且其中确定所述总体平均变体计数和所述总体分散参数包括:
基于所述读数层的所述分层测序深度,确定针对每个读数层的所述平均变体计数;
确定针对每个读数层的所述方差;
针对每个读数层的所述平均变体计数求和以确定所述总体平均变体计数;
组合每个读数层的所述方差以确定总体方差;以及
基于所述总体平均变体计数和所述总体方差,确定所述总体分散参数。
44.根据权利要求41所述的系统,其中组合所生成的所述噪声模型输出以产生所述组合结果包括:
确定每个读数层的观察到的分层变体计数;
在每个读数层中确定比每个读数层的观察到的所述分层变体计数更有可能的可能事件;
标识与比每个读数层的观察到的所述分层变体计数更高的发生可能性相关联的所述可能事件的组合;
对所标识的所述组合的概率求和以确定统计补数;以及
通过从1.0中减去所述统计补数来确定可能性值。
45.根据权利要求41所述的系统,其中所述步骤还包括:
基于所述组合结果,提供对具有变体的受试者的诊断。
46.根据权利要求45所述的系统,其中所述变体选自由以下组成的所述组: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、EPHA5、EPHB1、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、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和TNFRSF14。
47.根据权利要求45所述的系统,其中所述步骤还包括:
向被标识为具有所述变体的所述受试者提供施用治疗的指导。
48.根据权利要求47所述的系统,其中所述治疗包括施用选自由以下组成的所述组的药物:利妥昔单抗、曲妥珠单抗、西妥昔单抗、帕尼单抗、奥法木单抗、贝利尤单抗、伊匹单抗、帕妥珠单抗、曲美木单抗、纳武单抗、达西组单抗、乌瑞芦单抗、阿特珠单抗、派姆单抗、博纳吐单抗、CT-011、帕博利珠单抗、BMS-936559、MED14736、MSB0010718C、度伐鲁单抗、阿维鲁单抗和玛格妥昔单抗。
49.根据权利要求41所述的系统,其中所述可能性表示后续观察到的数据的总变体计数大于或等于所述多个已处理序列读数中观察到的总变体计数可归因于噪声。
CN202080077722.0A 2019-09-09 2020-09-08 用于分析dna数据的读数层特定噪声模型 Pending CN114746947A (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
US201962897923P 2019-09-09 2019-09-09
US62/897,923 2019-09-09
PCT/US2020/049751 WO2021050439A1 (en) 2019-09-09 2020-09-08 Read-tier specific noise models for analyzing dna data

Publications (1)

Publication Number Publication Date
CN114746947A true CN114746947A (zh) 2022-07-12

Family

ID=72615952

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202080077722.0A Pending CN114746947A (zh) 2019-09-09 2020-09-08 用于分析dna数据的读数层特定噪声模型

Country Status (7)

Country Link
US (1) US20220336044A1 (zh)
EP (1) EP4026130A1 (zh)
JP (1) JP2022546649A (zh)
CN (1) CN114746947A (zh)
CA (1) CA3150532A1 (zh)
IL (1) IL291145A (zh)
WO (1) WO2021050439A1 (zh)

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2019046804A1 (en) * 2017-09-01 2019-03-07 Grail, Inc. IDENTIFICATION OF FALSE POSITIVE VARIANTS USING A MODEL OF IMPORTANCE
CA3092352A1 (en) * 2018-02-27 2019-09-06 Cornell University Systems and methods for detection of residual disease

Also Published As

Publication number Publication date
EP4026130A1 (en) 2022-07-13
CA3150532A1 (en) 2021-03-18
WO2021050439A1 (en) 2021-03-18
US20220336044A1 (en) 2022-10-20
JP2022546649A (ja) 2022-11-04
IL291145A (en) 2022-05-01

Similar Documents

Publication Publication Date Title
Sammut et al. Multi-omic machine learning predictor of breast cancer therapy response
Ding et al. Expanding the computational toolbox for mining cancer genomes
JP7497084B2 (ja) がん治療の有効性を予測するためのシステムおよび方法
JP6987786B2 (ja) がんの進化の検出および診断
TWI814753B (zh) 用於標靶定序之模型
CN110305965A (zh) 一种预测非小细胞肺癌(nsclc)患者对免疫疗法的敏感性的方法
CN110958853A (zh) 用于鉴定或监测肺病的方法和系统
JP2003021630A (ja) 臨床診断サービスを提供するための方法
CN113096728B (zh) 一种微小残余病灶的检测方法、装置、存储介质及设备
Dlamini et al. AI and precision oncology in clinical cancer genomics: From prevention to targeted cancer therapies-an outcomes based patient care
McCabe et al. Investigating the suitability of in vitro cell lines as models for the major subtypes of epithelial ovarian cancer
TWI781230B (zh) 使用針對標靶定序的定點雜訊模型之方法、系統及電腦產品
JP6700376B2 (ja) 有意性が未知のバリアントに優先順位をつけるシステム及び方法
CN114746947A (zh) 用于分析dna数据的读数层特定噪声模型
US20220301654A1 (en) Systems and methods for predicting and monitoring treatment response from cell-free nucleic acids
WO2018209704A1 (zh) 基于dna测序数据的样本来源检测方法、装置和存储介质
Tao et al. Improving personalized prediction of cancer prognoses with clonal evolution models
JP2024512540A (ja) 人工知能基盤の無細胞DNAの腫瘍由来変異の検出方法及びこれを用いたがんの早期診断方法{Method for detecting tumor derived mutation from cell-free DNA based on artificial intelligence and Method for early diagnosis of cancer using the same}
US20200105374A1 (en) Mixture model for targeted sequencing
Subramanian et al. Novel multisample scheme for inferring phylogenetic markers from whole genome tumor profiles
Madjar Survival models with selection of genomic covariates in heterogeneous cancer studies
Dlamini et al. Informatics in Medicine Unlocked
Weber Applications of ctDNA Genomic Profiling to Metastatic Triple Negative Breast Cancer
Wang Statistical Models for Alternative Splicing with Applications to Heterogeneous Disease
Ogundijo Bayesian Inference for Genomic Data Analysis

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination