CN108603229A - 用于高保真测序的方法和系统 - Google Patents

用于高保真测序的方法和系统 Download PDF

Info

Publication number
CN108603229A
CN108603229A CN201780007584.7A CN201780007584A CN108603229A CN 108603229 A CN108603229 A CN 108603229A CN 201780007584 A CN201780007584 A CN 201780007584A CN 108603229 A CN108603229 A CN 108603229A
Authority
CN
China
Prior art keywords
sequencing
nucleic acid
sample
variant
assemblage
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
CN201780007584.7A
Other languages
English (en)
Inventor
奥利弗·克劳德·维恩
亚历山大·蒂洛·狄尔泰
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
Grier Co
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 Grier Co filed Critical Grier Co
Publication of CN108603229A publication Critical patent/CN108603229A/zh
Pending legal-status Critical Current

Links

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
    • G16B30/00ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • 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/6813Hybridisation assays
    • C12Q1/6827Hybridisation assays for detection of mutation or polymorphism
    • 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
    • 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
    • G16B20/20Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
    • 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/10Sequence alignment; Homology search
    • 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
    • G16B5/00ICT specially adapted for modelling or simulations in systems biology, e.g. gene-regulatory networks, protein interaction networks or metabolic networks
    • 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
    • G16B5/00ICT specially adapted for modelling or simulations in systems biology, e.g. gene-regulatory networks, protein interaction networks or metabolic networks
    • G16B5/10Boolean models
    • 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
    • G16B5/00ICT specially adapted for modelling or simulations in systems biology, e.g. gene-regulatory networks, protein interaction networks or metabolic networks
    • G16B5/20Probabilistic models
    • 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
    • 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
    • C12Q2535/00Reactions characterised by the assay type for determining the identity of a nucleotide base or a sequence of oligonucleotides
    • C12Q2535/122Massive parallel sequencing

Landscapes

  • Life Sciences & Earth Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Chemical & Material Sciences (AREA)
  • Biophysics (AREA)
  • Biotechnology (AREA)
  • General Health & Medical Sciences (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Evolutionary Biology (AREA)
  • Theoretical Computer Science (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Medical Informatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Molecular Biology (AREA)
  • Analytical Chemistry (AREA)
  • Organic Chemistry (AREA)
  • Genetics & Genomics (AREA)
  • Physiology (AREA)
  • Zoology (AREA)
  • Wood Science & Technology (AREA)
  • Immunology (AREA)
  • General Engineering & Computer Science (AREA)
  • Microbiology (AREA)
  • Biochemistry (AREA)
  • Probability & Statistics with Applications (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本文描述了用于高保真测序和鉴定样品中稀浓度的罕见突变的系统和方法。在多个方面,使用包括衔接子连接条件和杂合体捕获富集组的专门的文库制备技术与对照一起使用,以增加序列即用分子的产率并鉴定污染和错误且使污染和错误最小化。系统和方法还涉及使用系综和拟极大似然模型分析测序数据以区分真实变体和假阳性。

Description

用于高保真测序的方法和系统
相关申请的交叉引用
本申请要求于2016年1月22日提交的美国临时专利申请序列号62/286,110的优先权权益,该申请的公开内容通过引用以其整体并入本文。
发明领域
本发明涉及通过测定优化和数据分析用于样品中的稀变体(dilute variants)的高保真测序和鉴定的系统和方法。
背景
仅在美国每年诊断出数百万例癌症,并有数十万人死亡。包括癌症在内的许多疾病的根源是个体DNA中的遗传突变或变体。在癌症的情况下,这些突变可导致异常的细胞生长,这可能是无法控制的并导致死亡。这些疾病和潜在的突变的早期检测可能对成功治疗这些疾病至关重要。最近的进展允许从血液和其他体液分离肿瘤来源的核酸和其他无细胞核酸。这些发展允许患者突变的更便宜、非侵入性的检查和表征。不幸的是,感兴趣的突变,特别是在疾病发展的早期阶段,通常以低于标准测序错误率的频率发生。大多数无细胞核酸包括个体的正常基因组序列,而少得多的量为肿瘤起源的,呈现肿瘤特异性突变。现代测序和分析技术以在每1,000个位置读段中1个错误或99.9%的错误率进行,并且通常不足以在无细胞样品诸如血液或血浆中检测罕见肿瘤变体。问题在于,在测序过程中区分实际突变和通过错误引入的假阳性变得几乎不可能。因此,以低频率出现的疾病特异性突变的早期鉴定的诺言未被实现,并且丧失了早期干预的益处。
概述
本发明涉及用于罕见核酸变体的高保真测序和鉴定的方法和系统。本发明的系统和方法可用于鉴定无细胞核酸样品中的罕见变体,例如包含正常基因组核酸多数的样品中的肿瘤特异性突变。本发明的系统和方法允许样品中以低于1:10,000的频率发生的突变的可靠鉴定。这种罕见变体的鉴定来自测序过程中几个步骤的优化,然后是基于本文称为系综(ensembles)的对齐读段对的测序读段的分析。
本发明的系统和方法可以在罕见变体鉴定之外应用,诸如性能或灵敏度的期望水平的测序优化。通过使用本发明来对特定的应用定制测序程序,从业人员可以通过仅要求特定应用必需的确切数量的测序读段来避免另外的成本和时间。
本发明的各方面包括用于核酸测序的方法。方法的步骤可包括:获得核酸的测序读段,鉴定包含具有共享起始坐标和读段长度的两个或更多个测序读段的系综,确定系综所包括的测序分子的数量,鉴定系综中的候选变体,并且使用似然估计模型和所确定的测序分子的数量来确定候选变体是真实变体的可能性。在某些实施方案中,获得测序读段的步骤可进一步包括从核酸制备测序文库,扩增测试文库,并使用下一代测序(NGS)将测序文库测序。在某些实施方案中,衔接子可以在被配置为允许衔接子堆叠的条件下连接至核酸。测序文库的制备可包括使用约16小时的反应时间、在约16摄氏度的温度将衔接子连接至核酸。扩增步骤可以包括PCR扩增,并且本发明的方法可以进一步包括使用计算机模型选择所需的过度扩增因子和PCR循环数目以检测在样品中的指定浓度的变体。
在各种实施方案中,本发明的方法包括基于包括鸟嘌呤-胞嘧啶(GC)含量、靶群体中的突变频率和序列独特性的因素设计靶向基因组区域的杂合体捕获组,并在测序步骤之前使用杂合体捕获组捕获扩增的核酸。捕获步骤可以包括使用靶向靶基因座的有义链的第一杂合体捕获组和靶向靶基因座的反义链的第二杂合体捕获组。
在某些实施方案中,可以在扩增测序文库之前将合成的核酸对照(也称为对照序列、对照添加物(spike-in)或阳性对照)添加至核酸,并且错误率然后可以使用合成的核酸对照的测序读段确定。合成的核酸对照可以包含已知序列,所述已知序列在核酸所来源的物种中具有低多样性,并且具有与已知序列多于一个非天然存在的错配,并且在某些实施方案中,多于一个非天然存在的错配可以是4个。合成的核酸对照可包括代表杂合体捕获组的靶基因座的鸟嘌呤-胞嘧啶(GC)含量分布,或可包括包含与杂合体捕获组的下拉(pulldown)探针的不同重叠的多于一种核酸。错误率或候选变体频率可以使用合成的核酸对照的测序读段来确定。
在各种实施方案中,核酸可以包含无细胞核酸或可以从组织样品获得,其中获得测序读段还包括在制备步骤之前将核酸片段化。片段化可利用超声处理或酶促裂解生成。
本发明的方法可包括如果在核酸的有义链和反义链两者上均未鉴定出候选变体,则丢弃候选变体。
在某些方面,本发明包括用于鉴定核酸变体的系统。系统包括耦合到存储指令的有形的、非瞬时存储器的处理器,所述指令在由处理器执行时使系统执行各种步骤。本发明的系统可以是可操作的,以:鉴定包含具有共享起始坐标和读段长度的两个或更多个测序读段的系综,确定系综所包括的测序分子的数量,鉴定系综中的候选变体,并且使用似然估计模型和所确定的测序分子的数量来确定候选变体是真实变体的可能性。
在某些实施方案中,本发明的系统可以是可操作的,以:如果在核酸的有义链和反义链两者上均未鉴定出候选变体,则丢弃候选变体。本发明的系统可以是进一步可操作的,以基于包括鸟嘌呤-胞嘧啶(GC)含量、靶群体中的突变频率和序列独特性的因素对两种或更多种测序读段确定靶基因组区域。
附图简述
图1提供本发明方法的图解。
图2说明包括堆叠的衔接子的测序相容衔接子连接产物。
图3说明带有堆叠的衔接子的连接产物的PCR结果。
图4说明制备的无细胞DNA文库的分子长度的分布。
图5说明使用衔接子特异性引物的PCR扩增后无细胞DNA文库的分子长度的分布。
图6提供杂合体捕获组设计过程的图。
图7说明使用合成的DNA对照来鉴定无细胞DNA样品的污染。
图8说明本发明的计算机系统。
详细说明
本发明的系统和方法一般涉及使用优化的测序技术和测序读段分析的高保真测序和鉴定罕见核酸变体的方法。分子群体中低丰度突变(突变等位基因比例<5000-1)的检测和准确频率估计的必要条件是,在整个样品制备和文库制备过程中保持衍生的等位基因Nd(对应于体细胞变体)与祖先等位基因Na(对应于种系基因组)和来自其他来源的DNA N的比例。
衍生的等位基因的比例f可以由于在测序文库构建过程中经由损失消耗Nd、或通过污染增加分母而减小。因此,为了在样品、包括细胞中鉴定以低浓度水平在无细胞DNA中存在的突变或变体,人们必须在文库制备期间使污染最小化并使分子的损失最小化。本申请提出了用于实现这些目标的系统和方法,以及用于区分真实变体和假阳性的测序分析技术。通过优化文库制备和测序步骤,减少测序错误和包括变体验证步骤,本发明的系统和方法允许鉴定核酸样品中以1:10,000或更低的比存在的变体。罕见变体的鉴定有许多应用,包括鉴定主要由患者正常基因组DNA构成的无细胞DNA中的肿瘤、癌症、或疾病特异性突变。本发明的系统和方法利用与下一代NGS测序机的错误率相比高保真PCR酶的较低错误率,通过经由样品的PCR扩增增加待测序的分子数量,并与测序后分析组合来证实候选变体的有效性,来增加鉴定序列变体的灵敏度。
根据本发明的某些方面的系统和方法在图1中示出。步骤可包括测序文库制备101、测序文库扩增103和文库的测序105。本发明的系统和方法可以通过首先获得测序读段107来实现或者可以以核酸样品和上述步骤开始以产生测序读段。接下来,在测序读段中鉴定系综109,并且确定样品中在每个系综之下的原始分子的数量111。使用上述信息和参考序列,鉴定候选变体113,并且使用概率模型来确定候选变体是真实变体的可能性115。
样品制备
在某些实施方案中,核酸可以从患者样品获得。患者样品可以例如包括血液样品、全血、血浆、眼泪、乳头吸出物、血清、粪便、尿、唾液、循环细胞、组织、活组织检查样品或含有患者生物材料的其他样品。在优选的实施方案中,核酸是从患者血液或血浆分离的。血液样品在采集后迅速处理,以使由凋亡的有核细胞的DNA释放造成的污染最小化。
下面是用于从血液制备核酸的示例性程序。可以将血液收集在10ml EDTA管中(例如,可从Becton Dickinson获得)。Streck cfDNA管(Streck,Inc.,Omaha,Nebraska)可用于通过化学固定有核细胞使污染最小化,但是如优选实施方案中,当在2小时或更短时间内处理样品时,观察到来自基因组DNA的很少的污染。从血液样品开始,血浆可以通过在室温以3000rpm离心10分钟(减去制动)来提取。然后可将血浆以1ml等分试样转移至1.5ml试管,并在室温以7000rpm再次离心10分钟。然后可将上清液转移到新的1.5ml试管中。在这个阶段,样品可以保存在-80℃。在某些实施方案中,样品可以在血浆阶段储存用于后续处理,因为血浆可以比储存的提取无细胞DNA(cfDNA)更稳定。
然后使用商购可得的测定例如Qiagen QIAmp循环核酸试剂盒(Qiagen N.V.,Venlo Netherlands)从血液样品(例如血浆样品)提取核酸(例如DNA)。在某些实施方案中,可以使用以下修改的洗脱策略。可以使用Qiagen QIAmp循环核酸试剂盒按照制造商的说明提取DNA(每柱允许的最大血浆量为5ml)。如果从收集在Streck管中的血液的血浆提取cfDNA,则与蛋白酶K的反应时间可以从30分钟加倍到60分钟。优选地,应使用尽可能大的体积(即5mL)。在各种实施方案中,可以使用两步洗脱来使cfDNA产率最大化。首先,可以使用每个柱30μl缓冲液AVE洗脱DNA。为了增加cfDNA浓度,洗脱中可以使用完全覆盖膜必需的最少量缓冲液。通过用少量缓冲液减少稀释,可以避免样品的下游干燥,以防止双链DNA解链或物质损失。
随后,可以洗脱每个柱约30μl的缓冲液。在优选的实施方案中,可以使用第二洗脱来增加DNA产率。表1显示了使用上述方法中的第一和第二洗脱(其中两次洗脱体积为约30μl),来自六名黑素瘤患者的cfDNA样品观察到的DNA的量。另外洗脱的有用性可以通过平衡所获得的另外DNA与降低洗脱中最终DNA浓度来确定。然后可以使用商购可得的测定诸如Qubit DNA高灵敏度试剂盒(Thermo Fisher Scientific,Inc.,Cambridge,MA)将洗脱物组合并优选一式三份地定量DNA。
表1:洗脱物中的DNA浓度
样品ID 血浆体积(mL) 洗脱物1(ng) 洗脱物2(ng)
血浆009 3 12.63 5.22
血浆010 3 11.76 6.12
血浆045 3 21 4.14
血浆020 3 20.94 5.7
血浆062 3 17.1 5.88
血浆063 3 18.9 6.6
文库制备
尽管可以使用标准文库制备来生成文库,本发明的高产率方案相对于通常具有约40%的文库转换产率的标准方法提高了性能。本发明的方法提供约80%的文库转换。根据本发明,可以从核酸样品制备测序文库。商购可得的试剂盒可用于制备测序文库,诸如用于全基因组测序(WGS)的Illumina's TruSeq Nano试剂盒(Illumina,Inc.,San Diego,California)。可以修改试剂化学计量和孵育时间以增加通过该过程具有正确的测序衔接子连接的分子的数量(文库转换效率)。如果样品靶是样品中的cfDNA,则不需要片段化。在某些实施方案中,核酸可以从组织样品诸如肿瘤活组织检查获得。在这种情况下,核酸应该使用本领域已知的方法诸如超声处理或酶限制进行片段化。在实践中,未片段化的cfDNA群体的平均长度可以是约150-180个碱基并且因个体而异。在优选的实施方案中不使用固相可逆固定(SPRI)珠净化步骤,而是直接将样品进行末端修复以使cfDNA的损失最小化。这消除了携带乙醇进入PCR的风险;乙醇是PCR的抑制物,且在SPRI珠开始破裂之前去除所有乙醇液滴是挑战性的。避免SPRI净化步骤另外地减少了操作时间和成本。基于样品中DNA片段的估计数量,可以通过因子A调整试剂体积,以说明来自TruSeq Nano方案中规定的相对于超声处理的基因组DNA Ng的片段的不同数量的cfDNA片段Nf。这种调整可以应用于末端修复、3’末端腺苷酸化和衔接子连接步骤中使用的试剂。群体i中分子的数量Ni可以通过将群体的质量mi除以一个双脱氧核糖核苷酸的平均分子量(w=6.5E+11ng/摩尔)与每个分子中的平均碱基数Li的乘积,然后将此值乘以Avogadro常数来计算,如下所示:
然后调整因子A是Nf除以Ng的商:
在某些实施方案中,mg=100ng的输入DNA,并且指定的超声处理使片段长度Lg=350个碱基且对于给定样品,mf和Lf可以使用上述方程通过实验确定。然后可以使用已知的末端修复技术进一步处理核酸样品,以确保每个分子没有突出端,并含有5'磷酸和3'羟基,然后是3'腺苷酸化和衔接子连接。
在各种实施方案中,修改的衔接子连接程序可用于增加衔接子连接的cfDNA片段的产率。为了使具有连接的至少两个Y形Illumina测序衔接子(使用Illumina测序时)的cfDNA片段的数量最大化,衔接子连接反应时间可以增加到16小时和/或溶液中分子的动能可以使用16℃的较低孵育温度降低。在某些实施方案中,衔接子连接可以在鼓励衔接子连接并可导致衔接子的“堆叠”的条件下进行,诸如刚刚描述的那些,如图2所示。(203)。堆叠的衔接子在PCR扩增后分解(resolve),使得原始分子后代PCR产物不被阻止测序。图3说明了PCR过程期间堆叠的衔接子的分解(resolution)。空间位阻导致在扩增的PCR循环中最内侧的引物被选择。当最内侧引物在最外侧引物之前或与最外侧引物同时结合时,最外侧引物位点将在所得PCR产物中被消除。最内侧引物在最外侧引物之前退火的时间是几何分布的,成功概率为约.5,使得在4轮PCR扩增后,获得测序相容产物的概率为约15:16。
图4说明了来自肺癌患者的cfDNA文库的片段长度,其中平均分子长度为174个碱基,且每个衔接子为60个碱基。图5说明使用衔接子特异性引物PCR扩增后制备的文库。这些图说明,发生了衔接子堆叠,和通过PCR扩增有效分解了堆叠的衔接子,导致与配对末端测序(paired-end sequencing)相容的分子的更高产率。图4中的前三个峰对应于平均分子长度加上2、3和4个衔接子。
然后可以使用SPRI样品纯化珠以1:1.6、然后是1:1的样品:珠的比来净化扩增的样品,以去除游离的衔接子。然后可将样品洗脱至约27.5μl的体积。
根据某些实施方案,然后可以使用例如Bioanalyzer(Agilent Technologies,Santa Clara,California)或等效仪器来确定样品片段长度。可以输入约1μl的cfDNA来鉴定文库制备前后的平均片段长度。测序文库制备前cfDNA分子长度的分布可近似为正态分布取样,Xpre~N(μpre2),平均长度μ0约=150-180个碱基,和样品方差σ2。文库制备后分子长度的分布Xpost是以连接的测序衔接子的数量偏移的正态分布的叠加,每个测序衔接子具有固定的长度A,对于上述的Illumina平台通常为60个碱基(P5和P7衔接子)。可被测序(可测序)的分子具有连接到cfDNA片段的每一端的至少1个连接子,因此具有平均值μ0+kA,其中k≥2。如果文库被PCR扩增,如果连接的衔接子的数目k为至少2,则可以生成可测序的分子:
其中Yk是连接有k个衔接子的分子的贡献权重。在使用P5和P7PCR引物进行PCR扩增后,总群体应该以μpre+2A为主(如图3和4中所示)。
文库的质量可以使用Kapa文库定量试剂盒(Kapa Biosystems,Inc.Wilmington,Massachusetts)定量。文库可以使用任何已知的扩增方法(包括PCR扩增)扩增。为了进一步降低错误率,在优选的实施方案中,文库扩增可以使用Kapa HiFi Hotstart扩增(KapaBiosystems,Inc.Wilmington,Massachusetts KR0370-v5.13)进行。跨GC含量具有稳健性能的高保真PCR酶,诸如Kapa HiFi Hotstart,具有比Taq聚合酶的错误率低多达100倍的错误率。重复读段的水平可能影响所需的测序总量。模拟机可用于评估最佳过度扩增因子,以检测指定频率的变体、文库制备期间的共同地掺入损失、诱导误差和调用算法依赖性。在适用的情况下,模拟可考虑PCR扩增和杂合体捕获或其他下拉或富集技术中的损失。
系综中读段与潜在原始分子的比可以称为过度扩增因子。要计算一个测序运行中可以分析的样品的数目,可以应用以下公式:
这确保了每个测试运行的有效利用,同时确保在测序中有足够的读段用于被代表的系综。实现期望的冗余度所需的PCR循环数可使用与先前PCR运行拟合的模型来计算。首先,PCR效率可以通过将指数模型拟合到已知的cfDNA输入量来计算。然后,使用估计的参数,可以计算实现期望的过度扩增所需的扩增总数。
文库富集
在各种实施方案中,可以在测序之前使用文库富集以增加鉴定靶向区域中的变体的可能性。富集可以通过方法如靶向PCR(targeted PCR)或杂合体捕获组进行。靶向高通量测序可用于减少评估个体中指定基因座所需的测序读段总数。所需读段的减少是靶向序列长度除以基因组长度的商的函数,并且权重由靶向的和全基因组测序的分布测序读段覆盖深度(此后缩写为覆盖)确定。
增加的覆盖提高了灵敏度,因为含有靶等位基因的读段数目与真实变体比例(1-ε)x f和覆盖D大致呈二项式分布,其中ε是测序中的碱基错误率且f是分子群体中等位基因的频率。增加的覆盖可以通过使得能够跨越靶基因座聚合信息(整合出错误)来减少假阳性。由于测序中存在系统错误模式,例如均聚物中的错误,需要更复杂的错误模型。
选择基因组的哪些区域来靶向是设计靶向测序组的重要考虑因素。在使用遗传变体标志(signature)进行癌症检测的背景下,靶向组的统计功效(statistical power)是患者群体内跨这些基因座的变体重复出现(recurrence)的函数。在杂合体捕获设计中的另一个考虑因素是每个杂交探针的特异性和跨所有探针的灵敏度的均一性,两者都驱动在所需检测限值检测变体所需的测序读段的量。
本发明的系统和方法可以集中于选择多达总序列长度L的基因座的组合,其针对癌症患者中的最大组合重复出现负荷优化(组合驱动物和乘客遗传变体),考虑到影响杂合体捕获性能的决定因素诸如序列独特性和GC含量。此外,本发明可使用与cfDNA长度分布匹配并跨越靶区域上观察到的GC含量分布的合成的核酸添加物(spike-in)。基于指定的参考错配,添加物可以与cfDNA区分开来,选择错配的模式使得它们不可能从自然过程中观察到。这些添加物用于计算跨GC环境以及预测的杂合体捕获重叠的假阴性率估计。
本发明的杂合体捕获组可以通过以下来设计:鉴定为周期性体细胞突变(recurrently somatically mutated)(局部扩增、易位、倒位、单核苷酸变体、插入、缺失)和预先指定的基因座(诸如癌基因外显子)的区域,并选择提供最多信息的区域组合,直到指定的总组大小。可以通过以下考虑因素来设计杂合体捕获组:基因组长度、考虑的基因组改变和指定基因的强制包含;所考虑的肿瘤变异数据库和肿瘤类型、以及每个数据库的相对权重;对每种肿瘤类型的群体发病率的校正(以防止取样偏差);以及在外显子组或基因组水平的靶区域生成。
图6提供了根据某些实施方案的杂合体捕获组设计过程的图,包括数据转换。鼓形表示数据库,虚线框表示输入,菱形表示操作,且实心边框表示输出。
杂合体捕获组设计过程的输入可包括以碱基计的总允许组长度、预指定的将要靶向的区域、癌症类型的群体发病率的加权结果、留取(hold back)用于验证的样品的比例、对照添加物的数量、以及经验核酸长度分布。
参考数据库(DB)可包括以下:目标癌症类型的群体发病率、来自肿瘤测序的已知变体、诸如可以从基因组参考团体(http://www.ncbi.nlm.nih.gov/projects/genome/assembly/grc/human/)获得的人类参考基因组、来自健康群体的测序数据的已知变体、以及基因组独特性(例如,kmer比对可映射性和序列独特性)。数据库可以通过实验确定,并且可以通过应用本发明的方法将信息添加到数据库中。对数据库信息进行的操作可以包括在图6中在菱形内指出的那些操作。杂合体捕获组设计的输出可包括杂合体捕获目标集和以添加至样品或以其他方式用于评估跨鸟嘌呤-胞嘧啶(GC)含量分布的假阴性率的阳性对照。
为了鉴定基因组中提供最多信息的区域以靶向指定的总组长度,可评估记录在癌症突变数据库诸如COSMIC(癌症中的体细胞突变的目录http://cancer.sanger.ac.uk/cosmic)中的遗传改变的标志。优化可以使用前后优化(Forward-Backward optimization)或Greedy优化来进行。
杂合体捕获组设计可以使用交叉验证程序来验证,以说明通过从有限数量的样品构建组引起的潜在偏差。在设计癌症组时,交叉验证策略可以是重要的,因为样品中的遗传变异在肿瘤内(肿瘤内异质性)和患者之间(肿瘤间异质性)是异质的,并且受因素诸如遗传背景(例如,POLE突变状态)、环境暴露(例如,抽烟史、先前疗法)和肿瘤阶段影响。
对于前后优化,可以通过在前向和后向传递(forward and backward passes)之间交替直到从L基因座构建指定长度的组来鉴定基因座。基因座可以分为组中包含的那些(选择的位点)、以及未包括在组中的那些(可用的基因座)。对于每次迭代,在前向传递中,可以鉴定可用基因座中向组添加最大数量的体细胞突变的基因座f*。在后向传递中,f*可被包括在组集中,并且可鉴定包含的基因座中添加最少的体细胞重复出现的基因座b*。如果f*不等于b*,则可以排除b*。迭代可以重复。该方案可用于鉴定用于组合的体细胞重复出现的基因座的优化集。当达到组长度时,优化可以结束。
使用Greedy优化,该过程可以从添加最大体细胞突变负荷的基因座开始,将其添加到组,然后从剩余的基因座选择具有最大体细胞突变负荷的基因座。当组合的序列满足指定的组大小时,该过程可以终止。
交叉折叠验证(cross fold validation)可用于评估说明疾病数据库中结构的影响的所鉴定的组的稳定性。
在某些实施方案中,可以构建患者样品的两个相互排斥的集,其中集的基数由训练比例p确定。可以在第一集上生成具有基数p的组,基数p记录组上具有突变的患者的总数。然后可以在具有基数(1-p)的验证集中验证所提议的组,计算组上带有突变的患者的比例。如果患者比例在阈值T内,则该组可以被保留,并且如果比例不在T内,则该组可以被修改。
可以询问肿瘤活组织检查测序的数据库以获得遗传变异的样品,样品可以通过许多患者协变量诸如疾病类型、阶段、环境暴露和组织学来分层。然后可以去除在健康群体的群体测序中观察到的所有种系遗传变体,诸如1000Genomes数据库,以防止会使组设计混乱的癌症数据库中的假阳性变体(这一步骤仅当目标变体为疾病相关时有用,如在癌症诊断中)。存在使个体易患癌症的已知的种系突变,诸如BRCA1/2突变,其可通过这种方法消除,但是如果需要,可以将已知的感兴趣区域强制加入杂合物捕获组设计以克服这些遗漏。
为了说明杂交捕获中的差异性能,可以将关于人类基因组的序列特性的信息掺入到组选择过程中。在某些实施方案中,可将关于基因组中每个碱基的独特性的度量掺入在设计过程中,因为这驱动了杂合体捕获的特异性。例如,如果基因座与人类基因组中的99个其他基因座(例如,LINE元件)同源(相同),则捕获探针将仅下拉平均每100个中1个相关基因座。(使用的度量为1)。
可以通过使用从UCSC基因组浏览器数据库(https://genome.ucsc.edu/)可获得的基因组独特性的两个预先计算的汇总统计来掺入该信息。
可映射性,s,其定量kmer序列与基因组比对的独特性
独特性,u,在1个碱基滑动窗口中跨基因组的35个碱基窗口的独特性
其中x是准确的共享序列的数目。
可以组合这些映射,并且然后为人类参考基因组中的每个碱基生成字符(character)编码的独特性值。由此可以将参考基因组从核苷酸序列转换为通过杂交特异性得分f(s,u)注释的核苷酸序列。
一旦设计了杂合体捕获组,该组可用于使用其核苷酸序列富集样品的靶基因组区域。为了捕获分子,将双链DNA解链成单链DNA(例如,通过升高温度),然后加入杂交捕获探针(探针),并改变条件以促进链退火。探针与靶序列互补并具有使得分子能够被分离的选择标记(例如,生物素化的)。为了防止探针池中杂交探针之间的杂交,池中的所有探针都被设计为与靶基因座的有义序列或反义序列互补。因此,每个双链DNA分子仅捕获一条链。
在某些实施方案中,杂合体捕获组可以设计成特异性地靶向DNA的有义链和反义链二者。在杂交捕获之前对样品DNA进行PCR扩增的情况下,原始分子的两条链在有义和反义PCR重复群体中被代表。为清楚起见,考虑以下实例:x={x+,x-}是双链分子,α和β是长度为l的单链DNA分子,α的前n个连续碱基与β的最后n个连续碱基共有互补序列,其余序列是非互补的。因此,与β退火的α具有来自互补序列的双链DNA茎和来自非互补序列的单链DNA臂的叉状Y形结构。
接下来,可以使用已知的连接反应,例如钝末端连接产生分子,其中x的任一端侧翼为Y形{α,β}双链DNA:
αx-β
βx+α
然后可以使用与α和β互补的引物(分别表示为αc和βc)应用PCR扩增,以生成PCR重复家族:
将具有探针序列x-的杂交捕获方案应用于PCR产物将提取具有序列αcx+βc和βx+α的分子,其每个分别是反义链和有义链的后代。随后的PCR生成所有四个单链分子类,因此两条链被代表。然而,使用这种方法,一半的分子不被捕获。当前捕获方法仅使用有义探针或反义探针中的任一组,然而,链特异性分离可用于从原始取样的DNA生成两种相同地分布的样品。这种方法对于寻求在异质群体中检测低频率的分子的应用具有实用性,作为控制在取样的DNA的后续操作中引入的误差和丢失(dropout)的手段。例如,本发明的某些实施方案依赖于在有义链和反义链二者上被代表的候选变体用于验证。在这样的方法中,链特异性分离方法可以是特别有用的。可以使用以下步骤实现链特异性分离。
可以为感兴趣的基因座创建两个杂合体捕获组;一个有义(A)和一个反义(B)。然后可以将这些组连续地应用于DNA样品。可选择的探针可被应用于单链DNA,使用标准的杂交捕获方案将样品分离成分离物(被探针结合的DNA)分区和非分离物(未被探针结合的DNA)分区。可将组A应用于DNA群体。靶序列将被收集在分离物分区中。可以保留非分离物分区。可将组B应用于非分离物分区。可以在步骤2的分离物分区中收集组1靶序列的互补物。
可存在来自A的探针的一些携带污染,但如果优化了分离方法,这将是最小的。在可选择的实施中,样品可以分成两个等分试样,且A和B分别处理,从而避免在前一步骤中由探针携带引起的任何交叉杂交。
可以分别分析来自A和B的分离物,然后比较两个分析之间结果的一致性,这控制了在样品的下游处理中引入的伪相。这为分离物A和B之间的复制提供了机会,并通过分别评估A和B提高了灵敏度。
测序
样品最初可稀释至2nM,且然后在测序之前可稀释至600ul中19pM的最终浓度。合适的测序方法包括但不限于,通过杂交测序、SMRTTM(单分子实时)技术(PacificBiosciences)、真正的单分子测序(例如HeliScopeTM,Helicos Biosciences)、大规模并行下一代测序(例如SOLiDTM,AppliedBiosciences;Solexa和HiSeqTM,Illumina)、大规模并行半导体测序(例如IonTorrent)和焦磷酸测序技术(例如GS FLX和GS Junior Systems,Roche/454)。在优选的实施方案中,测序可以是通过合成测序技术(例如,HiSeqTM和SolexaTM,Illumina)。样品可以上样到HiSeq系统上。Illumina流通池上的读段簇的密度可以针对cfDNA进行优化,尤其是通过读段的长度分布驱动,并且簇密度可以通过测序各种上样浓度来通过实验优化。每个池可以上样的样品的数量可以通过计算每个测序运行的有效利用的分析公式来定义:这是可以同时运行的样品的最大数量,使得实现期望的过度扩增因子。上述浓度导致HiSeq2500上的最佳簇生成。然而,如果未获得在快速运行中850-1000K/mm2的所需簇生成,则上样浓度可以相应地改变。
分析
本发明的系统和方法是基于以下洞察:高准确度PCR酶比下一代测序机更不容易出错:如果目标是高保真测序,则创建每个单个分子的多个拷贝,对这些分别测序,且然后创建共有序列是个好主意,所述共有序列反映原始分子的序列并对在测序过程中产生的(大多数)错误取平均值。该方法的一个主要挑战是根据它们所源于的原始分子来对测序的分子进行分组。这可以通过在扩增之前用随机核苷酸序列对原始分子进行生物化学标记来实现,使得假定共享相同标记序列的所有测序分子来自相同的原始分子。在本发明的优选实施方案中,测序的分子可以在没有生化标记的情况下分组;相反,统计和生物信息学方法可用于鉴定每个原始分子的祖先。
这些概念可应用于鉴定可能包含突变的(低频)等位基因的BAM比对文件的列。BAM格式是用于存储序列数据的二进制格式。系综一致性检查的概念可以应用于通过寻找相容序列的系综中链平衡一致性来检查从文库构建的de Bruijn图中鉴定的推定变异。
根据本发明实施方案的系综是对齐的读段对的集合。在一些实施方案中,系综包括共享相同的起始坐标和终止坐标的对齐的读段对的集合。换句话说,对于每个读段对,存在读段对的碱基与其对齐的参考基因组坐标的一组坐标;每个这样的集具有最大值和最小值;系综是具有相同最大值和相同最小值的读段对的集。在一些实施方案中,系综包括具有近似相同的起始坐标和终止坐标的对齐的读段对的集合。忽略测序误差,单个系综包含源自参考基因组中具有相同或近似相同的起始/终止坐标的原始分子的PCR产物的读段。重要的是,原始分子的两条链应该被系综的成员代表,并且可以通过检查两个来源链是否是形成系综的“左侧”(意思是:较低的参考坐标)的第一或第二读段(在Illumina配对末端范例中)来区分它们。
上面讨论的过度扩增因子可被认为是根据来源于每个原始分子的读段的平均数。如果测序和PCR是完全的并且所有原始分子是独特的,则每个系综的读段数将等于过度扩增因子。
尽管可以通过实验确定过度扩增因子,在优选实施方案中,可以从输入BAM文件统计地估计过度扩增因子。估计程序可以基于以下洞察:大多数原始分子是独特的,并且大多数系综因此应包含与过度扩增因子类似数目的读段(即,过度扩增因子的第一近似值可以通过确定直方图的模式来计算,该直方图绘制x轴上每个系综的读段数与y轴上具有该读段数的系综的数量)。
为了将一组读段对对齐转换为系综列表(其中每个系综包含一组读段对对齐),可以使用上面给出的系综定义:具有相同的最大/最小坐标的所有读段对对齐变为同一系综的一部分。重要的是,该定义是基于完整对对齐的最大值/最小值,而不是基于2个单独读段的最大值/最小值(即,2个单独读段对齐的“内”端可以忽略)。在读段对齐的开始和结束处的测序误差(在对齐的坐标中,对应于由机器产生的两个单独成员读段的开始)将导致错误的读段形成它们自己的系综。另外,基于诸如以下标准,只有满足一系列一致性标准的读段对被考虑:
●两个成员读段被映射并在同一染色体上。
●读段1和读段2的相反链。
●两个对齐的最小值和最大值之间的总距离<常数(即假定的潜在分子的长度<常数;常数通常为~330)。
●由对准器(Burrows-Wheeler;Li H.和Durbin R.(2009)Fast and accurateshort read alignment with Burrows-Wheeler Transform.Bioinformatics,25:1754-60.[PMID:19451168],通过引用并入本文)设定的所有质量控制标志对于两个读段对齐OK;有一个“QC”标志和一个“正确对”标志。
●最小映射质量(诸如>0.95)。
●跨读段的错配的比例<常数(单独测量的)。
●所有读段成员都没有软剪切或填充。
●基因组一致性:基于原始分子来自相对正常的人类基因组的假设,要求与参考基因组的正链对齐的(读段对的)读段是另一个读段的“左侧”(通过其各自对齐的最小坐标测量)-且反之亦然。
通过检查系综的“左侧”读段(如上定义)是读段对的第一还是第二读段,可以确定来自原始分子系综成员的两条链中的哪一条。在一个优选实施方案中,使用其中一对的两个读段具有连续的对齐的对齐算法(例如,非分割读段对齐算法)。在一些实施方案中,使用分割读段对齐算法(例如,bwa mem)。
包括数据分析的本发明的方法可以由包括耦合到处理器的有形、非瞬时存储器的计算机来执行。从输入BAM文件开始,可以使用计算机进行以下一个或更多个分析步骤:
1.系综枚举:
鉴定存在于BAM中的所有系综,并且可以将它们的坐标(和协变量如长度、GC含量和成员读段的数量)写入文本文件(例如,clusters.txt)。输出文件后,可以从工作存储器中删除所有系综数据。
2.过度扩增的统计学估计:
考虑到协变量如GC含量、系综长度、与下拉探针的重叠,可以调用读取clusters.txt并估计过度扩增的统计模型的计算机脚本(例如,R脚本)。还估计了沿输入分子长度的分布和输入分子基因组覆盖。
3.确定性分析:
可以迭代BAM文件的所有列,并且鉴定可能包含突变的等位基因的那些列。列中的每个等位基因是簇的成员,且将等位基因按簇成员资格和它们来自原始分子的哪条链分组。用于鉴定具有可能突变的列的阈值考虑了来自统计学过度扩增模型的估计值。
概率分析:
对于候选列或对于所有列,可以应用完整的PCR扩增模型,其明确考虑扩增错误的不同情景(在PCR的不同循环,以及相对于原始分子的不同链)并且将它们的可能性与突变的输入等位基因的不同情景进行比较。
确定性和概率分析算法可以是基于列的,即,它们鉴定推定包含突变的等位基因的BAM对齐文件中的列。
可以对每个单独的读段等位基因分配全局有效系综ID,或者可以“在运行中(on-the-fly)”构建系综ID。“在运行中”生成的系综ID只能在每个BAM对齐列中被假定为独特/有效,并且它们对于“全局”系综列表没有确定的含义。这些函数可以基于回调:也就是说,它们将函数参考作为它们将为BAM对齐中的每个列调用的自变数(argument)。
它们也可以是多线程的(即,基于任何合适的并行化框架并行化(例如,使用openMP)),并行处理BAM文件的不同段。回调函数优选地不尝试访问全局变量,或使用受保护的存储器访问。回调函数也可以接收它们从其调用的线程数作为自变数,它也可以用于避免并发存储器访问的构建中(例如:如果有16个线程,则构建具有16个元素的矢量,并且每个线程只访问其相应的元素)。
如回调函数所见,列可以被建模为等位基因环境对象(allele contextobject)的矢量,其中每个等位基因环境对象代表对齐中的一个读段。通常,一个读段相当于一个碱基,但如果存在局部插入,等位基因环境对象也可以包含多于一个碱基。除了原始读段碱基之外,等位基因环境对象还可包含相关的碱基质量、关于对齐的另外信息(映射质量、读段中的位置、第一或第二读段等),以及,重要的是,指明读段属于哪个系综的系综ID(此ID是局部或全局独特的,见上文)。
用于构建这些列矢量的基本潜在算法是一个密集的过程,且可如下工作:
●对于参考基因组区域,可以为每个参考基因组位置构建空等位基因环境对象。
●对于同一区域,可以获得来自BAM文件的所有读段对齐。
●然后可以缩小所有对齐;也就是说,可以计算原始读段中的哪些碱基与哪些参考基因组位置对齐(该信息以CIGAR串或碱基长度序列和相关操作编码)。
●对于读段中的每个碱基,碱基(及潜在地关于其的另外信息)可被附加到其相应的参考基因组位置矢量。
可以在每列的基础上应用确定性算法并使用上述BAM访问功能。确定性算法的目的是鉴定推定含有突变的等位基因的混合物的列。分析算法可以如下起作用:
●低频等位基因可以在列中找到,其被处理为潜在的混合物等位基因。
●对于每种潜在的混合物等位基因:
●列中的所有等位基因可以根据其系综ID分组。
●对于每个系综:
o可以对来自系综中的读段的变体等位基因,分别对潜在分子的正和负链(即,其中读段对的第一个读段的对齐从系综的左侧开始的链),计算支持(即变体等位基因频率)。
o每个系综代表多个原始分子,其可以通过将系综中的总读段除以平均计算的过度扩增因子来估计。
o系综可归类为“推定地含有变体等位基因”,如果
■推定的变体等位基因的频率>=(1/'对潜在分子数的估计值')x因子,其中因子可以是系数诸如0.9。
●对于源自原始分子的正和负链的读段,可以要求分别满足该标准。
■此外,可需要观察来自两个单独链的最少数目的读段。在优选实施方案中,对于每个原始链可需要至少2个读段。
●如果有至少一个系综被归类为“推定地含有变体等位基因”,则列可归类为“推定地含有变体等位基因”。
概率算法也可以基于每个列应用。该算法的目的是计算列包含突变的等位基因的混合物的假设的证据强度。因此,在用确定性算法鉴定候选物之后,优选地采用其作为第二步骤(概率算法可以是计算上昂贵的,因此通过初始筛选来最小化其应用可以是期望的)。然而,该算法也可以单独使用,而无上文的确定性算法。在某些实施方案中,概率算法涉及确定候选变体是真实变体的可能性。概率算法可以使用任何已知的似然最大化模型,诸如例如,期望最大化、极大似然、拟极大似然、极大似然估计、M-估计、广义矩方法、最大后验、矩方法、支持方法、最小距离估计、限制极大似然估计或贝叶斯方法。
在优选实施方案中,概率算法可以如下应用:
●对于每个列,可以恢复存储的列数据。
●可以鉴定推定的突变(例如,通过找到低频变体等位基因,如在上面的确定性分析中)。
●等位基因可以通过系综ID在列中分组。
●对于每个假定的突变:
o观察的数据的似然可以在变体等位基因的频率为0的假设(H0)与指定非零变体等位基因频率的一系列假设下计算,其中特定频率在经验(列范围的)观察的变体频率周围生成。
o似然计算可在每个系综的基础上继续,其中假设系综是独立的(以指定的变体等位基因频率参数为条件)。为了获得在假设下观察的数据的似然,在该假设下可以乘以每系综的似然。
o然后可以选择具有最高似然的非零变体频率假设(类似于使变体频率的似然最大化,但是具有减少的计算应变)并且可以针对H0执行似然比检验以获得p值。
●对于每个列,可以报告具有最低p值的推定突变。
系综的似然可以在存在具有指定频率(其可以是0)的变体等位基因的假设下计算。由于列的似然被计算为每系综似然的乘积,这里描述的方法可形成概率分析方法的核心。每个系综源自未知数量的潜在分子。在系综中观察到的变体等位基因可以源自真正突变的潜在分子,或者它们可以由于测序和PCR错误而出现。真正突变的等位基因应该在源自原始分子的正和负链的读段上同等地被代表。根据PCR错误在其中发生的PCR循环,它们具有不同的结构(更早期的错误影响更多的分子)。测序误差被假设为随机发生(即,没有关于它们的特定结构)。
用于区分这些情景的统计模型可以基于完美PCR效率的假设,即,每轮PCR导致原始分子的加倍。这意味着,原始分子及其衍生分子的每个链可以表示为分叉树(即,对于每个原始双链分子有两个分叉树)——节点代表分子并且边缘代表PCR扩增过程。树中的级别数等于PCR轮次数+1(原始分子节点表示级别1)。可以假设错误模型作用于树的边缘,即每个边缘代表准确的扩增或错误。如果发生错误,其影响受影响的边缘下方的所有节点。错误在“非变体”和“变体”之间翻转分子的等位基因状态。树的尖端代表PCR扩增后的分子,即进入测序机的分子群体。由于每个系综源自未知数量的原始分子,每个系综可以与未知数量的分叉树相关联。
基于未知数量的未知分子、未知数量的错误等,存在无限数量的可能情景。为了限制所考虑的替代方案的空间,可以做出以下实用假设:每个系综的潜在原始分子的数量在1到8之间;PCR循环次数为4;并且一个原始分子在扩增期间的最大错误是2。这些假设可用于如下限制所考虑的情景数量:
●x=1-8之间的原始分子(即2-16个分叉树)
o其中y=0<=x可以“真实地”携带变体等位基因
■所有这些经历4个PCR循环
●对于完整系综,沿着所有树z=0-2个错误
o z个错误的每个落在一个定义的边缘上。
对于每个系综,总似然可以分成2个部分:系综中存在的读段总数、以及分别源自原始分子正和负链的读段中的变体等位基因频率。该因式分解可用于达到另一种简化。
可以假设,分别对于原始正和负链(“error_strand(错误_链)”),通过指定错误发生在树的哪个水平(“error_level(错误_水平)”),以及它是否影响携带变体等位基因的分子的祖先(“error_variant(错误_变体)”),可以对PCR错误对变体等位基因跨树的尖端节点(即,测序的分子)的频率的影响进行建模。“场景”的正式定义可以作为x、y和z值(在上面指定的边界内)加上对于z个错误的每个的(error_strand(错误_链),error_level(错误_水平),error_variant(错误_变体))集的组合提供。对于完整的概率评价,可以计算在所有场景下的数据的似然。
分别针对正和负链衍生分子,每个场景在所包含的树的尖端水平具有关联的变体等位基因频率,以x和y以及错误集为条件。可以使用计算机来如下处理此信息:
●F_mutatedAllele_plus(F_突变的等位基因_正)可定义为突变的等位基因跨起源于原始分子的正链的系综成员的频率(在所考虑的情景是真的的假设下),且F_mutatedAllele_minus(F_突变的等位基因_负)可定义为突变的等位基因跨起源于原始分子的负链的系综成员的频率。
●F_mutatedAllele_plus(F_突变的等位基因_正):=F_mutatedAllele_minus(F_突变的等位基因_负):=y/x然后可以初始化。
●对于定义为(error_strand(错误_链),error_level(错误_水平),error_variant(错误_变体))的z个错误的每个:
o levels_downstream_affected(水平_下游_受影响):=roundsPCR(PCR轮)-error_level(错误_水平);(基于1的水平索引,即第一轮PCR中的错误具有水平1)。
o如果error_strand(错误_链)=“+”:
■如果error_variant(错误_变体)=“non_variant(非_变体)”:
F_mutatedAllele_plus(F_突变的等位基因_正)=F_mutatedAllele_plus(F_突变的等位基因_正)+oneMutation_effect(一个突变_影响)
■如果error_variant(错误_变体)=“variant(变体)”:
F_mutatedAllele_plus(F_突变的等位基因_正)=F_mutatedAllele_plus(F_突变的等位基因_正)-oneMutation_effect(一个突变_影响)
o如果error_strand(错误_链)=“-“:
■如果error_variant(错误_变体)=“non_variant(非_变体)”:
F_mutatedAllele_minus(F_突变的等位基因_负)=F_mutatedAllele_minus(F_突变的等位基因_负)+oneMutation_effect(一个突变_影响)
■如果error_variant(错误_变体)=“变体”:
F_mutatedAllele_minus(F_突变的等位基因_负)=F_mutatedAllele_minus(F_突变的等位基因_负)-oneMutation_effect(一个突变_影响)
●F_mutatedAllele_plus(F_突变的等位基因_正)和F_mutatedAllele_minus(F_突变的等位基因_负)可被限制于0和1的边界。
在各种实施方案中,对于z个错误中的每个,程序可以任选地仅指定a.它是否影响携带变体等位基因的分子的祖先(“error_variant(错误_变体)”);b.它是否影响正或负链原始分子的祖先(“error_strand(错误_链)”);和/或c.错误的树水平(“error_level(错误_水平)”)。在某些实施方案中,该程序可以指定a.错误影响了1..X个分子(+祖先)中的哪一个;b.其是否影响原始正或负链的祖先;和/或c.精确地在相应树的哪个边缘发生了错误。
为了在每个考虑的情景下计算数据的似然,可以获得先前情景似然并且乘以该情景下的数据的似然。可以如下给出每个场景的先验概率:X可具有来自过度扩增计算机脚本的统计估计的输出的概率分布,考虑到原始分子基因组覆盖,条件是系综的长度(例如,较长的系综具有仅来源于一个原始分子的更高的机会)。y可以具有(泊松)概率分布,由假设的变体等位基因的频率参数化。z,错误的总数,可以具有(泊松)概率分布(来自根据边缘数量缩放的PCR酶的实验估计的错误频率),并且假设每个边缘同等可能地被错误命中(即,携带变体等位基因和非变体的原始分子的祖先以与情景中的这些分子的数量(变量x和y)成比例的概率被命中)。在此情景中仅有的追踪的考虑因素是错误是否命中变体/非变体分子祖先树、它是否命中正/负链树、以及它命中了哪个水平(如上所述)。
系综的数据可以基于情景给出似然。可以注意到,系综数据由具有相关的质量值(通常是FASTQ碱基质量)的等位基因组成,并且每个等位基因与变体等位基因相同或不同(‘非变体’)。此外,对于每个考虑的情景,在树的尖端水平处的变体等位基因的频率可以代表原始变体和非变体分子的正和负链的祖先。
使用这些频率,观察到的系综数据可以被建模为伯努利分布(分别对于正和负链祖先),在单独等位基因碱基质量上整合。
对于给定变体等位基因频率,类别likelihoodTree<int roundsPCR、intmaximumUnderlyingMolecules、int maximumErrors>代表所有情景的集合。也就是说,为了进行完全分析,可能需要将H0(变体等位基因频率=0)与多个假设、多个likelihoodTree对象进行比较。基本情景参数,诸如PCR的轮次、最大潜在分子和每个系综的最大错误数量,可以表示为模板自变数,从而实现高效的编译器优化。
类别likelihoodBranch<int roundsPCR,int maximumErrors代表单独的情景,其由以下信息组成:
●潜在分子的总数
●这些潜在分子中的多少携带变体等位基因
●存在多少错误:
o在代表正链非变体原始分子的祖先的树的每个水平上
o在代表正链变体原始分子的祖先的树的每个水平上
o在代表负链非变体原始分子的祖先的树的每个水平上
o在代表负链变体原始分子的祖先的树的每个水平上
方法likelihoodBranch::likelihood_data(..)(似然分支::似然_数据(..))可计算likelihoodBranch对象所代表的情景下的一个系综的似然。likelihoodTree对象需要被所有一致的likelihoodBranch对象填充。函数likelihoodTree::computeErrorConfigurations(..)(似然树::计算错误配置(..))计算所有一致的情景,其然后(在构建者likelihoodTree中)转换为likelihoodBranch对象。每个情景的先验概率也可以在likelihoodTree构建者中计算。
对于指定长度、GC含量等的观察到的系综以及以特定数量的读段,R组分可帮助确定沿潜在分子数的概率分布。为了回答这个问题,可以推导出以下数量的估计值:
●以系综长度为条件,推断沿着潜在分子数的先验概率分布。这种分布受潜在分子的总基因组覆盖及其长度分布(因此需要被估计)的影响。
以假设存在特定数量的潜在分子为条件,推断沿着由这些潜在分子生成的读段的概率分布。这种分布受到过度扩增过程的性质的影响,过度扩增过程被假定为独立作用于原始分子并且假定为遵循泊松分布。
对于每个单独的原始分子,泊松的平均值可通过用截距(Mu)和以下的系数的线性函数(的指数)来参数化
o系综的长度(Length(长度))。
o系综的GC比例从0.5的偏差(GCm50)。
o如果已经应用了下拉捕获富集,系综与最近的下拉探针具有少于90个碱基重叠的程度(PulldownLess90)。可以对其他富集方法如前面所讨论的那些构建类似的度量。
上述数量估计可使用每个系综的潜在分子数的概率分布来进行。
这个概率分布可以形成矩阵,其中系综在行中且可能的潜在分子数为列,其中每行总和为1。这个概率分布可通过考虑沿着每个系综的读段的直方图来初始化:在来自血浆的cfDNA测序的应用中,大多数分子可以被认为是独特的(如使用从全基因组无PCR无细胞DNA测序获得的测序数据的分子长度分布的计算机模拟所示),因此,大多数系综可具有等同于其实现的过度扩增因子的数目的读段。为了考虑影响协变量,可以通过协变量值(以多维分位数(quantile))对系综数据进行分层,并且然后可以分别对每个分位数进行程序。这为每个系综提供了第一猜测过度扩增因子。可以通过假设观察到的读段计数遵循泊松分布来填充矩阵,其中平均值等于number_underlying_molecules(数量_潜在_分子)xover-amplification_factor_of_ensemble(系综_的_过度扩增_因子)。矩阵可以以逐行方式用获得的似然填充,并按行归一化。这为每个系综提供了沿着潜在分子的概率分布的第一近似。
分布可通过采用期望最大化(EM)样程序来完善,以完善概率矩阵。在此过程中可以进行一些简化的独立假设。
对于EM算法,可以保持观察到的读段计数遵循泊松,具有meannumber_underlying_molecules(平均值_潜在_分子数)xover-amplification_factor_of_ensemble(系综_的_过度扩增_因子)的假设,但over-amplification_factor_of_ensemble(系综_的_过度扩增_因子)可被替换为exp(过度扩增(Mu,Length,GCm50,PulldownLess90)),其中过度扩增(Mu,Length,GCm50,PulldownLess90)是单独分子的过度扩增因子的线性预测物。可以针对每个系综单独计算过度扩增(Mu,长度,GCm50,PulldownLess90),考虑到全局系数以及系综的GC含量、下拉重叠等的单独值。
对于EM部分,可以在矩阵的列上引入先验概率,条件是系综长度(即,每个系综具有其自己的逐列先验)。这些先验概率取决于在基因组的每个位置处原始分子的起始率(覆盖)和分子长度分布,这些是也可以估计的量-并且假设独立于过度扩增协变量,条件是固定的每个系综的潜在分子数量概率分布。下文更详细地描述估计程序。
EM样算法可以如下构造:
1.初始化E=clusterData_P_underlying(簇数据_P_潜在)
2.(M步骤)保持E固定,估计潜在分子基因组覆盖、长度分布和沿着潜在分子数的先验分布,以系综大小为条件。
3.(M步骤)保持E固定,估计Mu、Length、GCm50、PulldownLess90。
4.(E步骤)保持Mu、Length、GCm50、PulldownLess90和潜在分子先验分布固定,从观察到的数据(每个系综的读段数)估计E。
5.在E和所有参数估计值下测量观察到的数据的似然;如果有足够的改进,转到步骤2,如果没有足够的改进则中止。
估计基因组覆盖和潜在分子的长度分布和对每个系综的潜在分子数的先验概率可以使用填充的矩阵来完成,该矩阵指定沿着每个系综的潜在分子数的概率分布。每个位置的潜在分子的起始率可以估算,然后可以估算长度分布,然后可以以系综长度为条件估算先验分布。
起始率/覆盖估计:
可以鉴定用于测量覆盖的第一位置。在某些实施方案中,可以仅测量显示出与下拉探针的足够重叠的位置处的覆盖(或更精确地:从这些位置开始的假设cfDNA分子与下拉探针的重叠需要足够)。如果鉴定出过多的位置,则系综数据可被向下取样,以仅包括从位置的子集开始的系综(即:移除不在这些位置之一处开始的所有系综)。该子采样可以在进入算法的EM部分之前执行一次,并影响估计程序的所有步骤,包括Mu、Length、GCm50、PulldownLess90的估计。分子的起始率的估计可以通过鉴定在选择的位置之一处开始的所有系综并在其预期的潜在分子数上求和来推导。然后可以将该数除以所考虑的位置数。如果需要,随后可以通过乘以平均分子长度来获得覆盖。
长度分布估计:
对于每个系综,可以推断潜在分子的期望值。然后可以计算系综长度的加权平均值(通过每个系综的潜在分子估计进行加权)。可以内插缺失值(例如,由“Coverage”部分期间的子采样引起的)。
在指定长度上每系综潜在分子数的先验分布的估计:
每个位置的起始率和长度的分布可以使得能够计算在指定长度的系综潜在的原始分子的数量的先验概率。首先人们可以在(与长度无关的)一个位置迭代可能的起始分子的数目x,且然后(从长度分布)计算这些分子的y<=x=1,2,3…等于我们系综的指定长度的概率。然后可以在可能的(x,y)值上对概率分布进行归一化,并对y进行边际化。
根据某些实施方案,本发明的系统和方法可包括模拟器。模拟器功能可以采用输入,其指定参数诸如覆盖、突变的等位基因混合物和所选择的箱。两个最重要的参数是PCR前“原始cfDNA”产物的覆盖和设想的测序数据覆盖。(在我们感兴趣的区域测量,见下文)。PCR前“原始cfDNA”产物的覆盖包含来自突变的亚克隆(见下文)的分子以及非突变的分子。两个参数之间的延伸可用于确定过度扩增因子。在某些实施方案中,模拟过程可以通过以下性质表征:
●模拟的基因组区域可限于由下拉组捕获的区域。
●许多突变沿着指定区域延伸。每个突变具有相关的混合物频率(其存在于我们的模拟的cfDNA中的频率)。每个混合物频率可以被视为单独的亚克隆,且因此属于一个箱的所有突变被一起模拟(即,如果它们彼此足够接近,它们将形成单体型)。
●可以创建代表总cfDNA产物的分子池(即,包括突变的和非突变的片段)。该池可以通过分别模拟来源于非突变的参考基因组和来自指定的亚克隆(即:来自指定的混合物频率)的分子来填充。如果分子来源于非参考亚克隆,则它(如果它重叠)携带与其来源亚克隆/混合物频率相关的突变。通过在不同的亚克隆(具有指定的混合物比例)和非突变的参考基因组(接收剩余的、非混合比例)上延展PCR前产物的总期望覆盖,可以确定模拟程序的每个部分的总覆盖。
●以下是覆盖如何在亚克隆和非突变的参考基因组之间延展的两个例子:
o如果指定了原始分子的所需的总覆盖为1,000x,并且如果有一个亚克隆/混合物箱具有10%频率,则获得以下覆盖:900x为“非突变的”,且100x为“10%混合物”亚克隆。
o如果添加了具有1%频率的另外的混合物箱,则发现以下分子箱:890x为非突变的,100x为10%突变的,10x为1%突变的。
●具有预定义序列的对照序列也可被添加到分子池中(作为创建池之后的第一个步骤)。每个对照序列可被多个相同的分子代表,并且每个对照序列的相同分子的数量可以从泊松分布中得出(其平均值可以是用户指定的,并且对于不同的对照序列可以是不同的)。
●填充分子池后,可以模拟P5和P7衔接子的连接和PCR扩增(分别对于正和负链并保留连接的P5/P7分子的方向)。模拟可以在池上进行,即,池中分子的数量随着每个模拟轮次而增长。PCR过程模拟可包括模拟和测序错误以及不完全扩增。不完全扩增的概率可以对池中的每个分子单独计算,并且取决于分子的GC含量。可以从期望的测序读段覆盖和指定的测序效率计算PCR循环的数量。PCR后分子池中所需的覆盖可以通过将所需的测序覆盖乘以1/指定的测序效率来计算。然后,考虑到(PCR前池中的分子的)平均扩增效率,人们可计算需要多少个PCR循环以将池中的覆盖从PCR前水平带到期望的PCR后水平。
o为这一计算提供实例,如果设想的总覆盖是160,000,并且总原始分子估计为20,000(即原始分子在感兴趣区域上的覆盖为20,000x),且PCR效率为100%,并且指定的测序效率为0.5:需要PCR后池中320,000x的覆盖,且这要求4个PCR循环。
●最后,人们可以从池中的PCR后分子取样(以测序效率率(sequencingefficiency rate)),并为每个选择的分子生成配对末端测序读段。P5/P7连接方向确定分子的哪一端生成第一读段。测序读段的生成可以包括测序错误的模拟。
模拟器可以保持追踪许多重要事件,例如,PCR错误的位置和时机(哪轮PCR)。这些数据可以作为文本文件存储在模拟输出目录中。
在测序读段的模拟完成后,可将模拟的读段映射到参考基因组。映射完成后,可以对数据进行分析并用于产生有多少模拟的突变被调用以及有多少假阳性的分析。该输出可被发送到输入/输出设备,例如打印机或显示器。
在优选实施方案中,测序数据的分析可以以BAM文件作为输入数据开始,其中输出是一个或更多个文本文件。
对照
在某些方面,本发明的系统和方法涉及使用样品中的体细胞改变来估计测序错误和非均匀覆盖对变体等位基因频率估计的影响。为此,可以鉴定具有与种系基因组不同的N个连续碱基(N>1)的体细胞变体,并由矢量V={a(1),a(2),…,a(n-1),a(n)}表示,其中元素a(i)代表变体中在位置i处的不同碱基。这种变体可以通过体细胞改变生成:易位、倒位、插入、缺失、扩增。
对于变体中的每个碱基a(i),对支持该碱基的等位基因的总数进行计数,这生成了V的频率f(V)的n个估计。所有观察到的频率f(a(i))应该等于f(V),但是由于覆盖的变化和测序错误,情况可能并非如此。然后,可以使用已知的统计方法来定量在测序期间产生的频率估计中的离差。然后,这可用于校正频率估计。一个例子将是使用样品平均值和方差来使用适当的采样分布来估计置信区间。
在二倍体生物中,等位基因在杂合位点处的比应为1/2。在人类群体中存在SNP分离的大数据库。对于给定的个体,可以询问这些位点,并将杂合位点鉴定为具有大致相等的等位基因频率的两个等位基因的基因座。然后可以从在杂合位点处的第二等位基因的观察到的频率构建等位基因频率的经验分布。如果杂合位点的数量足够大,可构建每个等位基因组合(A>C、A>G、…、T>G)的频率估计。然后可以使用分布来校正样品数据中在体细胞变体位点处的频率估计。
在某些实施方案中,可以将具有与患者不同的序列的已知输入量的DNA添加到样品中。这些是样品中变体等位基因的阳性对照。为了生成可鉴定的添加物,可生成在人类群体中不可能观察到的序列。这可通过以下来完成:1)选择在群体测序数据库中具有低报告多样性的区域,2)向序列引入不反映自然突变过程的变化(例如序列(相同)n,{变化,相同,变化,相同,变化},(相同)n)。对照序列可被进一步区分,因为添加物的长度(120个碱基)是已知的,且所引入的变化的位置也是已知的。
已知,杂交捕获可受捕获探针和靶DNA之间的错配数量的影响。在某些实施方案中,然后将4个突变引入每个对照中。也可以构建添加物以便1)GC含量和2)探针-靶重叠的影响可以通过以下观察:1)选择跨靶向区域与已知GC含量分布具有不同的GC百分比的序列,和2)改变120碱基长的对照DNA与其相应的下拉探针的重叠百分比。
可以在血液抽取之前将添加物加入到血液收集真空采血管中,以便a)样品可以从它们的测序中鉴定,允许在测序中鉴定样品混合,b)以便可以估计来自有核白血细胞凋亡的污染,并且c)以便可以检测假阴性。
除了绝大部分来自人的正常(通常是健康的)基因组的分子以外,来自人血浆的无细胞循环DNA(pDNA)还包含癌症患者中的肿瘤DNA的片段和孕妇中胎儿DNA的片段。调查肿瘤或胎儿DNA的混合部分本质上具有挑战性,因为癌症/胎儿衍生分子的混合比例可以低至5000个分子中之1。
任何给定的未加工的血液样品(通常但不总是储存在EDTA管或不同类型的血液收集容器中)将含有一定比例的无细胞DNA以及白血细胞和红血细胞(WBC和RBC)。经过一段时间后(并受环境因素诸如温度的影响),所含的WBC将经历细胞死亡并开始将所含的DNA片段释放到循环中。由于该过程,血液样品中包含的任何肿瘤或胎儿衍生的无细胞DNA将进一步稀释,使得它们的检测和表征甚至更具挑战性。
存在技术解决方案(例如Streck管)以防止所含的WBC破裂并释放其DNA,但是这些解决方案并不完美并且稀释问题仍然存在,特别是如果血液储存较长时间段或当运输血液样品时。
对于基于调查肿瘤或胎儿衍生DNA的存在或特征的任何诊断方法,因此期望测量和控制潜在的污染。
在本发明的某些实施方案中,可将合成的扰乱DNA掺入收集容器中以追踪污染。可确定人类基因组中的段或区域,即:a)在人类群体的绝大多数中为纯合的,即具有已知和/或可确定的频率阈值(或在期望的靶群体的绝大多数中为纯合的)和b)基因组复杂性高,即使用用于读段对齐的标准算法方法,对于从该区域衍生的分子建立基因组起源是明确且无挑战的。通常,该段的长度在50和150个碱基之间变化,但是此处描述的方法可使用更长和更短的区域。然后可以通过用不同核苷酸取代一定数量的核苷酸或引入或缺失一定数量的核苷酸来扰乱段或区域的序列。通常,该步骤将包括用不同的核苷酸取代位于序列中央的一个或两个核苷酸。接下来,可以证实扰乱的序列不存在于正常人类群体中。存在多种标准方法来实现这一点,例如基因组比对或与由群体测序数据,诸如1000Genomes Project生成的de Bruijn图比较。如果此验证失败,可重复步骤2或1。
然后可以使用DNA合成方法合成扰乱序列以产生(近似或精确)n个拷贝的如此扰乱的序列。通常选择数量n,使得当n个分子被引入患者样品(在这种情况下为血液)时,n与所抽取的血液体积中人类基因组的预期拷贝数之间的比类似于肿瘤/胎儿衍生片段与正常基因组片段之间预期的/所需的最小比。(例如,如果人们预计1000个循环片段中的1个是肿瘤起源的,并且如果1ml血液通常包含约1000个拷贝的人类基因组,并且人们抽取5ml血液,则每个管n=5将是明智的选择)。
扰乱序列的合成拷贝可以在收集之前存在于收集容器中,或者可以在收集后添加到样品中。将合成的扰乱DNA在时间X与样品接触。在样品分析期间,可以通过离心提取无细胞循环DNA,并且可以从提取的DNA制备DNA文库。可以使用将在样品的下游解释中使用的技术(例如,基于数字PCR的方法或基于测序的方法,使用全基因组测序方法或靶向测序方法)来测量观察到的扰乱序列的频率(fP)和非扰乱序列的频率(fn)。
可以如下分析观察到的频率:fP/(fP+fn)是样品中以n个拷贝原始地(即,在由于WBC破裂的稀释开始之前)存在的肿瘤或胎儿来源的等位基因的稀释后频率的估计。依据所采用的下游解释技术的特性,如果fP/(fP+fn)为0或低于指定的阈值,则应拒绝或不解释样品。将稀释后数据中假定的肿瘤或胎儿衍生等位基因的观察数频率乘以([(fP+fn/fP)]x n将得出该等位基因的稀释前绝对计数的估计值。肿瘤等位基因计数及其随时间的发展已被证明是疾病状态和进展的重要指标。使用预接种的收集容器的上述过程在图7中举例说明。
上述程序可用于不同的基因组基因座和不同的n值,以提供另外的优点,诸如控制GC含量偏倚和实现稀释总量(更准确的)估计(在稀释衍生的分子片段中测量)和因此血液样品中DNA片段的稀释前数量。
如本文提及的计算机通常包括通过总线耦合到存储器和输入输出(I/O)机制的处理器。存储器可包括RAM或ROM,并且优选地包括至少一个有形的、非暂时性介质,其存储可被执行以使系统执行本文所述的功能的指令。如本领域技术人员将认识到对于执行本发明的方法所必要的或最适合的,本发明的系统包括一个或更多个处理器(例如,中央处理单元(CPU)、图形处理单元(GPU)等)、计算机可读存储设备(例如,主存储器、静态存储器等)或其组合,其经由总线彼此通信。
处理器可以是本领域已知的任何合适的处理器,例如由Intel(Santa Clara,CA)以商标XEON E7销售的处理器或由AMD(Sunnyvale,CA)以商标OPTERON 6200销售的处理器。
根据本发明的输入/输出设备可以包括视频显示单元(例如,液晶显示器(LCD)或阴极射线管(CRT)监视器)、字母数字输入设备(例如,键盘)、光标控制设备(例如,鼠标或触控板)、磁盘驱动器单元、信号生成设备(例如,扬声器)、触摸屏、加速度计、麦克风、蜂窝射频天线(cellular radio frequency antenna)和网络接口设备,其可以是,例如,网络接口卡(NIC)、Wi-Fi卡或蜂窝调制解调器。
图8中描绘了本发明的示例性系统501。计算机901包括耦合到处理器309的输入/输出设备305和有形的、非瞬时存储器307。在某些实施方案中,计算机901可以通过网络517与服务器511通信。服务器511还可包括耦合到处理器309的I/O设备305和存储器307。服务器可以存储一个或更多个数据库385,其能够存储在如上所述的本发明的方法中有用的记录399。
本发明的方面包括算法和实现协议,如本文所述。SENTRYSEQ技术是基于以下洞察:高准确度PCR酶比下一代测序机更不容易出错:如果目标是高保真测序,则创建每个单独的分子的多个拷贝,对这些分别测序,且然后创建反映原始分子的序列的共有序列,并对在测序过程中产生的(大多数)错误取平均值是个好主意。
本主题方法的方面包括鉴定可能包含突变的(低频)等位基因的BAM比对文件的列。系综一致性检查的概念可以应用于通过寻找系综链平衡对于相容序列的一致性来检查从SENTRYSEQ文库构建的de Bruijn图中鉴定的推定变异。
系综
系综是共享相同的起始和终止坐标的对齐的读段对的集合(精确定义:对于每个读段对,存在读段对的碱基与其对齐的参考基因组坐标的一组坐标;每个这样的集具有最大值和最小值;系综是具有相同最大值和相同最小值的读段对的集)。
忽略测序误差,单个系综包含源自参考基因组中具有相同的起始/终止坐标的原始分子的PCR产物的读段。重要的是,原始分子的两条链应该被系综的成员代表,并且可以通过检查两个来源链是形成系综的“左侧”(意思是:较低的参考坐标)的第一读段还是第二读段(在Illumina配对末端范例中)来区分它们。
过度扩增因子
过度扩增因子是来源于每个原始分子的读段的平均数;如果测序和PCR是完全的并且所有原始分子是独特的,则每个系综的读段数将等于过度扩增因子。
尽管可以通过实验测量过度扩增因子,在本范例中,从输入BAM文件统计地估计过度扩增因子。估计程序是基于以下洞察:大多数原始分子是独特的,并且大多数系综因此应包含与过度扩增因子类似的数目的读段(即,过度扩增因子的第一近似值可以通过确定直方图的模式来计算,该直方图绘制x轴上每个系综的读段数与y轴上具有该读段数的系综的数量)。
有效系综
真实的测序数据包含测序错误,且并非所有读段可被完全映射。为了将一组读段对对齐转换为系综列表(其中每个系综包含一组读段对对齐),使用上面给出的定义:具有相同的最大/最小坐标的所有读段对对齐变为同一系综的一部分。重要的是,该定义是基于完整对对齐的最大值/最小值,而不是基于2个单独读段的最大值/最小值(即,忽略2个单独读段对齐的“内”端)。
在读段对齐的开始和结束处的测序误差(在对齐的坐标中,对应于由机器产生的两个单独成员读段的开始)将导致错误的读段形成它们自己的系综。另外,只有满足一系列一致性标准的读段对被考虑。这些可包括:
●两个成员读段被映射并在同一染色体上。
●读段1和读段2的相反链。
●两个对齐的最小值和最大值之间的总距离<常数(即假定的潜在分子的长度<常数;常数通常为~330)。
●由对准器(BWA)设定的所有质量控制标志对于两个读段对齐OK;有一个“QC”标志和一个“正确对”标志。
●最小映射质量(当前>0.95)。
●跨读段的错配的比例<常数(单独测量的;当前常数=2,即,标签不活跃)。
●所有读段成员都没有软剪切或填充(当前不活跃)。
●基因组一致性:基于原始分子来自相对正常的人类基因组的假设,要求(读段对的)读段与参考基因组的正链对齐在另一个读段的“左侧”(通过其各自对齐的最小坐标测量)-且反之亦然。
区分来自正和负链的系综成员
通过检查系综的“左侧”读段(如上定义)是读段对的第一读段还是第二读段,可以区分来自原始分子系综成员的两条链中的哪一条。
重要的技术考虑
如上定义的标准要求对的两个读段都具有连续的对齐;由支持分段读段对齐的对齐算法产生的BAM,如BWA-mem,是有问题的。
分析过程的概述
提供有输入BAM后,SENTRYSEQ进行以下步骤:
1.系综枚举
●寻找存在于BAM中的所有系综,并且将它们的坐标(和协变量如长度、GC含量、成员读段数量…)写入文本文件clusters.txt。
●输出文件后,从工作存储器中删除所有系综数据。
●主要函数:clusterGenerator/clusterGenerator.cpp中的clusterGenerator::enumerateClustersInBAM(..)。
2.过度扩增的统计学估计
●考虑到协变量如GC含量、系综长度、与下拉探针的重叠,调用读取clusters.txt并估计过度扩增的统计模型的R脚本。还估计了沿输入分子长度的分布和输入分子基因组覆盖。
●主要文件:R/analyeSENTRY.R。
3.确定性分析:
●迭代BAM文件的所有列,并且鉴定可能包含突变的等位基因的那些列。列中的每个等位基因是簇的成员,且将等位基因按簇成员资格和它们来自原始分子的哪条链分组。鉴定具有可能的突变的列的阈值考虑了来自统计学过度扩增模型的估计值。
●主要函数:analysis/deterministic/deterministicAnalysis.cpp中的deterministicAnalysis::kickOff()。
4.概率分析(不一直活跃)
●对于候选列或对于所有列,应用完整的扩增模型,其明确考虑扩增错误的不同情景(在PCR的不同循环)并且将它们的可能性与突变的输入等位基因的不同的情景进行比较。
●主要函数:analysis/probabilistic/probabilisticAnalysis.cpp中的probabilisticAnalysis::kickOff()。
本发明的方面包括在以下进一步详细描述的高保真测序方法和方案。分子群体中低丰度体细胞突变(突变等位基因比例<5000-1)的检测和准确频率估计的必要条件是,在整个样品制备和文库制备过程中保持衍生的等位基因Nd(对应于体细胞变体)与祖先等位基因Na(对应于种系基因组)和来自其他来源的DNA N的比。
衍生的等位基因的比例f可以由于(a)在测序文库构建过程中经由损失消耗Nd,或(b)通过污染增加分母而减少。
对于来自无细胞DNA(cfDNA)样品的测序循环肿瘤DNA(ctDNA)测序的应用,必须采取措施以控制(a)通过使血液抽取期间和/或之后由凋亡细胞释放的细胞核DNA污染最小化,和控制(b)必须采取措施,以是文库制备期间分子的损失最小化。
检测低频等位基因(f<〖10〗^(-3))中的一个挑战是,高通量测序具有约O(1个错误/1000个碱基)的测序错误率。存在Illumina测序错误的已知协变量,例如,读段中的位置、碱基、均聚物长度等。为了控制错误率,生成原始分子的PCR重复,且然后使用统计学模型来评估在所鉴定的重复(这被称为系综)上聚集的每个检测的变体处真实变异与错误的证据。通过扫描共享对齐和读段长度来从头构建系综,以鉴定由潜在PCR重复引起的读段,说明在原始群体中可以存在多个相同分子的事实(相同原始分子的数量是cfDNA浓度和cfDNA长度分布的函数)。每个原始分子的平均重复数被称为过度扩增因子。
通过使用统计学模型传播覆盖潜在候选变体的序列读段中的不确定性并说明推断的潜在分子数来最小化过度扩增因子。与其他方法相比,这具有减少所需测序(主要成本组分)的影响。因此,本文描述的文库制备方案已经与用于鉴定变体及其相关统计显著性的统计模型联合优化。
本发明的方面包括用于从无细胞DNA(cfDNA)制备用在Illumina测序平台上的测序文库的方法,除文库制备外,该方法可应用于任何鸟枪测序仪上的任何片段化的DNA。例如,这意味着,通过使DNA片段化(使用例如限制酶或超声处理)且然后应用相同的系综生成策略,可以在细胞群体中检测少数细胞群体。
图2示出了Illumina衔接子连接产物。方案修改导致衔接子堆叠。这样做是为了使测序相容产物的数量最大化(参见图3,用于堆叠的衔接子的PCR分辨)。图3显示通过引物结合竞争对堆叠的衔接子和所得PCR产物的分辨。如果最内侧的引物在最外侧的PCR引物退火位点之前或同时结合,则结果是从PCR产物中消除最外侧的引物。由于最内侧首先结合的等待时间是几何分布的,因此在4轮PCR之后,不获得与测序相容的产物的机会为仅1/16。
图4-5示出了来自肺癌患者的cfDNA文库的实例。使用该方法观察到大概加倍的可测序产物。在图4中,观察到四个峰,前三个峰涉及平均分子长度加上2、3和4个衔接子。在PCR(图5)之后,模式转换到平均分子长度加上2个测序衔接子。还观察到两个较长的片段群体。
本发明的方面包括用于使用杂交捕获技术将DNA样品分离成两个群体分区的方法。杂交捕获是一种基于DNA分子的核苷酸序列从群体中分离特定DNA分子的方法。为了捕获分子,将双链DNA解链成单链DNA(例如,通过升高温度),然后加入杂交捕获探针(探针),并改变条件以促进链退火。探针与靶序列互补并具有使得分子能够被分离的选择标志物(例如,生物素)。为了防止探针池中杂交探针之间的杂交,池中的所有探针都被设计为与靶基因座的有义序列或反义序列互补。因此,每个双链DNA分子仅捕获一条链。
通常,在杂交捕获之前对样品DNA进行PCR扩增,这导致原始分子的两条链在有义和反义PCR重复群体中被代表。为清楚起见,考虑以下toy实例:x={x+,x-}是双链分子,α和β是长度为l的单链DNA分子,α的前n个连续碱基与β的最后n个连续碱基共有互补序列,其余序列是非互补的。因此,与β退火的α具有来自互补序列的双链DNA茎、和来自非互补序列的单链DNA臂的叉状Y形结构。
现在,使用已知的连接反应,例如钝端连接产生分子,其中x的任一端侧翼为Y形{α,β}双链DNA:
αx-β
βx+α
然后使用与α和β互补的引物(分别表示为α_c和β_c)应用PCR,生成PCR重复家族:
将具有探针序列〖x〗_(-)的杂交捕获方案应用于PCR产物,将提取具有序列α_c x_(+)β_c和β〖x〗_+α的分子,其每个分别是反义和有义链的后代。随后的PCR生成所有四个单链分子,因此两条链被代表。然而,使用这种方法,一半的分子不被捕获。目前捕获方法仅使用一组有义探针或反义探针。
可以使用链特异性分离来从原始取样的DNA生成两种相同分布的样品。这对于寻求在异质群体中检测低频率的分子的应用作为控制在取样的DNA的后续操作中引入的误差和丢失(dropout)的手段是有用的。提议以下两步骤过程:
为感兴趣的基因座分别制造两个杂合体捕获组;一组有义,一组反义。有义、反义的应用顺序无关紧要,因此称为组A和B。然后将组连续地应用于DNA样品,如下所示。亲和力选择继续:将可选择的探针应用于单链DNA,使用标准的杂合体捕获方案将样品分离成分离物(被探针结合的DNA)分区和非分离物(未被探针结合的DNA)分区。
步骤1:将A应用于DNA群体。靶序列将被收集在分离物分区中。保留非分离物分区。
步骤2:将B应用于非分离物分区。将在步骤2的分离物分区中收集组1靶序列的互补序列。
可存在来自A的探针的一些携带污染,但如果优化了分离方法,这将是最小的。在另一种实施中,样品可以分成两个等分试样,且A和B分别应用,从而避免在前一步骤中由探针携带引起的任何交叉杂交。
分别分析来自A和B的分离物,然后寻找两个实验之间结果的一致性,这控制了在样品的下游处理中引入的伪相。这为分离物A和B之间的复制提供了机会,并通过分别评估A和B提高了灵敏度。
本发明的方面包括用于进行杂合体捕获区域选择程序的方法。靶向的高通量测序通过减少评估个体中指定基因座所需的测序读段总数而被鼓励。所需读段的减少是靶向序列长度除以基因组长度的商的函数,并且权重由靶向的和全基因组测序的分布测序读段覆盖深度(此后缩写为覆盖)确定。
增加的覆盖提高了灵敏度,因为含有靶等位基因的读段数目与真实变体比例(1-ε)×f和覆盖D大致呈二项式分布,其中ε是测序中的碱基错误率且f是分子群体中等位基因的频率。覆盖与测序错误之间的关系是复杂的,但假设无系统错误,覆盖可以通过使得跨越靶基因座的读段信息聚合(整合出错误)来减少假阳性。由于测序中存在系统错误模式,例如均聚物中的错误,需要更复杂的错误模型。
选择基因组的哪些区域来靶向是设计靶向测序组中的重要考虑因素。在使用遗传变体标志进行癌症检测的背景下,靶向组的统计功效是患者群体内跨这些基因座的变体重复出现的函数。在杂合体捕获设计中的另一个考虑因素是,每个杂交探针的特异性和跨所有探针的灵敏度的均匀性,两者都驱动在所需检测限值检测变体所需的测序读段的量。
这里描述的方法是选择多达总序列长度L的基因座的组合,其针对癌症患者中的最大组合重复出现负荷优化(组合驱动物和乘客遗传变体),说明影响杂合体捕获性能的决定因素,例如序列独特性和GC含量。此外,设计了与cfDNA长度分布匹配并跨越跨目标区域上观察到的GC含量分布的合成的DNA添加物。基于指定的参考错配,添加物可以与cfDNA区分开来,选择错配的模式使得它们不可能从自然过程中观察到。这些添加物用于计算跨GC环境以及预测的杂合体捕获重叠的假阴性率估计。
模型概述
模型鉴定周期性体细胞突变(局部扩增、易位、倒位、单核苷酸变体、插入、缺失)的区域和预先指定的基因座(如癌基因外显子),并选择提供最多信息的区域组合,直到指定的总组大小。
●指定基因组长度、所考虑的基因组改变并强制包含指定基因。
●指定所考虑的肿瘤变异数据库和肿瘤类型、以及每个数据库的相对权重。
●指定每种肿瘤类型的群体发病率是否正确(以防止采样偏差)。
●指定区域是否应对外显子或在基因组水平生成。(这些区域已针对参考基因组中的独特性进行了校正)。
图6提供了杂合体捕获组设计过程的示意表示,包括数据转换。鼓形表示数据库,虚线框表示输入,菱形表示操作,且实心边框表示输出。
区域选择优化
为了鉴定基因组中提供最多信息的区域以靶向指定的总组长度,评估记录在癌症突变数据库例如COSMIC中的遗传改变的标志。优化使用前后优化或Greedy优化来完成。
然后使用交叉验证程序来验证设计,以说明通过从有限数量的样品构建组引起的潜在偏倚。在设计癌症组时,交叉验证策略是重要的,因为样品中的遗传变异在肿瘤内(肿瘤内异质性)和患者之间(肿瘤间异质性)是异质的,并且受因素诸如遗传背景(例如,POLE突变状态)、环境暴露(例如,抽烟史、先前疗法)和肿瘤阶段影响。因此,潜在群体的结构可以影响组设计,交叉验证是一种防范这种结构的众所周知的策略。
前向后向
通过在前向和后向传递(forward and backward passes)之间交替来鉴定基因座,直到从L基因座构建指定长度的组。基因座分层为组中包含的那些(选择的位点)、以及未包括在组中的那些(可用的基因座)。
前向迭代传递:鉴定可用基因座中向组添加最大数量的体细胞突变的基因座,f*。
后向:将f*包括到组集中,然后鉴定所包括的基因座中添加最少的体细胞重复出现的基因座,b*。
如果f*不等于b*,则排除b*。开始下一个迭代。
该方案鉴定用于组合的体细胞重复出现的基因座的最佳集。当达到组长度时,优化结束。
Greedy优化
从添加最大体细胞突变负荷的基因座开始,将其添加到组,然后从剩余的基因座选择具有最大体细胞突变负荷的基因座。当组合的序列满足指定的组大小时,终止。该算法不保证鉴定全局最优。
交叉折叠验证
交叉折叠验证用于评估所鉴定的组的稳定性,说明疾病数据库中结构的影响。
构建患者样品的两个相互排斥的集,其中由训练比例p确定集的基数。在第一集上生成具有基数p的组,基数p记录组上具有突变的患者的总数。在具有基数(1-p)的验证集中验证所提议的组,计算组上带有突变的患者的比例。如果患者比例在阈值T内,保留该组。否则修改。
数据库询问
为了获得遗传变异的样品,询问了肿瘤活组织检查测序的数据库,它由许多患者的协变量如疾病类型、阶段、环境暴露、组织学等分层。为了防止假阳性变体在癌症数据库中出现,将健康群体的群体测序诸如1000Genomes数据库中观察到的所有种系遗传变体去除,假阳性变体会使组设计(鉴定在癌症中、而不是在健康个体中发现的突变为目标)混乱。存在已知的种系突变,诸如使个体易患癌症的BRCA1/2突变,其可通过该方法消除。然而,允许感兴趣的区域被强制进入设计以掺入此信息。
数据转换
为了说明杂交捕获中的差异性能,将关于人类基因组的序列特性的信息掺入到组选择过程中。
具体地,掺入关于基因组中每个碱基的独特性的量度,因为这驱动了杂合体捕获的特异性。例如,如果基因座与人类基因组中的99个其他基因座(例如LINE元件)同源(相同),则捕获探针将仅下拉每100个中平均1个相关基因座。使用的量度为1)。
通过使用可从UCSC基因组浏览器数据库获得的基因组独特性的两个预先计算的汇总统计来掺入该信息。
1.可映射性,s,其定量kmer序列与基因组比对的独特性
2.独特性,u,在1个碱基滑动窗口中跨基因组的35个碱基窗口的独特性
组合这两种映射,并且然后为人类参考基因组中的每个碱基生成特征编码的独特性值(character encoded uniqueness)。因此,将参考基因组从核苷酸序列转换为由杂交特异性得分f(s,u)注释的核苷酸序列。
软件描述
●createEndcodeReferenceGenomes.pl
输入:
人类参考基因组的BED文件(基于零的坐标)。
wgEncodeDukeMapabilityUniqueness35bp.bed—kmer序列如何独特地与参考基因组对齐,s,其中s=1/((“基因组中的匹配数”)),例如,对于基因组中的一个匹配,s=1,对于两个匹配,s=0.5。
wgEncodeCrgMapabilityAlign36mer.bed—每个序列在特定碱基处开始和为特定长度(此处为36)的正链上如何独特,u,其中对于>=4个匹配u=0,对于4个匹配u=0.25,对于3个匹配u=0.33,对于2个匹配u=0.5,和对于独特匹配u=1。
输出:
FASTA格式文件*.refGen。基因组中的每个碱基编码有编码根据以下的参考基因组独特性/可映射性的字符
"chr"(65+"int"(20*V))
其中V是如输入中描述的s或u。
●explore[MUSIC,COSMIC]samplesIDs.pl
输入:
*Uniq*.refGen from createEncodedReferenceGenomes.pl
*Map*.refGen from createEncodedReferenceGenomes.pl
变异数据库
TCGA somatic_mafs_cleaned_filtered/*_cleaned_filtered.maf
COSMIC…
PARAMS:
Ignore1000G<BOOL>排除在1000G中观察到的变体
checkMappabilty<BOOL>
mappabilityThreshold<DOUBLE>用于碱基的阈值>=阈值
输出:
Exons.txt<gene-exon#,length[bp],gene,exon,chromosome,start,end>
Bins.txt<chromosome-start-stop,chromosome,start bp>
Mutations_inBins.txt<TCGA_tumour-v-TCGA_normal,chromosome-start-stop,mutation count>
Mutations.txt<TCGA_tumour-v-TCGA_normal,gene-exon#,count>
Kernel.txt<chromosome,postion,mutation count,mutation count*prevalence of disease>
Samples.txt<TCGA_tumour-v-TCGA_normal,disease type,mutationcount>
allPositions_preQC.txt
操作:
加载1000G数据并从分析中排除所有位点
排除具有reference Uniqueness(参考独特性)<=mappabilityThreshold(可映射性阈值)的所有基因座
排除具有reference Mappability(参考可映射性)<=referenceMappability(可映射性阈值)的所有基因座
以200个碱基的间隔跨基因组生成箱
排除具有>=1000个基因组变体的所有样品
排除具有>=肿瘤位点的所有样品
对于COSMIC:排除所有TCGA样品,保留基因组范围的、非编码的插入。
计算每个样品id、每个外显子和每个箱的突变。
本发明的方面包括用于估计测序错误以校正变体频率估计的方法。已经观察到,循环肿瘤DNA(ctDNA)比例(fraction)与肿瘤大小、阶段、治疗响应和预后相关。成像的肿瘤大小用于追踪治疗响应和缓解。已经表明,追踪ctDNA变体与成像的肿瘤直径具有高相关性(>90%,Pearson相关性(其他研究已显示使用追踪肿瘤鉴定的突变的类似结果)。因此,从ctDNA准确估计体细胞突变具有为患者的临床决策提供依据的潜力。
使用患者中的多个核苷酸体细胞变体
在这里,描述了使用样品中体细胞改变估计测序错误和非均匀覆盖对变体等位基因频率估计的影响的方法。为此,鉴定与种系基因组具有不同的N个连续碱基(N>1)的体细胞变体,使其由矢量V={a(1),a(2),…,a(n-1),a(n)}表示,其中元素a(i)代表变体中在位置i的不同碱基。这种变体可以通过体细胞改变生成:易位、倒位、插入、缺失、扩增或突变。在一些实施方案中,一种或更多种考虑的碱基不需要含有体细胞改变,条件是这种考虑的碱基彼此足够接近(例如,在彼此约1、2、3、4、5、6、7、8、9、10、11、12、13、14或15个碱基内)。
对于变体中的每个碱基a(i),计数支持该碱基的等位基因的总数,这生成了对V的频率f(V)的n个估计。所有观察到的频率f(a(i))应该等于f(V),但是由于覆盖的变化和测序错误,情况可能并非如此。然后,可以使用已知的统计学方法来定量在测序期间产生的频率估计中的离差。然后,这可用于校正频率估计。一个例子将是使用样品平均值和方差来使用适当的采样分布来估计置信区间。
使用患者中的杂合种系变体
在二倍体生物中,等位基因在杂合位点处的比应为1/2。在人类群体中存在SNP分离的大数据库。对于给定的个体,可以询问这些位点,并将杂合位点鉴定为具有大致相等的等位基因频率的两个等位基因的基因座。然后可以从在杂合位点处观察到的第二等位基因的频率构建等位基因频率的经验分布。如果杂合位点的数量足够大,可构建每个等位基因组合(A>C、A>G、…、T>G)的频率估计。然后可以使用分布来校正在体细胞变体位点处的频率估计。
掺入DNA对照
将具有与患者不同的序列的已知输入量的DNA添加到样品中。这些是样品中变体等位基因的阳性对照。
为了生成可鉴定的添加物,生成在人类群体中不可能观察到的序列。这通过以下来完成:1)选择在群体测序数据库中具有低报告多样性的区域,2)向序列引入不反映自然突变过程的变化(例如序列(相同)n,{变化,相同,变化,相同,变化},(相同)n)。对照序列被进一步区分,因为添加物的长度(120个碱基)是已知的,且所引入的变化的位置也是已知的。
已知,杂合体捕获可受捕获探针和靶DNA之间的错配数量的影响。将四个突变引入每个对照。构建添加物以便1)GC含量和2)探针-靶重叠的影响可以通过以下观察:1)选择跨靶向区域具有与已知GC含量分布不同的GC百分比的序列,和2)改变120碱基长的对照DNA与其相应的下拉探针的重叠百分比。
在血液抽取之前将添加物加入到血液收集真空采血管中,以便a)样品可以从它们的测序中鉴定,允许在测序中鉴定样品混合,b)以便可以估计来自有核白细胞凋亡的污染(这在本文进一步描述),并且c)以便可以检测到假阴性。
本发明的方面包括用于检测人类中的无细胞循环DNA的污染的方法。图7提供了根据本发明的实施方案的一种方法的示意视图。
来自人血浆的无细胞循环DNA(pDNA)除了多数比例的来自人的正常(通常是健康的)基因组的分子外,还包含癌症患者中的肿瘤DNA片段和孕妇中胎儿DNA的片段。调查肿瘤或胎儿DNA的混合部分本质上具有挑战性,因为癌症/胎儿衍生分子的混合比例可以低至5000个分子中之1。
任何给定的未加工的血液样品(通常但不总是储存在EDTA管或不同类型的血液收集容器中)将含有一定比例的无细胞DNA以及白细胞和红细胞(WBC和RBC)。经过一段时间后(并受环境因素诸如温度的影响),所含的WBC将经历细胞死亡并开始将所含的DNA片段释放到循环中。由于该过程,血液样品中包含的任何肿瘤或胎儿衍生的无细胞DNA将进一步稀释,使得它们的检测和表征甚至更具挑战性。
存在技术解决方案(例如Streck管)以防止所含的WBC破裂并释放其DNA,但是这些解决方案并不完美并且稀释问题仍然存在,特别是如果血液储存较长时间段或运输。
对于基于调查肿瘤或胎儿衍生DNA的存在或特征的任何诊断方法,因此期望测量和控制潜在的污染。
潜在使用情形包括:
1.建立样品特异性的检测下限,
2.对样品进行质量控制,包括拒绝带有过多污染量的样品或对具有过多污染物的样品不产生诊断读出,
3.通过估计的污染物分子数重新缩放检测到的异常/肿瘤衍生/胎儿衍生分子的数量,
4.使用重新缩放的检测分子数作为比检测到的分子的原始数量将允许的疾病状态或进展的更准确表示。
5.估计稀释的绝对量。
在一些实施方案中,方法包括以下步骤中的一个或更多个:
1.鉴定人类基因组中的段或区域,即:a)在人类群体的绝大多数中为纯合的(或在期望的靶群体的绝大多数中为纯合的)和b)基因组复杂性高,即使用用于读段对齐的标准算法方法对从该区域衍生的分子建立基因组起源是明确且无挑战的。通常,该段的长度在50和150个碱基之间变化,但是本文描述的方法可使用更长和更短的区域。
2.通过用不同核苷酸取代许多核苷酸或引入或缺失许多核苷酸来扰乱段或区域的序列。通常,该步骤将包括用不同的核苷酸取代位于序列中央的一个或两个核苷酸。
3.证实,如此扰乱的序列不存在于正常人类群体中。存在多种标准方法来实现这一点,例如基因组比对或与由群体测序数据,诸如1000Genomes Project生成的de Bruijn图比较。如果此验证失败,需要重复步骤2或1。
4.使用DNA合成方法生物化学合成(近似或精确地)n个拷贝的如此扰乱的序列。通常选择数量n,使得当n个分子被引入步骤6中抽取的血液时,n与所抽取的血液体积中人类基因组的预期拷贝数之间的比类似于肿瘤/胎儿衍生片段与正常基因组片段之间预期的/所需的最小比。(例如,如果人们预计1000个循环片段中的1个是肿瘤起源的,并且如果1ml血液通常包含约1000个拷贝的人类基因组,并且人们抽取5ml血液,则每个管n=5将是明智的选择)。
5.以下步骤之一:
a.产生包含如此扰乱的序列的n个合成拷贝的血液收集容器。这可以是标准的真空管,专门设计用于防止WBC破裂的容器,或任何其他种类的血液收集容器。
b.产生添加剂组分,包含如此扰乱的序列的n个合成拷贝。添加剂组分在与血液接触后会溶解,且然后释放扰乱的序列的n个拷贝。
6.在时间X用人类血液填充所用的血液收集容器。如果使用步骤5.b,在加入血液之后立即或之前加入添加剂组分。通常,血液容器现在将被运输到处理设施,例如使用快递服务。
7.通过离心提取无细胞循环DNA并从提取的DNA制备DNA文库。
8.使用将在样品的下游解释中使用的技术测量扰乱序列的频率(fP)和非扰乱序列的频率(fn)。通常,使用基于数字PCR的方法或基于测序的方法来测量频率,使用全基因组测序方法或靶向测序方法。
9.观察到的频率的解释:
a.fP/(fP+fn)是样品中以n拷贝原始地(即,在由于WBC破裂的稀释开始之前)存在的肿瘤或胎儿来源的等位基因的稀释后频率的估计。
b.依据所采用的下游解释技术的特性,如果fP/(fP+fn)为0或小,则应拒绝或不解释样品。
c.将稀释后数据中假定的肿瘤或胎儿衍生等位基因的观察频率乘以([(fP+fn)/fP]x n将得出该等位基因的稀释前绝对计数的估计值。肿瘤等位基因计数及其随时间的发展已被证明是疾病状态和进展的重要指标。
对不同的基因组基因座和不同的n值进行上述程序赋予了重要的另外优点,包括但不限于以下:
●利用不同的基因组基因座将减少统计变异性,并且如果基因座选自一系列局部序列背景,则可用于控制GC含量偏倚。
●使用不同的n值将使得能够(更准确地)估计稀释总量(在稀释衍生的分子片段中测量)并因此血液样品中DNA片段的稀释前数目。
为了估计由c表示的污染分子(即,从血液样品中的凋亡有核细胞释放的那些DNA分子)的数量,可以使用两步取样方法。注意,c随着时间单调增加。扰乱序列(如上鉴定和合成的)被称为基准序列。将基因组中在该位置的取样的pDNA分子数用d表示。
在血液采集时,b1立即添加到样品中(实际上,收集容器甚至可以预先用第一基准序列分子接种,参见上文)。因此,第一基准的频率是
f1=b1/(d+c(t=0))。
然后将样品运输到收集设施。在时间T,在从样品分离pDNA之前,进行基准序列频率的第二次测量。观察样品频率f(1)和f(2),然后计算观察频率的差异以确定污染分子的数量。
通过引用并入
在本公开内容中,已经对其他文献,诸如专利、专利申请、专利出版物、期刊、书籍、论文、网页内容进行了参考和引用。所有此类文件为了所有目的通过引用整体并入本文。
等同物
除了本文所示和所述的那些之外,本发明的各种修改及其许多其他实施方案对于本领域技术人员来说将从本文件的全部内容(包括对本文引用的科学和专利文献的参考文献)变得明显。本文的主题包含重要的信息、示例和指导,其可以适应于本发明在其各种实施方案及其等同物中的实践。

Claims (24)

1.一种用于核酸测序的方法,所述方法包括:
获得样品中的核酸的多个测序读段;
鉴定包含具有共享起始坐标和读段长度的两个或更多个测序读段的系综;
确定所述样品中存在的对应于所述系综测序读段的原始输入分子的数量;
鉴定所述系综中的候选变体;和
使用概率模型和所确定的原始输入分子的数量来确定所述候选变体是真实变体的可能性。
2.根据权利要求1所述的方法,其中获得测序读段包括:
从所述样品制备测序文库;
扩增所述测序文库;和
使用下一代测序(NGS)将所述测序文库测序。
3.根据权利要求2所述的方法,其中所述制备步骤包括使用约16小时的反应时间、在约16摄氏度的温度将衔接子连接至所述核酸。
4.根据权利要求2所述的方法,其中所述扩增步骤包括PCR扩增,并且所述方法还包括使用计算机模型选择所需的过度扩增因子和PCR循环数目以检测样品中的指定浓度的变体。
5.根据权利要求2所述的方法,还包括:
基于包括鸟嘌呤-胞嘧啶(GC)含量、靶群体中的突变频率和序列独特性的因素设计靶向基因组区域的杂合体捕获组,和
在所述测序步骤之前使用所述杂合体捕获组捕获扩增的核酸。
6.根据权利要求5所述的方法,其中所述捕获步骤包括使用靶向靶基因座的有义链的第一杂合体捕获组和靶向所述靶基因座的反义链的第二杂合体捕获组。
7.根据权利要求2所述的方法,还包括在扩增所述测序文库之前将合成的核酸对照添加至所述样品,并且使用所述合成的核酸对照的测序读段确定错误率。
8.根据权利要求7所述的方法,其中所述合成的核酸对照包含已知序列,所述已知序列在所述核酸所来源的物种中具有低多样性,并且具有与所述已知序列的多于一个非天然存在的错配。
9.根据权利要求8所述的方法,其中所述多于一个非天然存在的错配是4个。
10.根据权利要求7所述的方法,其中所述合成的核酸对照包括代表所述杂合体捕获组的靶基因座的鸟嘌呤-胞嘧啶(GC)含量分布。
11.根据权利要求7所述的方法,其中所述合成的核酸对照包括包含与所述杂合体捕获组的下拉探针的不同重叠的多于一种核酸。
12.根据权利要求7所述的方法,还包括使用所述合成的核酸对照的测序读段确定错误率。
13.根据权利要求7所述的方法,还包括确定候选变体频率。
14.根据权利要求1所述的方法,其中所述核酸包括无细胞核酸。
15.根据权利要求2所述的方法,其中所述样品包括组织样品,其中获得测序读段还包括在所述制备步骤之前将所述核酸片段化。
16.根据权利要求15所述的方法,其中所述片段化步骤包括超声处理或酶促裂解。
17.根据权利要求1所述的方法,还包括在应用所述概率模型之前向所述候选变体应用确定性模型。
18.根据权利要求17所述的方法,其中所述确定性模型包括当在所述核酸的有义链和反义链两者上均未鉴定出候选变体时丢弃所述候选变体。
19.根据权利要求1所述的方法,其中所述概率模型是似然估计模型。
20.一种用于鉴定核酸变体的系统,所述系统包括耦合到存储指令的有形的、非瞬时存储器的处理器,所述指令在由所述处理器执行时使所述系统:
鉴定包含来自样品的核酸的两个或更多个测序读段的系综,所述测序读段具有共享起始坐标和读段长度;
确定所述样品中存在的对应于所述系综测序读段的原始输入分子的数量;
鉴定所述系综中的候选变体;和
使用概率模型和所确定的原始输入分子的数量来确定所述候选变体是真实变体的可能性。
21.根据权利要求20所述的方法,所述方法是进一步可操作的,以在应用所述概率模型之前向所述候选变体应用确定性模型。
22.根据权利要求21所述的系统,其中所述确定性模型包括当在所述核酸的有义链和反义链两者上均未鉴定出候选变体时丢弃所述候选变体。
23.根据权利要求20所述的系统,所述系统是进一步可操作的,以基于包括鸟嘌呤-胞嘧啶(GC)含量、靶群体中的突变频率和序列独特性的因素对所述两种或更多种测序读段确定靶基因组区域。
24.根据权利要求20所述的系统,其中所述概率模型是似然估计模型。
CN201780007584.7A 2016-01-22 2017-01-20 用于高保真测序的方法和系统 Pending CN108603229A (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
US201662286110P 2016-01-22 2016-01-22
US62/286,110 2016-01-22
PCT/US2017/014426 WO2017127741A1 (en) 2016-01-22 2017-01-20 Methods and systems for high fidelity sequencing

Publications (1)

Publication Number Publication Date
CN108603229A true CN108603229A (zh) 2018-09-28

Family

ID=59362079

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201780007584.7A Pending CN108603229A (zh) 2016-01-22 2017-01-20 用于高保真测序的方法和系统

Country Status (4)

Country Link
US (1) US20190338349A1 (zh)
EP (1) EP3405573A4 (zh)
CN (1) CN108603229A (zh)
WO (1) WO2017127741A1 (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113628683A (zh) * 2021-08-24 2021-11-09 慧算医疗科技(上海)有限公司 一种高通量测序突变检测方法、设备、装置及可读存储介质

Families Citing this family (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110016499B (zh) 2011-04-15 2023-11-14 约翰·霍普金斯大学 安全测序系统
CN104956225B (zh) 2012-10-29 2018-10-02 约翰·霍普金斯大学 卵巢和子宫内膜癌的帕帕尼科拉乌测试
US11286531B2 (en) 2015-08-11 2022-03-29 The Johns Hopkins University Assaying ovarian cyst fluid
US11332784B2 (en) 2015-12-08 2022-05-17 Twinstrand Biosciences, Inc. Adapters, methods, and compositions for duplex sequencing
KR102326769B1 (ko) * 2016-03-25 2021-11-17 카리우스, 인코포레이티드 합성 핵산 스파이크-인
MX2020001575A (es) 2017-08-07 2020-11-18 Univ Johns Hopkins Materiales y métodos para evaluar y tratar el cáncer.
EP3775198A4 (en) * 2018-04-02 2022-01-05 Grail, Inc. METHYLATION MARKERS AND TARGETED METHYLATION PROBE PANELS
CN109097458A (zh) * 2018-09-12 2018-12-28 山东省农作物种质资源中心 基于ngs读段搜索实现序列延伸的虚拟pcr方法
EP3856903A4 (en) 2018-09-27 2022-07-27 Grail, LLC METHYLATION MARKER AND TARGETED METHYLATION PROBE PANEL
WO2020264565A1 (en) * 2019-06-25 2020-12-30 Board Of Regents, The University Of Texas System Methods for duplex sequencing of cell-free dna and applications thereof

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150044687A1 (en) * 2012-03-20 2015-02-12 University Of Washington Through Its Center For Commercialization Methods of lowering the error rate of massively parallel dna sequencing using duplex consensus sequencing
WO2015083004A1 (en) * 2013-12-02 2015-06-11 Population Genetics Technologies Ltd. Method for evaluating minority variants in a sample
US20150324519A1 (en) * 2014-05-12 2015-11-12 Roche Molecular System, Inc. Rare variant calls in ultra-deep sequencing
US20150368708A1 (en) * 2012-09-04 2015-12-24 Gaurdant Health, Inc. Systems and methods to detect rare mutations and copy number variation

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6312892B1 (en) * 1996-07-19 2001-11-06 Cornell Research Foundation, Inc. High fidelity detection of nucleic acid differences by ligase detection reaction
EP1112378A1 (en) * 1998-07-17 2001-07-04 GeneTag Technology, Inc. Methods for detecting and mapping genes, mutations and variant polynucleotide sequences
US8055034B2 (en) * 2006-09-13 2011-11-08 Fluidigm Corporation Methods and systems for image processing of microfluidic devices
WO2011143231A2 (en) * 2010-05-10 2011-11-17 The Broad Institute High throughput paired-end sequencing of large-insert clone libraries
US20130173177A1 (en) * 2010-08-24 2013-07-04 Mayo Foundation For Medical Education And Research Nucleic acid sequence analysis
US20160040229A1 (en) * 2013-08-16 2016-02-11 Guardant Health, Inc. Systems and methods to detect rare mutations and copy number variation
EP3122894A4 (en) * 2014-03-28 2017-11-08 GE Healthcare Bio-Sciences Corp. Accurate detection of rare genetic variants in next generation sequencing

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150044687A1 (en) * 2012-03-20 2015-02-12 University Of Washington Through Its Center For Commercialization Methods of lowering the error rate of massively parallel dna sequencing using duplex consensus sequencing
US20150368708A1 (en) * 2012-09-04 2015-12-24 Gaurdant Health, Inc. Systems and methods to detect rare mutations and copy number variation
WO2015083004A1 (en) * 2013-12-02 2015-06-11 Population Genetics Technologies Ltd. Method for evaluating minority variants in a sample
US20150324519A1 (en) * 2014-05-12 2015-11-12 Roche Molecular System, Inc. Rare variant calls in ultra-deep sequencing

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
JUSTIN M. ZOOK等: "Synthetic Spike-in Standards Improve Run-Specific Systematic Error Analysis for DNA and RNA Sequencing", 《PLOS ONE》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113628683A (zh) * 2021-08-24 2021-11-09 慧算医疗科技(上海)有限公司 一种高通量测序突变检测方法、设备、装置及可读存储介质
CN113628683B (zh) * 2021-08-24 2024-04-09 慧算医疗科技(上海)有限公司 一种高通量测序突变检测方法、设备、装置及可读存储介质

Also Published As

Publication number Publication date
WO2017127741A1 (en) 2017-07-27
EP3405573A1 (en) 2018-11-28
US20190338349A1 (en) 2019-11-07
EP3405573A4 (en) 2019-09-18

Similar Documents

Publication Publication Date Title
CN108603229A (zh) 用于高保真测序的方法和系统
US11932910B2 (en) Combinatorial DNA screening
US11367508B2 (en) Systems and methods for detecting cellular pathway dysregulation in cancer specimens
CN110800063B (zh) 使用无细胞dna片段大小检测肿瘤相关变体
CN102171565B (zh) 等位基因调用和倍性调用的方法
KR102028375B1 (ko) 희귀 돌연변이 및 카피수 변이를 검출하기 위한 시스템 및 방법
CN110770838B (zh) 用于确定体细胞突变克隆性的方法和系统
EP3924502A1 (en) An integrated machine-learning framework to estimate homologous recombination deficiency
EP3571615B1 (en) Methods for non-invasive assessment of genetic alterations
AU2020398913A1 (en) Systems and methods for predicting homologous recombination deficiency status of a specimen
US20210130900A1 (en) Multiplexed parallel analysis of targeted genomic regions for non-invasive prenatal testing
US20210102262A1 (en) Systems and methods for diagnosing a disease condition using on-target and off-target sequencing data
CN112218957A (zh) 用于确定在无细胞核酸中的肿瘤分数的系统及方法
CN104346539A (zh) 从目标测序面板中寻找变异的方法
WO2019025004A1 (en) METHOD FOR NON-INVASIVE PRENATAL DETECTION OF FETUS SEX CHROMOSOMAL ABNORMALITY AND FETUS SEX DETERMINATION FOR SINGLE PREGNANCY AND GEEMELLAR PREGNANCY
CN114207727A (zh) 用于从变体识别数据确定起源细胞的系统和方法
US11869630B2 (en) Screening system and method for determining a presence and an assessment score of cell-free DNA fragments
Mayrink et al. A Bayesian hidden Markov mixture model to detect overexpressed chromosome regions
WO2023031485A1 (en) Method for the diagnosis and/or classification of a disease in a subject
Wilson Jr Statistical Methods for the Estimation of Cell-type Composition and Cell-type Specific Association Studies
Gevaert A Bayesian network integration framework for modeling biomedical data

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
TA01 Transfer of patent application right
TA01 Transfer of patent application right

Effective date of registration: 20220126

Address after: California, USA

Applicant after: Greer Co.,Ltd.

Address before: California, USA

Applicant before: Grail, Inc.