CN113053460A - 用于基因组和基因分析的系统和方法 - Google Patents

Abstract

本发明涉及一种对人类核酸样品进行基因组和/或基因分析的方法，所述方法包括以下步骤：提供一组人类参照基因组；测试所述人类核酸样品的性别和/或血统；基于所述性别和/或血统测试的结果，从所述一组人类参照基因组中选择一个或多个种群特异性的人类参照基因组(population‑specific human reference genomes,PHREG)；将所述人类核酸样品与选择的PHREG进行比对；以及相对于选择的PHREG进行变体识别。本发明还提供了相应的计算机系统和计算机程序。

Description

用于基因组和基因分析的系统和方法

本发明涉及一种用于对人类核酸样品进行基因组和基因分析的系统和方法。

背景技术

下一代测序(NGS)

下一代测序，又名高通量测序，是本领域技术人员所熟知的一种对核酸片段进行高通量并行测序的常规方法。可从各种供应商处商购获得(参见，例如，www.illumina.com)下一代测序的设备和方法学。

下一代测序是用来描述一些不同的现代测序技术的总称，包括：

·Illumina(Solexa)测序；

·Ion torrent：质子/PGM测序；

·SOLiD测序。

NGS技术产生高质量的DNA序列(“读长”)。这些读长与基于毛细管的Sanger测序技术所产生的读长(650-1000bp)相比大大缩短。Sanger测序技术由Frederick Sanger及其同事于1977年开发，是近30年来最为广泛使用的测序方法。Sanger的读长以低通量高成本的形式产生，而NGS方法则以适中的成本产生短得多的读长(25-500碱基)，但一次NGS运行进行测序的碱基对的总数要高出几个数量级。这两个因素造成许多新的信息学挑战，包括处理这些数百万或甚至数十亿的短NGS读长的能力。通常以下列两种方式之一处理这些读长：将它们映射回其在已有骨架/参照序列中的正确位置，建立一条与骨架序列相似但不一定相同的序列(称为“读长映射(Read Mapping)”)，或者将它们构建成一条新的序列(称作“从头组装(De Novo Assembly)”)。

相对于从头组装而言，读长映射至参照基因组的主要优点在于，其大大简化了基因组推断的过程。从头组装需要去发现整条基因组序列且会产生许多模糊性，而基于参照序列的再测序只需要去发现样品与所述参照序列的差别。就复杂性和时间要求而言，从头组装比映射组装要慢几个数量级且占用更多内存。

读长映射是NGS分析管线中的第一步且是最基本的一步，旨在去发现新测序的人类基因组(或其片段，如一小部分靶向基因或外显子组)相对于先前所测的人类参照基因组的变体。

也可以用读长映射来比对这些数百万或数十亿的短NGS读长，从而检测覆盖度(在特定位置/位点的读长数目)，这是NGS实验以及所有将来所得结论的关键质量参数。

人类参照基因组(HRG)

2001年2月，美国联邦政府的人类基因组计划(Human Genome Project)和私人公司Celera Genomics一起成功完成了整个人类基因组的草图，并随后进行了数次修改[Lander et al.,2001；Venter et al.,2001；Church et al.,2001]。多年来，基因组组装已稳步改进，并发布了新版本(“build”)，因此当前的基因组参照合作体(GenomeReference Consortium,GRC)人类基因组组装GRCh38[Scheneider et al.,2017]可以说是目前存在的最佳组装的哺乳动物基因组，其只剩875个组装间隙和少于1.6亿个未指定的“N核苷酸(截止至GRCh38.p8)，而第一个版本具有约150,000个间隙[Editorial(October2010)."E Pluribus Unum".Nature Methods.7(5):331.Doi:10.1038/nmeth0510-331]。

HRG是现今人类遗传学和基因组学中使用的最重要的单一资源，因此，其作为通用的坐标系统，也是描述了注释(基因、启动子等)和基因变体的空间[Harrow et al.,2012；ENCODE,2012；1000 Genomes Project Consortium,2012]。HRG还可以作为下一代测序分析管线中读长比对步骤中的参照。在该映射的下游，它可用于功能测定和变体识别管线[LiH&Durbin 2009；DePristo et al.,2011]。

HRG的初始版本由13名主要有欧洲血统的匿名DNA捐献者在纽约水牛城志愿捐献的一小群DNA序列组成[Snyder et al.]。1997年3月23日，星期天，通过在水牛城新闻上的广告招募捐献者。前十名男性和前十名女性志愿者受邀与该计划的基因咨询师预约并捐献了可从中提取DNA的血液。这些DNA样品的处理方式导致大约80％的参照基因组来自八个人，而一位标记为RP11的男性占了总共的66％。

为了识别和解决较大的组装问题，例如，包含大规模重复和结构变体的复杂区域、由新的基因组映射技术测得的序列数据和源自新供体的单倍型单一资源已被引入到最新的内部版本中。在提交本申请时，GRCh38包含了来自大约50位不同个体的序列，参见http://www.bio-itworld.com/2013/4/22/church-on-reference-genomes-past-present-future.html.

HRG的局限性

1.HRG是线性的

人类DNA被装在称为染色体的物理上独立的单元中。人类是二倍体生物，包含两组遗传信息，一组遗传自母亲，另一组遗传自父亲。因此，每个体细胞具有22对常染色体(每对中的成员各来自一个亲本)和2条性染色体(男性中的X和Y染色体以及女性中的两条X染色体)。每条染色体包含一条非常长的线性DNA分子。最小的人类染色体中的DNA分子由大约5千万对核苷酸组成；最大的染色体包含大约2.5亿对核苷酸。

由上可知，二倍体人类基因组由24种不同类型的46个单个DNA分子组成。因为人类染色体成对存在，且每对中的两个DNA分子几乎完全相同，故获取完整的代表性人类基因组的信息只需对30亿对核苷酸(单倍体基因组)进行测序。所以常说人类基因组包含30亿对核苷酸，尽管大多数人类细胞包含60亿对核苷酸。单倍型人类基因组由22条常染色体以及Y和X染色体组成。

每条染色体代表一个单条DNA分子，即数百万个核苷酸碱基的序列。这些分子是线性的，故可预期每条染色体应该由一单条的、连续/线性的核酸序列来表示。然而事实并非如此，主要有两个原因：1)由于基因组DNA的性质和测序方法的局限性，基因组中的有些部分还没有被测序，及2)基因组的某些区域在不同个体之间的变化很大，以至于它们不能用单条连续序列来表示。而HRG由24条由正常碱基(A、C、T或G)组成的线性DNA序列来表示，其中由一系列的“N”来表示间隙，以清楚表明该组装内间隙的位置。

人类基因组计划的主要目的是，即便存在不确定的区域，也要为每条物理染色体生成单条代表性序列，即单条“框架(scaffold)”。尽管它也包括了少数表示等位基因变体的替换框架(存在于SNP位点的DNA碱基的不同版本称为等位基因)，但它们与主要的框架没有正式关系(formalized relationship)。由于认识到基因组的某些高度多态区域特别难以用单条参照序列表示，所以从GRCh37[Church et al.,2011]开始添加了正式模型来引入具有代表性的高度可变区域替换版本。相对于“主要”(单倍体)组装，锚定于沿着主要框架的位置，描述了从千碱基(kilobase)至数兆碱基(multi-megabase)的“替代位点框架”形式的序列。在提交本申请时的组装(GRCh38.p9)包括178个区域和总共261条线性序列[Patenet al.,2011]。

另一个复杂因素是HRG是由最初国际基因组测序计划的多名匿名个体的一群DNA推导而得。因此，最终得到的HRG确实是一个随机混合的复合体(conglomerate)、不同DNA序列的单倍型拼合(mosaic)，以至于在一些情况下，可能无法以单条线性序列来进行正确地表示。

2.HRG显然不是无疾病的

Chen&Butte(2011)在HRG中识别了3,556个疾病易感变体，其中15个为罕见变体(主要等位基因频率<1％)。作者们使用精选高质的定量人类疾病-SNP关联数据库，评估了参照基因组针对104种疾病相对于健康人群风险增加的可能性比率，并发现1型糖尿病、高血压和其他疾病具有高风险。这清楚地证明了HRG不能代表一个普通人，显然也不是没有疾病的。尽管HRG极大地加速了个人基因组测序的分析工作，但仅关注相对于该参照基因组的变体有可能会错过许多致病变体，包括罕见变体[Chen&Buttle，2011]。

3.参照等位基因偏向欧洲血统

现有技术的NGS分析管线中，使用HRG组装的主要问题是，其实际上来自相对较少的匿名捐献者的DNA样品，而这些匿名捐献者偏向欧洲血统，因此只代表小部分人类基因变体的样品。

尽管对大部分基因组而言，以该参照基因组作为坐标系具有相对的有效性和通用性，但人们越来越担心，用HRG作为镜头(lens)来研究所有其他人类基因组，会排除大量常见的人类变体并引入普遍的参照等位基因偏倚[Petrovski et al.,2016，Paten et al.,2017]。参照等位基因偏倚倾向于过度报告参照基因组中存在的等位基因，而过低报告基础DNA与参照等位基因不匹配的其他等位基因[Degner et al.,2009,Brandt et al.,2015]。

这种偏倚主要产生于测序实验中的读长映射和比对步骤。为了正确地映射，产生读长的基因组序列既要在参照序列中有表示，又要与参照序列足够相似以使其能被识别为相同的基因组元素。当这些条件都不能满足时，映射的错误会对真实序列产生系统盲目性[Paten et al.,2017]。基于参照基因组在每个位点偏倚的血统历史，参照等位基因偏倚还可能对某些遗传亚群和基因组某些区域的影响要比其他基因更大[Petrovski et al.,2016,Paten et al.,2017]。高度多态的区域(如HLA基因)特别容易受到参照等位基因偏倚的影响[Nielsen et al.,2011]，特别是在用单个参照基因组作为NGS读长比对的指标时。在这种情况下，很多真实的变体无法被识别，因为单倍体中的这些真实的变体不同于作为指标的基因组，所以无法比对从这些区域产生的读长，导致信息丢失[Brandt et al.,2015]。

由上所述，参照偏倚是人类基因组测序中用HRG来进行变体识别的已知问题，对参照的修改可以提升识别准确度和可译性[Fakhro et al.2016]。缓解此问题的一种方法是，在基因组-解释过程的早期通过修改参照基因组来修改变体普遍性信息，以使得基因组中发现的变体是该种群中的次要等位基因[Dewey et al.,2011]。对参照的修改可简化分析工作流程，因为它可以减少假阳性的数量且需要解释的变体更少[Fakhro et al.,2016]。

未来：基于图形的参照结构/基因组图

人们越来越认识到，单一的单倍型参照基因组对于人类遗传学和基因组学来说是一个不良的通用参照结构，因为它仅代表了人类变体的一小部分：还有相对于参照基因组无法轻易描述的变体和注释[Horton et al.,2008,Pei et al.,2012]。而且，如上所述，作为读长映射和解释的靶标，它引入了参照等位基因偏倚。为了缓解这些问题，参照基因组组装的最新版本(例如，在提交本申请时的人类基因组组装(GRCh38.p9))已包括“替代位点”序列(alternate locus sequence,“alts”)：被认为是高度多态化的人类基因组区域的其他多种序列代表，其末端固定在“主要”(单倍型)参照组装内的位置。这种包含多个部分重叠的序列路径的结构，可以被认为是数学图形的一种形式——基因组图[Novak et al.,2017]。

图形在生物序列分析中具有很长的历史，它们经常被用来在其中紧凑地表示可能的序列集合。通常，序列本身被隐式编码为图中的游动。这使得图形很自然地适于表示参照集合，因为它们本质上是相关序列的集合[Paten et al.,2017]。该图形不仅包含样品的近似序列，也包含了很多它们特定的变体。

基因组图有望改进读长映射、变体识别(variant calling)和单倍型确定。可以预期，基于图形的参照将取代人类和其他具有被测序个体集合的应用中的线性参照[Novaket al.,2017]。很多项目正在构建和应用这些基因组图。现在可以从常见变体库构建基因组图，尽管仍处于实验阶段，但一些工具说明了基于图形的方法有巨大潜力。

尽管具有这些理论上的优势，但使用基因组图谱进行变体识别的研究仍然相对较新。还有很多问题有待解决。如何表示重复(duplication)和重复(repeat)？如何最好地映射至图谱？应该如何解析同源性不清楚的短变异体？如何使用图谱来实现更全面的变体分类法？这些问题都为将来的研究开辟了道路。

为了在实践中有用，基因组图谱必须能够将其允诺的参照偏倚的减小转化为在变体识别中与已有方法相比具有可测量的改进。相应地，为基因组图谱开发变体识别算法是当前一个重要的研究前沿。

卡塔尔基因组(QTRG)

卡塔尔是波斯湾上的一个小半岛，总人口约30万卡塔尔人。卡塔尔是世界上近亲结婚率最高的国家之一，并且仍在增加，而卡塔尔的内婚率接近100％。所有这些因素以及庞大的家庭规模是造成本土遗传疾病发病率高的主要原因，这意味着对卡塔尔预算的财政负担。这些因素触发了卡塔尔政府寻求方法来保护其国民以免受遗传疾病的威胁[Zayed2016]。

2013年，政府官员决定启动卡塔尔基因组计划(Qatar Genome Project,QGP)(http://www.gulf-times.com/story/374345/Qatarlaunches-genome-project)。该计划的目的在于对每个卡塔尔公民的基因组进行测序，以通过绘制致病变体/罕见变体的图谱并建立一个卡塔尔参照基因组作为通向个性化医疗的路径，来保护卡塔尔人免受本土遗传病高发病率的侵害。该计划的最终目标是将信息应用于临床实践，并使这种方法成为卡塔尔医疗体系的常规组成部分[Zayed 2016]。为了实现QGP预期的临床应用前景，必须解决几个严峻的挑战，包括达到较高的变体识别敏感度和准确度[Koboldt 2010]。

为了促进中东和北非地区精准医疗的发展，通过整合来自1,161名卡塔尔人(代表了0.4％的人群)的全基因组测序的等位基因频率数据，构建了种群特异性的基因组，该基因组专门针对卡塔尔本土阿拉伯人群的疾病研究而定制。在卡塔尔总共观察到2090万个单核苷酸多态性(single-nucleotide polymorphism,SNP)和310万个插入和删除(insertions and deletion,inDel)，每个基因组平均包含1.79％新的变体[Fakhro etal.,2016]。

1000基因组计划(1kG)

1000基因组计划(1000 Genomes Project)成立于2008年，旨在对全球至少1000人的基因组进行人类遗传变体(相对于HRG GRCh37)和单倍型的测序并生成目录(因此被称为1000基因组计划)。该计划当前的第3阶段的分析包含来自26个种群的2,504名个体，并定义了5个所谓的超级种群，每个超级种群由4至7个种群的联合体组成[1000 Genomes ProjectConsortium et al.,2015]。这种更精细的单倍型资源将有助于在基因组水平和地理水平上了解遗传变异[Baye,2011]。

发明目的

NGS技术的最新进展使得DNA和RNA测序变得更为快速且价格低廉，从而彻底改变了基因组学和分子生物学的研究。针对健康人群和疾病人群的基因组测序计划已经检测出许多功能或疾病相关的基因组变体，这些变体可以为新型临床应用提供治疗靶标或基因组标记的线索。

基因的变体识别主要基于将原始序列的读长相对于参照基因组进行比对(读长映射)。这种基于比对的方法有许多局限性，包括基因组组装的不完整性[Meyer,L.R et al.,2013]、正常个体基因组中存在的结构变体[Sudmant et al.,2015]、读长中的测序误差以及单核苷酸多态性(SNP)对读长映射的干扰[Iqbal,Z.et al.,2012]。

目前，在提交本申请时，相对于线性HRG的读长映射是唯一的标准方法，也将会在临床NGS分析管线和人类个体测序中持续作为标准，因为HRG作为多数基因组的坐标系具有相对的有效性和普遍性。而且(不同于用基因组图谱进行的新兴基因组推断)，已经公开了很多成功使用线性参照基因组进行变体识别的方法学[Nielsen et al.,2011]。

然而，如上所述，HRG的主要问题是其偏倚，忽视了关于物种内遗传变异的先前信息。目前，通常通过修改参照基因组来解决该问题，这样，相对于修改后的参照基因组识别到的变体是种群中的少数等位基因。

使用NGS技术的临床基因组学的成功需要对个人基因组变体进行准确一致的识别。这些目标的前提是精确的读长映射(比对)以及随后的变体识别。

本发明的一个目的是检测新的生物标记物，特别是基因变体，例如单核苷酸变体(SNV)、插入和删除(inDel)、拷贝数目变体(CNV)和结构变体(SV)(例如，染色体易位(translocation)、倒位(inversion)、重复(duplication)、大的插入和删除(inDel))，以便在人类基因组研究中使用下一代测序技术。

另一个目的是提高已有的基于NGS的生物标记物的准确性和可信度(例如，用于癌症治疗，其中该技术用于分析肿瘤细胞及其损坏的DNA)。

发明内容

根据本发明的第一方面，本发明提供了一种对人类核酸样品进行基因组和/或基因分析的方法，包括以下步骤：

a)提供一组人类参照基因组；

b)测试所述人类核酸样品的性别和/或血统(ancestry)；

c)基于步骤b)中的性别和/或血统测试的结果，从所述一组人类参照基因组中选择一个或多个种群特异性的人类参照基因组(PHREG)；以及

d)将所述人类核酸样品与步骤c)中选择的PHREG进行比对。

下文将种群特异性的人类参照基因组(Population-specific Human ReferenceGenomes,PHREG)理解为血统特异性的参照基因组和性别特异性的参照基因组。PHREG极大限度地减小了参照偏倚并提高了比对准确度，而且，如果之后进行变体识别的话，还能提高变体识别的准确度。有利地，本发明不仅提高了比对的精确度，还提高了计算速度、正确比对的读长数目和比对的计算步骤数目。使用PHREG对人类核酸样品进行基因组和/或基因分析还可以提高读长覆盖深度，并可以用变体识别的敏感度改善来评估使用PHREG的好处。

在本发明的内容中，术语“人类核酸样品”通常是指从人类样品中分离出来的任何核酸样品。特别地，人类核酸样品可能包括NGS读长，其在下文中有更详细的定义。

人类核酸样品通常可以包含来自各种适于准备人类核酸的标准生物化学、分子和/或细胞生物学操作的样品。这样的操作包括穿刺、活检、液体活检、无细胞DNA分离试剂盒等。人类核酸样品可以是或衍生自各种合适的来源，包括但不限于，体液、粘膜、组织、组织提取物或细胞或其任何组合。人类核酸样品也可以是来自各种合适来源的对照样品。人类核酸样品可以包含例如血液样品、血浆样品、尿液样品、肿瘤样品，其可能包括组织处理程序FFPE(formalin-fixed paraffin-embedded tissue or formaldehyde-fixedparaffin-embedded tissue)中的固定过程导致的不希望的错误(undesired artefact)。

特别地，人类核酸样品可能包含DNA、RNA和/或大小分割(size fractionated)的总DNA或RNA。提供来自目标样品的DNA可能包括一个或多个生物化学纯化步骤，例如，离心、裂解和/或分层(fractionation)，通过机械或化学破坏步骤进行的细胞裂解，包括但不限于多次的冰冻和/或解冻循环、盐处理、苯酚-氯仿萃取、十二烷基硫酸钠(sodium dodecylsulfate,SDS)处理和蛋白酶K消化。任选地，提供来自目标样品的DNA可进一步包括通过在存在聚乙烯或盐的情况下进行沉淀来除去大型RNA(例如，丰富的核糖体rRNA)的步骤，或通过在存在盐(优选地，氯化钾溶液)的情况下进行沉淀来除去干扰性的十二烷基硫酸钠(SDS)的步骤。从细胞和/或组织中纯化总DNA或RNA的方法对本领域人员而言是熟知的，所述方法包括，例如，标准操作，例如使用硫氰酸胍-酸性苯酚-氯仿萃取(例如，

Invitrogen,美国)。但是，同样优选的是，还可以不以任何在此描述的生物化学沉淀和/或纯化步骤的形式提供目标样品中的DNA。

在本发明的内容中，术语“核酸”通常指任何种类的单链或双链寡核苷酸分子，其由脱氧核糖核苷酸或核糖核苷酸或两者组成，包括基因组DNA、核DNA、体细胞DNA、种系DNA、合成设计的和/或制备的DNA，合成设计的和/或制备的DNA包括但不限于，来自信使RNA谱的体外产生的DNA，优选的是cDNA形式。术语“核酸”通常指具有相同或相似长度的单链或双链寡核苷酸分子，例如，由相同或相似数目的核苷酸组成。

人类核酸样品可以包含基因组序列，其可用于在基因组、转录或转录后的水平评估、分析、比对、指标化和/或绘制特定突变。因此，根据本发明的人类核酸样品可以指以及包括但不限于任何编码区域、非编码区域、外显子、内含子、染色体和/或染色体内区域、启动子区域、增强子区域、编码小和/或长调控RNA的区域、活性转录区和/或非转录区、转座子、热点突变区域、移码突变的区域等。

“一组人类参照基因组”包含至少两个人类参照基因组，优选地，多个人类参照基因组。步骤b)中的性别和/或血统测试允许在步骤c)中从一组人类参照基因组中选择一个或多个最适合的人类参照基因组。在优选的情况下，步骤b)中的性别和/或血统测试使得性别和/或血统会自动分类，并允许从一组人类参照基因组中选择一个PHREG以进行后续步骤d)中的比对，但是为后续的分析选择一个或多个额外的PHREG也是可能的。

优选地，步骤b)中的性别和/或血统测试是基于血统和/或性别特异性的序列变体的子类，所述序列变体与从专业数据库(curated database)中提取出来的性别和/或血统相关。优选地，这些序列变体是单核苷酸多态(SNP)和/或单核苷酸变体(SNV)。这种用于性别和/或血统测试的序列变体的子类也被称为种群依赖性人类血统和性别模式(Population dependent Human Ancestry and Sex Patterns，PHASPs)。优选地，该专业数据库包含所有种群中所有已知的序列变体。PHASP数据集是该专业数据库的摘录，其数据集远小于PHREG数据集，且为分类中最有区别的子集。用于生成PHASP的技术是来自机器学习(包括特征缩减，其中所述特征是基因型)的计算机方法。可以用标准分类的结果来比较和测试这些学习。

优选地，性别和/或血统测试包括一个初始的比对步骤，来检测所述样品的单个序列变体模式，其中将所述人类核酸样品和一个单人类参照基因组(例如，GRCh37或GRCh38)进行比对。这种用于步骤b)中测试的单人类参照基因组不是血统或性别特异性的。通过将样品的所述序列变体模式与PHASP数据集进行比对，就确定了患者的血统和性别。

根据一个实施方案，所述测试可以包括性别测试。根据另一个实施方案，所述测试可以包括血统测试。根据另一个实施方案，所述测试可以包括性别测试和血统测试。

在一个示例性实施方案中，一组人类参照基因组包括男性和女性参照基因组。如果步骤b)中的性别测试确定了所述人类核酸样品是男性或女性参照基因组，那么在步骤c)中，各自的男性或女性参照基因组或基因组们将被选为代表性的PHREG来进行后续的比对步骤d)。

由于性别染色体包含同源性序列，用一个性别调整过的参照基因组(有X和Y染色体是男性，没有Y染色体是女性)可以防止读长的错误比对。因此，使用性别特异性的参照基因组减少了后续变体识别的假阳性和假阴性。

在另一个示例性实施方案中，一组人类参照基因组包含一些血统特异性的参照基因组。步骤b)中的血统测试确定了一些血统特异性的参照基因组中最好的一个或一些参照基因组。然后在步骤c)中，最接近的一个或一些参照基因组将被选为PHREG或PHREGs来进行后续的比对步骤d)。

选择了错误的血统会导致大量的假阳性和假阴性的变体识别。用血统特异性的参照基因组可以有效地增加正确比对的读长的数目，减少假阳性和假阴性。

同样地，当一组人类参照基因组包含了血统特异性的男性参照基因组和血统特异性的女性参照基因组时，性别和血统测试的组合是具有决定性的。

步骤b)中的术语“测试”应被理解为包括人类核酸样品的至少一种基因和/或基因组测试。基因和/或基因组测试比任何来自“自我报告(self-reporting)”的信息更为可靠。自我报告的和研究者指派的血统通常依赖于对一些遗传或非遗传信息复杂组合的主观解释，这些信息包括行为、文化、社会准则、肤色和其他影响因素。研究参与者或患者无差错地报告其种族是很罕见的。自我报告种族的差错可以由许多原因造成；有些人可能不能完全意识到他们真实的血统或只知道最近的血统(或他们的地理上的来源)，而其他人可能只认识到一种种族，尽管他们具有混合的背景[Mersha&Abebe 2015]。文献证实了自我报告的血统和性别经常是不正确的[Ainsworth,2015；Mersha&Abebe,2015]。事实上，Ainsworth甚至解释了100个人中有一个人会被性别发育紊乱症所影响，导致其外观与其基因组不符。

有利地，该方法还可以根据性别和血统来检测样品是否被互换，以此作为额外的质检。自我声明和测序运行预测的性别和血统之间的不匹配可以揭示，例如，样品移位和实验室处理中的其他错误。

术语“比对”通常指一种计算机步骤，其中被测序的样品与参照序列进行比较和拟合。为此，需要为生成的测序数据中的每条读长找到其在该序列中的相应部分。换言之，比对或读长映射是为观察到的核酸测序读长决定其在基因组序列中最有可能的来源的过程。在典型的实施方案中，读长是NGS读长，但应该要明白的是，来自其他测序技术的读长也包括在本发明的教导中。

来自人类核酸样品的比对的读长可以被显示、储存、打印、通过通信网络发送、或其它进一步处理。特别地，对所比对的人类核酸样品的进一步的应用或使用可以包括下列中的一种或多种：

1)在插入和删除(inDel)周围进行局部的再比对

术语“inDel”表示在基因组中插入或删除碱基对，通常包括长度为1至1000bp的小基因变体。在inDel周围的再比对改善了后续的数据分析，特别是后续的变体识别。

2)碱基质量分数再校准(Base quality score recalibration,BQSR)

术语“碱基质量分数”描述了每个碱基的误差估计，其代表了测序仪器识别碱基的可信度。该分数可以用来例如对后续的变体识别的证据进行加权。BQSR允许通过考虑系统技术误差来调整质量分数，该系统技术误差由如何进行测序的物理学或化学造成。

3)通过机器学习从而将真实的分离的变体与下一代测序技术中常见的机器假象(machine artifact)分开。

4)变体发现和基因分型从而发现所有潜在的变体，在此也被称为变体识别。

变体发现可包括发现SNP/SNV、InDel、CNV和SV(染色体易位、倒位、复制、大型inDel)。

5)进化分析研究

进化分析研究可包含测量核苷酸多样性、种群差异、连锁不平衡以及一个或多个种群突变频谱的工具。进化分析通常可以包括用于计算进化序列统计的计算工具。所述计算工具可以适于在跨染色体或框架的滑动窗口中执行分析。所述计算工具可以例如产生样品的系统发育树。

可以用例如“POPBAM”软件进行这种进化分析，具体描述可见，例如，https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3767577/.

6)测试野生型生物标记物

进一步，可以测试比对的人类基因组样品中是否存在野生型生物标记物，即，那些由于被包括在PHREG内而不会在变体识别中被检测到的生物标记物。因此，比对之后的计算步骤可以包括测试每个已知的生物标记物，该测试表明该生物标记物是否存在于比对的人类基因组样品中，不管PHREG在这个位置的信息是什么。

根据一个实施方案，所述方法包括对与选择的PHREG比对的人类核酸样品进行变体识别的额外步骤。有利地，本发明通过引入初始性别和/或血统测试来确定正确的PHREG，用于后续比对和变体识别步骤，从而提高了变体识别的准确度。

因此，比对的人类核酸样品，更具体地，比对的来自人类核酸样品的NGS读长还可以由一种或多种所谓的变体识别器进行进一步处理，这些变体识别器是包括不同的变体识别算法的计算机模块，所述算法可以检测任何种类的变体(例如，SNV、InDel、拷贝数变化和结构变体)。后续的方法步骤可以包括变体解释。变体识别和/或变体解释的结果可以被显示、储存、打印、通过通信网络发送或其它方式进行进一步处理。有利地，通过从使用的参照基因组中去除偏倚，该方法能够检测到之前(例如在癌症或其它疾病中)没有发现的生物标记物。特别地，根据本发明的方法允许区分多种基因突变，包括但不限于SNV、多核苷酸变体(MNV)、复杂事件(complex event)和大结构变体，特别是热点突变、移框突变(frame-shiftmutation)、非沉默突变、终止密码子突变、核苷酸插入、核苷酸缺失、拷贝数变异(copynumber variation)、拷贝数改变(copy number alteration)和/或剪接位点。

人类核酸样品的供体可以是患者，即，患有疾病或怀疑患有疾病的人。但该方法的应用不应理解为仅限于患者。

变体识别或解释可以包含指示是否存在某种疾病的基因组序列分析。基于变体解释，可以将患者分为包括未指示某种治疗的第一组患者，和包括指示了某种治疗的第二组患者。因此，本发明可以有利地用作疾病筛查程序的一部分，以评估患者中疾病的存在与否。

另外或或者，该方法可以包括检索得到(retrieve)与人类核酸样品有关或相关的疾病指征的步骤。该疾病指征可以从例如电子健康记录中检索得到，或者由患者本人或主治医生通过计算机设备的输入方法来手动加入。可以根据疾病分类数据库(diseaseontology)来识别疾病指征，例如ICD-10、MeSH或MedDRA。对于某些类别的指征，可能还会有专门的分类数据库，它们具有诸如更精确的指征分类的优势。在肿瘤学中，使用ICD-O-3和/或TNM分期系统可能是有益的。

基于变体识别和解释的结果，以及考虑到患者的疾病，该方法可能涉及为患者提供一种治疗计划。在本文中，治疗计划可以特别是针对患者的个性化治疗计划，其中所述个性化治疗计划包括针对患者的遗传数据，特别是针对他/她的临床、分子和/或遗传状况而定制的治疗选项。

为了为患者确定有希望的治疗方法，该方法可以包括检查是否有任何变体(例如，在患者中(例如在患者的肿瘤或者正常对照组织中)发现的突变)能指示患者在任何治疗下的结果。该方法还可以包括确定所有与发现的任何变体相关的治疗。该方法可以包括为确定的治疗方法打分并根据分数进行排名，以便为患者提供一种疗法选项或治疗禁忌的优先次序。

在本发明的内容中，术语“治疗”包括以预防、改善或治疗伴随指征的症状的药物组合物形式给予治疗有效的药物或药物活性化合物。术语“治疗”还包括任何种类的手术、放射疗法和/或化学疗法或其任何组合。

对于这两种选择，即在筛查方法或个性化治疗计划的背景下，本发明可以使医生的诊断能力得到提高，例如，由于比对和变体识别的准确度的提升，使得医生的治疗决策得到改善。

根据一个实施方案，在多数等位基因水平上与PHREG进行比对。该多数等位基因水平上使用PHREG中独特的核苷酸代码(A、C、G、T)将参照序列调整至一个种群。选择在该种群中特定位点最常观察到的单个核苷酸。在等位基因频率相关的情况下，可以使用存在于基础参考序列(例如，GRCh37或GRCh38)中的等位基因。

根据另一个实施方案，在非罕见等位基因水平上与PHREG进行比对。该非罕见等位基因水平使用了根据已建立的IUPAC命名法的模糊核苷酸代码(ambiguity nucleotidecodes)[Cornish-Bowden,1985]，例如，“R”代表“A”或“G”。非罕见等位基因水平可编码最多两个或三个，优选地，两个在种群中频率很高的等位基因。频率很高可以定义为超过或等于30％、20％、15％、10％、5％、3％、1％或0.1％，特别地超过或等于5％。由于PHREG的每个基因组位置整合了一个以上的变体等位基因，因此可以预期读长比对会更精确。在一个实施方案中，在非罕见等位基因水平上只考虑了单核苷酸变体(SNV)。在其他实施方案中，也考虑了插入和删除(inDel)以及其他结构变体。

根据一个实施方案，在多数等位基因水平上相对于PHREG进行变体识别。在一些实施方案中，在非罕见等位基因水平上进行比对并在多数等位基因水平上进行变体识别。或者，在非罕见等位基因水平上进行变体识别。

根据一个实施方案，步骤a)中的人类参照基因组是公开的人类参照基因组。特别地，公开的人类参照基因组可以包括HRG的内部版本，特别地GRCh37和GRCh38的内部版本。另外或或者，公开的人类参照基因组可以包括QTRG。另外或或者，公开的人类参照基因组可以包括来自1000基因组(1kG)计划的基因组。对于1kG计划，可以从1kG FTP站点ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/下载和使用最新版本的所有染色体的VCF文件。如果有更多个体和种族的数据集被公开(例如，研究阿拉伯联合酋长国人群的1000个阿拉伯基因组计划[Al-Ali,M.et al.,2018])，这些也可被用于本发明的方法。

另外或或者，步骤a)中的人类参照基因组衍生自公开的人类参照基因组。术语“衍生自”可能特别地包含错误校正和/或将人类参照基因组调整到多数等位基因编码水平或非罕见等位基因水平。

可以执行错误校正，使得在给定种群的零个个体中观察到的参考核苷酸被相应的多数核苷酸取代。

在一个实施方案中，步骤a)包括将人类参照基因组调整至编码水平，所述编码水平包含独特的核苷酸代码或模糊的核苷酸代码。特别地，包含独特的核苷酸代码的所述编码水平可以用于在多数等位基因水平上定义PHREG。特别地，包含模糊的核苷酸代码的所述编码水平可以用于在非罕见等位基因水平上定义PHREG。

在一个实施方案中，对于调整至编码水平，考虑了单核苷酸变体。对每个种群(或超级种群)，使用了所有报道的SNV及其等位基因频率。在其他实施方案中，也考虑inDel、CNV和/或SV。

根据一个实施方案，提出了将参照序列调整至种群的四种不同的水平，其中两种水平限制在独特的核苷酸代码(A、C、G、T)，另外两种水平利用了根据IUPAC命名的模糊核苷酸编码[Cornish-Bowden,1985]，例如，“R”代表“A”或“G”。这些PHREG编码水平的定义如下：

1.最大程度的保守错误校正(Maximally conservative error correction)：在种群的零个个体中观察到的参照核苷酸被相应的多数核苷酸(例如，相应的多数1kG核苷酸)取代。

2.多数等位基因：选择在种群给定位点最常观察到的单核苷酸(如果与等位基因频率有关联，则使用基础参照序列(例如GRCh37或GRCh38)中存在的等位基因)。

3.非罕见等位基因：编码种群中频率很高(例如>＝5％)的至多两个等位基因，必要的时候可以使用IUPAC代码。

4.对观察到的等位基因进行完整建模：在每个位置编码所有(至多四个)等位基因，其至少在种群中的一位个体中有报道。

然而，第4水平PHREG中1kG变体的完整表示是以不成比例的大量基因组修饰为代价，这些修饰引入了模糊性(ambiguity)，从而可能严重阻碍读长映射器对种子的发现。因此，在一个实施方案中，用第3水平进行比对，其使用了IUPAC模糊感知比对算法(IUPACambiguity-aware alignment algorithm)。由于目前表现最佳的变体识别器不能用来处理模糊代码，因此后续的变体识别使用第2水平PHREG，除非有更好的IUPAC模糊感知比对算法。

因此有利地，该方法可以根据目标种群和下游的分析，允许PHREG对种群遗传变异进行用户定义水平的调整。

根据一个实施方案，步骤a)中的人类参照基因组是PHREG。因此，步骤a)可以包括，例如，从公共资源下载PHREG。

如上述定义，首先，PHREG被理解为血统特异性的参照基因组和/或性别特异性的参照基因组。在一个实施方案中，步骤a)中提供的人类参照基因组已经是种群特异性的了，因为他们包含指示其血统和/或性别的元数据(meta data)。例如，在提交本发明时，目前1kG计划第3阶段的分析包含来自26个种群的2,504位个体和5个所谓的超级种群，每个超级种群为由4至7个种群组成的联合体。可以在http://www.internationalgenome.org/faq/which-populations-are-part-your-study找到来自1kG研究第3阶段的26个种群和它们相关的5个超级种群(AFR，非洲人；AMR，Ad混血美洲人；EAS，东亚人；EUR，欧洲人；SAS，南亚人)。

在一个实施方案中，使用来自1kG计划的数据，为31个(超级)种群中的每一个以及另一个包含所有其他种群的超级种群构建优化的种群特异性的人类参照基因组。

当步骤a)中提供的人类参照基因组是PHREG时，可能也提供了(例如，通过从公共资源下载)PHREG的公开元数据。该元数据可能为该方法提供质量控制。如果元数据与性别和血统分类数据重合，则可以认为质量控制是成功的。如果没有重合，则软件可能会产生一个警告或者警报展示给用户，另外或或者，软件可能会例如在比对步骤之前停止程序。

根据一个实施方案，性别测试包括至少一个下列步骤：测试X染色体和/或Y染色体上的性别特异性基因中的至少一个位置；利用X染色体和/或Y染色体上的人类基因组样品的比对差异；细胞遗传学测试；FISH分析；CGH分析，或任何其它可以直接或间接确定人类核酸样品的性别的实验方法。

因此，性别测试也可以是人类核酸样品的FISH分析(荧光原位杂交分析)[GallJ.G.1969]的结果或副产物。因此，性别测试也可以是人类核酸样品的CGH分析(比较基因组杂交)[Kallioniemi A.et al.,1992]的结果或副产物。

性别测试可以有效且可靠地区分男性或女性人类核酸样品。

由于来自一个血统或种族的个体共享许多SNP，从而使他们与其他血统或种族区分开，因此可以通过检查一系列血统确定的SNP来确定最适合读长比对和变体识别的PHREG。因此，可以基于血统测试的结果从一组人类参照基因组中选择PHREG。

可以在上游基因组分析管线步骤中使用不同的实验设置，以便在进行比对之前确定个体的血统，从而确定最佳匹配的PHREG参照并避免错误。

1)血统测定可以基于用在人类核酸样品上的机器学习算法，或者利用血统特异性变体的另一个分类方案。特别地，血统测试方法可以基于机器学习，其利用了外显子位置的基因型，例如，超过100、500、1000、2000或优选地超过5000个外显子位置的基因型。

2)可以基于NGS数据或者另一个实验方法来确定相关的基因型，例如SNP阵列，如其在法医研究中所做的那样[Fondevila et al.,2013]。这里，使用非编码SNP可以帮助确定种族。

3)可以将与选项2)的法医SNP阵列中测试的相同的非编码SNP(加上侧翼区域)添加到现有的靶向NGS小段(panel)中，以确定相关的基因型。

特别地，血统测试可以包括使用至少一个基因组位置的基因型。

在一个特别的实施方案中，血统测试包括测试至少一个选自附录中的序列方案的基因。为产生准确的结果，展示了来自附录中的序列方案的249个基因。

另外或或者，血统测试可以包括测试SNP阵列和/或SNP芯片和/或测试来自Sanger测序或质谱的标记，或任何其它用于确定相关基因型的实验方法。

在一个特别的实施方案中，血统测试包括测试选自下组基因的至少一个基因：ABL2、ATP1A3、CIC、CYP2C8、CYP2C9、EPHA3、EPHA7、ERBB3、ERG、ETV1、F2、FAS、HFE、IL11RA、IL2RA、ITGB6、KIF11、KIT、KLK3、LRP6、MDM4、NAT2、NTRK2、PDGFB、PIK3R1、PLA2G3、PLAU、PRKCB、RICTOR、SLC7A11、STAT3、T、TSC1、VCAM1、VDR、VEGFB、ACVRL1、AXL、CA9、CALCR、CASP9、ENG、EPHB1、ERBB4、ESR1、FGFR2、HPSE、HSP90AA1、ITK、MRE11A、PLK1、PTPRC、SERPINE1、SMC4、TERT、TLR3、WISP3、WT1、XRCC1、ANGPT2、ARID2、BARD1、CBR3、CDH2、CYP1B1、DDR2、DNMT3A、EPCAM、ERCC2、FANCG、FANCL、GSTP1、IRS2、ITGB1、JAK3、LHCGR、MSH6、NCF2、RNF43、SLC5A5、TMPRSS2、TNFRSF8、AKT1、CD248、CD4、ESR2、EZH2、IGF1R、ITGAV、ITGB2、KLHL6、MAP3K1、MET、MLL、MTHFR、NFKB1、NUP93、PARP8、RB1、RPE65、TSHR、ABL1、BLM、CYP19A1、DPP4、EPHA6、ERBB2、EWSR1、FOXP4、ITGAM、KDM5A、LPA、LTK、MLH1、PBRM1、PHLPP2、SF3B1、TNFRSF10A、ABCG2、ACPP、ADAM15、DPYD、EPHA5、EPHB6、FOLH1、KDR、MSH3、MST1R、NTRK1、ROCK2、SLC6A2、TET2、TGM2、TH、ABCB1、CD22、CD40、CD44、CDH20、CYP11B2、ERCC5、GPR124、IL7R、ITGB3、ITGB5、NCL、NOD2、NR4A1、PGR、PLCG1、PPP2R1A、PRAME、PTCH2、RET、SETD2、XPC、ASXL1、EPHB4、PLA2G6、SYK、TET1、EP300、FLT1、ITGA1、LOXL2、PDGFRB、PIK3CD、SSTR5、TEC、APC、ATR、CLU、CREBBP、CYP2D6、EML4、MMP2、PARP2、PDGFRA、TRPM8、CSF1R、DOT1L、FGFR3、FGFR4、GLP2R、IKBKE、JAK1、NOTCH2、SPEN、SPG7、BRCA1、CYP11B1、GNAS、ITGA5、LTF、NRP2、PTK2B、TNKS、ABCC1、CEACAM5、CYP4B1、EGFR、FLT3、INSR、PTCH1、SMARCA4、ZNF217、BCR、EEF2、SELP、SLCO1B1、ABCC2、FLT4、MTR、IL4R、MTOR、RPTOR、TEK、ATM、CARD11、FANCD2、MEFV、NF1、TP73、BRCA2、CD109、PTPRD、ABCC6、IGF2R、P2RX7、ROS1、ACE、PARP1、PRKDC、CENPE、TSC2、ALK、NOTCH1、TNC、NOTCH3、POLE、MLL2、MYH11、POLD1、GRIN3B、F5、FANCA、LRP1B、LRP2、VWF。

在一个更特别的实施方案中，血统测试包括测试至少一个从附录1列出的基因组坐标中选择的基因组坐标。附录1描述了用于血统分类器(ancestry classifier)的特征的基于GRCh37的基因组坐标。根据BED文件标准(https://www.ensembl.org/info/website/upload/bed.html)设置前3列的格式，并(从左至右)对应于染色体、特征的从0开始的起点以及特征的从0开始的终点(即，在特征的结尾后的第一个位置)。第4列显示了与这个位置的分类器有关的碱基，第5列显示了对应的基因名称。

基因名称由HUGO基因命名委员会(Gene Nomenclature Committee，HGNC，https://www.genenames.org/)批准。HGNC负责批准人类基因位点的唯一符号和名称，包括蛋白质编码基因、ncRNA基因和假基因(pseudogenes)，以便进行明确的科学交流。本申请中使用的基因名称是于2013年8月检索得到。

在另一个特别的实施方案中，血统测试包括至少一个列在附录2中的SNP[Fondevila et al.,2013]。附录2表明了SNP所位于的染色体号码(左列)、确切的染色体位置(中间列)以及相应的rs号(右列)。rs号是NCBI(国家生物技术信息中心，NationalCenter for Biotechnology Information)在其SNP数据库(dbSNP，https://www.ncbi.nlm.nih.gov/projects/SNP/)中分配的登录号(accession number)，并且被广泛用于指整个基因组数据库中的特定SNP。当研究人员识别出SNP时，他们向dbSNP数据库发送一份报告(该报告包括紧邻SNP的序列)。如果发送了重叠的报告，它们将被合并到相同的非冗余参照SNP簇中，该簇被分配了唯一的rsid。有关更多信息，请参见以下URL http://www.ncbi.nlm.nih.gov/sites/books/NBK44406/。

这样的血统测试可以包括基因和/或基因组测试，来区分血统类别。根据1kG计划，此类血统类别可以被定义为AFR、AMR、EAS、EUR、SAS。然而，该方法不限于1kG计划数据，例如，如果出现具有更多个体/种族的更全面的数据集，则可以将这些数据集用于相同的目的。

根据一个实施方案，人类核酸样品包含一组来自下一代测序程序的读长，并且其中所述比对包括将所述读长映射至选择的PHREG的步骤。另外或或者，人类核酸样品包括一组来自目标测序程序，例如，来自小段测序(panel sequencing)的读长。

有利地，该方法可以无缝地集成到任何现有的基于对HRG进行读长映射的NGS分析工作流程中。

通过将读长映射至选择的PHREG来比对人类核酸样品和选择的PHREG，这样可以预先假定已准备好测序文库，其通过对DNA或cDNA样品的随机断裂(fragmentation)，然后进行5'和3'-衔接子连接(ligation)制备。在一些实施方案中，断裂和连接反应被组合为一个步骤，然后对衔接子连接的片段进行PCR扩增。

通过将读长映射至选择的PHREG来比对人类核酸样品和选择的PHREG，这样可以预先假定这组DNA片段已经过测序，其结果导致得到长度为大约在28个碱基对(bp)和1000个碱基对(bp)之间的读长[Goodwin S.et al.,2016]。该组包括了足够多的读长以达到预定的目标区域的覆盖度(通常在几倍到几千倍之间)，从而适于回答所问的实验问题。

在一个实施方案中，下一代测序程序涉及全外显子组测序。在另一个实施方案，下一代测序程序涉及全基因组测序。术语“全外显子组测序”通常是指对一个基因组中所有编码蛋白质的基因(称为外显子)进行测序的技术。它包括首先仅选择编码蛋白质的DNA子集(称为外显子)，然后使用任何高通量DNA测序技术对该DNA进行测序。人类具有约18万个外显子，约占人类基因组的1.5％，或约3000万个碱基对。特别地，可以通过下一代测序来进行外显子测序。“全基因组测序(whole genome sequencing)”(也称为WGS、全基因组测序(full genome sequencing)、完整基因组测序(complete genome sequencing)或整条基因组测序(entire genome sequencing))是一种实验室过程，可一次确定生物体基因组的完整DNA序列。这需要对生物体的所有染色体DNA以及线粒体中包含的DNA进行测序。

根据本发明的另一个方面，一种对人类核酸样品进行基因分析的计算机系统包括：

a)第一模块，所述第一模块包含用于提供一组人类参照基因组的计算机指令；

b)第二模块，所述第二模块用于测试人类核酸样品的性别和/或血统；

c)第三模块，所述第三模块包含用于基于所述性别和/或血统测试的结果，从所述一组人类参照基因组中选择一个或多个种群特异性的人类参照基因组(PHREG)的计算机指令；以及

d)第四模块，所述第四模块包含用于将所述人类核酸样品与所述确定的PHREG进行比对的计算机指令。

特别地，计算机系统可以适合于或可以被配置用于执行以上公开的任何一种方法。因此，应当理解，已在方法的上下文中描述的特征也适用于计算机系统，反之亦然，将在计算机系统的上下文中描述的特征也适用于以上方法。

这些模块可以是存储在诸如永久或可重写存储装置之类的机器可读存储介质上，或分配给计算机装置的存储介质(诸如移动存储介质(例如CD-ROM、DVD、蓝光光盘、存储棒或存储卡))上的软件模块、软件例程或软件子例程。另外或或者，可以在诸如服务器或云服务器之类的计算机装置上提供模块，以用于下载，例如通过诸如因特网之类的数据网络或通过诸如电话线或无线线路之类的通信线路来下载。

本文公开的任何模块可以是不必在物理上彼此分离的功能单元。例如，如果在软件包中实现多个功能，则模块的几个单元可以以单个物理单元的形式实现。

本文公开的计算机模块不一定是集成系统的一部分，而是可以分布在通过通信网络彼此交互的几个单独的系统上。

根据一个实施方案，测试人类核酸样品的性别和/或血统的第二模块是包括计算机指令的计算机模块。另外或或者，第二模块可以包括湿实验室实验，例如进行FISH测试的实验。可以通过电子或目测方式分析FISH测试的结果，以确定样品的性别。

根据本发明的另外一个方面，计算机程序包含下述指令：当计算机执行所述程序时，所述程序使得所述计算机执行根据由上所述的任何一种方法的步骤a)、b)、c)和d)。

根据本发明的另外一个方面，计算机可读存储介质包含下述指令：当计算机执行所述指令时，所述指令使得所述计算机执行根据由上所述的任何一种方法的步骤a)、b)、c)和d)。

如上面已经讨论的，本发明的方法特别适合于鉴定患者基因组中的变化，该变化对某一给定疾病具有指示性或对患者对给定疗法的敏感性具体指示性。

在本文中，术语“疾病”包括以一种或多种基因组变化为特征的任何疾病。这包括癌症、自身免疫性疾病、心血管疾病和任何遗传性疾病。该患者可以是任何物种，但优选是哺乳动物，更优选是人类。

取决于个体疾病和疗法，本领域技术人员将能够选择对患者有益的个体治疗模式。

所以，在本发明的进一个方面，本发明涉及一种诊断患者的疾病的方法，包括：

获取所述患者的疾病指征的识别信息；

获得所述患者的核酸样品；并且

基于对本文所述的人类核酸样品进行基因组和/或基因分析的方法对该核酸样品进行基因组和/或基因分析，从而确定患者的疾病状态。

可以通过本领域中已知的任何方法来检索得到疾病指征的识别信息，例如，用户输入、电子健康记录或电子病历或包含病历的患者数据库。

在本发明的这一方面的上下文中，术语“疾病状态”在一个实施方案中是指患者的疾病得到确认。在另一个实施方案中，该术语意味着疾病得到更准确地诊断，即，疾病的个体亚型得到确认。

本发明还涉及一种治疗患者疾病的方法，包括：

获取所述患者的疾病指征的识别信息；

获得所述患者的核酸样品；并且

基于本发明所述的方法对所述核酸样品进行基因组和/或基因分析，从而确定患者的疾病状态，并治疗患者。

在本发明的另一个方面，本发明涉及一种确定患者是否对某种药物的治疗具有易感性的方法，包括：

获取所述患者的疾病指征的识别信息；

获得所述患者的核酸样品；

基于本文所述的对人类核酸样品进行基因组和/或基因分析的方法对所述核酸样品进行基因组和/或基因分析；

获取所述患者的所述疾病指征的可能疗法；

进行变体识别和解释；以及

基于变体解释对检索到的可能疗法进行分类，其中疗法被分类为对患者适用或对患者禁用。

用这种方法，有可能为患者决定有哪种疗法或者哪种疗法有利。例如，有可能决定患者是否对一种给定疗法易感或者是否可以预期一种给定疗法只有可接受的副作用。

可以再次通过本领域已知的任何方法来检索疾病指征的识别信息，例如，用户输入、电子健康记录或电子病历或包含病历的患者数据库。

可以通过本领域中已知的任何方法，例如从数据库，来检索用于患者疾病指征的可能疗法。

本发明也涉及一种治疗患者的方法，包括：

获取所述患者的疾病指征的识别信息，

获得所述患者的核酸样品，

基于本文所述的对人类核酸样品进行基因组和/或基因分析的方法对所述核酸样品进行基因组和/或基因分析，

获取所述患者的所述疾病指征的可能疗法，

进行变体识别和解释，

基于变体解释对检索到的可能治疗方法进行分类，其中疗法被分类为对患者适用或对患者禁用，

选择一种指示的疗法，以及

根据所选的疗法治疗患者。

可以再次通过本领域已知的任何方法来检索疾病指征的识别信息，例如，用户输入、电子健康记录或电子病历或包含病历的患者数据库。

可以再次通过本领域已知的任何方法，例如从数据库中，检索用于患者的疾病指征的可能疗法。

附图说明

通过参考以下结合附图的描述，本公开的前述和其他目的、方面、特征和优点将变得更加明显和更好地理解，其中：

图1是描述根据本发明的用于人类核酸样品的基因组和/或基因分析的方法的流程图；

图2是描述根据本发明的数据分析方法的流程图；

图3表示读长映射步骤；

图4是描述根据本发明的用于人类核酸样品的基因组和/或基因分析的方法的流程图；

图5是表示在MH面板数据(MH Panel data)上计算出的为性别分类器(sexclassifier)选择的特征的按类分布图；和

图6是两个Ansextry分类器(性别分类器和血统分类器)和EthSEQ的内存使用情况和运行时间的箱线图。

附图详述

图1示出了用于人类核酸样品的基因组和/或基因分析的一般工作流程，其包括提取人类核酸样品、制备测序文库、测序和随后的数据分析的过程。在本发明的上下文中，提取人类核酸样品、制备测序文库和测序的过程可能涉及众所周知的标准过程，因此将不作更详细的解释。图2更详细地示出了本发明的数据分析部分。

图2示出了图1中的数据分析步骤，包括性别和血统测试步骤(第一步)、比对(或读长映射)步骤、变体识别步骤和注释步骤。读长映射计算模块的输入文件是原始序列数据(例如，以FASTQ文件的形式)。读长映射计算模块的输出文件是，例如，BAM文件，其作为变体识别计算模块的输入文件。变体识别计算模块的输出文件是，例如，VCF文件。随后的注释计算模块可能注释来自VCF文件的数据，并将其以所需的格式(例如，PDF、HTML或其它)进行输出。文件格式只是示例性的，也可以是不同的格式，例如，不是BAM，而可能是SAM或CRAM文件或其它。图2中的数据分析管线也可能包括一个将输入或输出文件从一种文件格式转换成另一种文件格式的计算机模块。

图2也比较了现有技术的情况和本发明的情况。现有技术的方法(由图2中的“A”表示)没有提供性别和血统测试，所以其比对和变体识别是相对于标准HRG进行的。根据本发明的方法(由图2中的“B”表示)提供了性别和血统测试，允许选择一种或多种确定的PHREG。随后的比对和变体识别是相对于确定的PHREG进行的。

图3是示例性读长映射步骤的示意图。在该示例中，NGS读长带有血统特异性的SNP“A”。血统特异性的SNP“A”位于先前未发现的生物标记物变体“G”的附近。所述附近可以远至读长的长度。

在比对过程中，当NGS读长相对于标准HRG进行比对时会产生2个不匹配，即，血统特异性SNP和生物标记物变体。但是，在比对过程中，当NGS读长相对于相应的PHREG进行比对时仅会产生1个不匹配，即，生物标记物变体，因为该PHREG已经在血统特异性的位置进行了修改以使其与该血统特异性的SNP一致。

比对算法使用了评分系统，该评分系统涉及对每个不匹配的惩罚(penalties)和/或测序的读长与选择的参照基因组之间的间隙。读长随后被比对至得分最佳的位置，或由于总分很低或太多基因组位置有相同的比对得分而根本不被比对。由于比对算法中的不匹配惩罚，与PHREG相比，读长相对于HRG比对的可能性较小，尤其是如果之后发现的变体位于读长的长度内。因此，该读长会被丢弃，或更糟糕的情况，甚至在HRG的错误位置进行比对。

因此，PHREG具有这样的效果：其可以挽救源自血统特异性变体区域的读长，尤其是如果它们带有除了血统特异性变体的其他变体(例如，致病变体)时。这使得能够检测到先前未被发现的生物标记物。

图4示出了流程图，其描绘了根据本发明的用于人类核酸样品的基因组和/或基因分析的方法。

第一步中，向包含处理单元的系统提供一组人类参照基因组。为此，该系统的第一计算机模块可以从一个远程设施(例如，互联网数据库)下载参照基因组。处理单元可以是任何可编程计算机装置，其基本上包括至少一个处理器，该处理器具有内部存储器，例如RAM(随机存取存储器)，其允许存储和执行指令。处理单元可以访问非易失性存储装置(non-volatile storage means)，其可以存储数据集和计算机文件(例如，人类参照基因组以及患者的临床数据和遗传概况)。该系统可以访问通信网络，例如LAN或互联网。

第二步中，由系统的一个计算机模块将人类参照基因组调整至编码水平，优选地由系统用户来操作。该编码水平可以包括独特的核苷酸代码或模糊的核苷酸代码。在一些实施方案中，提出了四种不同的水平将人类参照基因组调整至一个种群，其中两个水平限制在独特核苷酸代码(A、C、G、T)，另两个水平利用根据IUPAC命名的模糊核苷酸编码，特别地，最大程度的保守错误校正(Maximally conservative error correction)、多数等位基因水平、非罕见等位基因水平和对所有观察到的等位基因进行完整建模。

第三步中，提供患者的人类核酸样品。为此，系统的另一个计算机模块可以从测序实验室下载原始序列数据(例如，以FASTQ文件的形式)，该测序实验室对目标样品在远程平台上进行测序。在一个替代的实施方案中，在本地进行测序且内部进行结果转移。在第三步的上下文中，系统还可以从其他输入源(例如，有关患者所患疾病的信息、有关当前治疗的信息等)接收患者的其他临床数据。临床患者数据可以例如是直接从患者那里收到、例如可以在键盘上键入或可以从在键盘上键入的自由文本推导出来，也可以从GUI中的多项选择元素接收。也可以从电子健康记录(electronic health record,EHR)或电子病历(electronic medical record,EMR)中检索得到临床患者数据，其可能在芯片卡上或在可通过通信网络检索的数据库中。

第四步中，对人类核酸样品进行性别和/血统测试。同样，可以在本地进行测试，或者系统的另一个计算机模块可以通过通信网络从外部服务提供者得到测试结果。可以由第二计算模块或者另一个湿实验室试验(wet lab experiment)执行性别和/或血统测试。

第五步中，基于性别和/或血统测试的结果，从一组人类参照基因组中选择一个或多个PHREG。可以由第三计算模块执行该选择。

第六步中，将人类核酸样品与选择的PHREG进行比对。所述比对包括将来自NGS程序的读长组映射至选择的PHREG。可以由第四计算模块执行比对，且输出文件可以是BAM文件。

第七步中，相对于选择的PHREG对比对的人类核酸样品进行变体识别。在变体识别前，系统的一个计算模块可以将人类参照基因组再次调整至编码水平，优选地由系统用户来操作。编码水平可以包括独特的核苷酸代码或模糊的核苷酸代码，且与比对步骤中使用的编码水平不同。可以使用最合适的最新算法来识别变体。可以由第五计算模块执行变体识别，并且输出的内容可以包括变体识别格式(variant call format,VCF文件)的序列数据，该序列数据以相对于PHREG的变体形式存在。

第八步中，执行变体解释。因此，系统可以包括另一个适于执行分析所识别的变体的后处理计算模块。在一个实施方案中，后处理模块可以分析能暗示患者中是否存在疾病的一组基因和/或变体。另外或或者，后处理模块可以为患者的疾病确定一组疗法，考虑到其他临床患者数据，可以基于患者的遗传数据，特别是基于识别的遗传变体为患者确定最适合的个体化治疗方案。在另一个实施方案中，后处理模块可以执行统计学分析并确定来自识别的变体的突变记载、核苷酸替换率和热点突变。

发现的变体还可以作为分类器的输入，所述分类器可以预测治疗效果或治疗安全性或用于诊断和/或治疗目的。

第九步中，可能会产生并提供诊断和/或治疗意义。为此，该系统可以包括与第三、第四、第五和后处理模块中的任一个功能连接的输出接口，使得可以输出它们的结果。输出接口可以耦合到任何显示装置或打印机，使得可以呈现由处理单元计算的信息。此外，可以有用于内联网和/或互联网的通信系统的链接，例如用于经由输出接口实现发送和接收电子邮件的程序。

图5示出了根据MH面板数据(MH Panel data)计算出的按性别分类的所选特征的分布图(F：女性；M：男性。彩色竖线表示类别中位数。顶部：chrX/chrY对齐的读长比率；中：对于chrX上的500个常见SNP位置，bin 0.8-1.0中多数等位基因频率的分数；底部：chrY上正确配对的读长的百分比。应该在以下描述的实施例的上下文中查看图5。

图6显示了300个TCGA全外显子样品上两个Ansextry分类器和EthSEQ的内存使用情况(顶部，以GB为单位)和运行时间(底部，以分钟为单位)的箱线图。应该在以下描述的实施例的上下文中查看图6。

实施例1

AnSextry是基于机器学习的工具，其使用来自全外显子组测序数据的读长比对从样品中得到性别和血统。已经知道对这两种特征的自我声明是不可靠的，AnSextry的预测在样品交换检测和无偏见的基因组变体解释中都是有用的，尤其是在大型队列研究(largecohort study)中。AnSextry对1,300多个样品的基准测试表现出高精度、用时少和较低的内存要求。

1.引言

在过去十年间，随着成本的急剧下降，大型队列的下一代测序已变得司空见惯[Cancer Genome Atlas Research Network et al.,2013；Rand et al.,2016]，全外显子组方法在大型研究中起着重要作用，尤其是在精准医学或疾病的全面表征领域。在这种情况下，对样品的血统和性别的可靠了解会带来多种好处。首先，它可以用作简单的质量控制，以帮助识别由于样品处理中涉及的复杂方案和手动操作而引起的样品交换。其次，血统对于变体影响的解释、规避大多数基因组研究和人类参照基因组中存在的强烈欧洲偏倚、以及改善对具有多种血统的人的临床护理至关重要[Petrovski et al.,2016；Mersha etal.,2015；Fakhro et al.,2016]。最后，血统被广泛用于遗传关联研究中，以避免由于人口分层而与疾病产生虚假关联[Wu et al.,2011]。自我声明性别和血统常常是不可靠的[Mersha et al.,2015；Ainsworth，2015]，这也呼吁使用基因组信息对性别和血统进行鉴定。

开发出了一种基于逻辑回归的机器学习方法，AnSextry，以快速可靠地从全外显子组测序配对末端读长比对(whole-exome sequencing paired-end read alignment)中表征性别和血统。该算法依靠标准文件格式，且容易被整合至现有的下一代测序分析工作流程。它提供了一个现成的模型，只需要输入一个简单的BAM文件。此外，AnSextry的低内存要求使其可以在台式计算机上运行。基准测试程序表明，AnSextry在精确性、运行时间和内存使用方面与EthSEQ[Romanel et al.,2017](这是唯一已知的另一种基于BAM文件的全外显子组血统推断工具)相比均具有优势。迄今为止，尚无其他公开的性别预测方法。

2.方法

2.1算法

基于全外显子组测序的配对末端读长比对，准备了两个分类器，它们可以推断一个个体最可能的性别和血统。该工具利用了读长映射和各个基因型之间的差异来进行预测。

性别和血统分类器基于使用Python和scilit-learn机器学习库的逻辑回归。两者的特征(features)均来自一个输入的BAM文件。配对末端读长的比对用了BWA0.7.15默认的比对设置，而没有类似局部重新比对或去除重复的后处理步骤。用GRCh37作为参照基因组，没有非染色体超重叠群(supercontig)，但有屏蔽的伪常染色体区域PAR1和PAR2，以防止在X和Y染色体上发生比对变形。在本发明的内容中，术语“超重叠群”通常被理解为一组有序的重叠群，即碱基顺序为高置信度已知的连续长度的基因组序列。

性别分类器使用L1正则化(L1-regularization)与两类(two-class)逻辑回归一起运行，并返回得到每个类(class)的概率。5倍交叉验证(5-fold cross-validation)被用来确定合适的正则化长度。选择在训练数据的精确召回曲线下面积最大的模型作为最佳模型，并在测试数据集上进行评估。

血统分类器基于使用L2正则化和主要成分分析(Principle ComponentAnalysis,PCA)的多项逻辑回归，并返回得到在1000基因组计划中定义的每一个五大洲血统的概率：非洲人(AFR)、Ad混血美洲人(AMR)、中东人(EAS)、欧洲人(EUR)、南亚人(SAS)(The 1000 Genomes Project Consortium et al.,2015)。5倍交叉验证(5-fold cross-validation)被用来确定合适的参数。选择训练数据的F1得分最高的模型作为最佳模型，并在测试数据集上进行评估。

2.2特征

性别分类器的特征是基于X和Y染色体之间的比对差异(图5)。使用了chrX和chrY的读长率以及chrY上正确配对的读长的百分比。而且，复合了chrX上的500个常见外显子SNP位置的多数等位基因频率。为了避免种群偏倚，选择了不同的主要血统间经常出现的SNP。

对于血统分类器，从第2.3节描述的1000基因组数据确定了所有常染色体SNP的基因型，这些常染色体SNP处于常用安捷伦全外显子试剂盒(Agilent All Exon kits)(V5、V6、V6+COSMIC)的靶标区域与Molecular Health Pan-Cancer基因小段(panel)(靶标大小为2.9Mbp)的交集。使用特征选择来保留有意义的SNP，这些SNP显示了不同血统之间的变体，其产生对应于5,040个基因组位置的10,000个基因型，可被用作该分类器的特征。可在附录1中找到对应的BED文件，其可被用来确定与任何靶标测序试剂盒的重叠。

2.3数据

为了从不同的血统获取数据，使用了来自1000基因组计划第3阶段的1735名个体的基因型数据，以训练和测试血统分类器。用于分类的各个洲的血统(AFR、AMR、EAS、EUR、SAS)和个体被随机选取以获得平衡的类(class)。694名个体是测试组的一部分。

从TCGA(cancergenome.nih.gov)下载了300名自我报告种族和性别的个体的主要全外显子对照数据(primary whole-exome control data)作为测试组，对应于三种癌症(尿路上皮膀胱癌、肺腺癌/鳞状上皮细胞肺癌、胃腺癌)。所有样品均由安捷伦SureSelect人类全外显子50Mb试剂盒(Agilent SureSelect Human All Exon 50Mb kit)进行测序。随机选择记录以达到与TCGA类别相对应的平衡的类的规模：150名男性和150名女性个体，以及100名白人、100名亚洲人和100名黑人或非裔美国人。

使用了来自988名癌症患者自我报告性别的靶标序列数据来训练和测试性别分类器，该序列数据用Molecular Health Pan-Cancer基因组合(gene panel)来测序得到。随机选择个体以达到男女类(class)平衡。随机选择396例作为性别分类器的测试数据。上述300个TCGA案例用作附加测试集。

3.方法

3.1性别分类器

使用Molecular Health Pan-Cancer基因组合(gene panel)测得的592组数据集来训练性别分类器。比对了配对末端读长并根据方法小节中描述的内容计算了特征。在用交叉验证(cross-validation)调整方法后，对两组测试数据集(用上述基因小段(panel)测序的396名个体，以及具有全外显子数据的300名TCGA个体)的表现进行了评估。

在panel的测试数据上，性别分类器达到了97.5％的平均准确度，其中10名个体(5名男性和5名女性)分类错误(参见表1)。错误分类与较低的覆盖度没有关联。

性别	预测的性别	真实性别的预测概率[％]	平均覆盖度(两次重复)
				F	M	39.0	2579
F	M	15.3	2099
				F	M	33.8	1656
F	M	17.1	1787
				F	M	0.8	1797
M	F	0.0	6016
				M	F	28.0	2401
M	F	0.3	3603
				M	F	0.0	1606
M	F	0.0	1705

表1：MH面板(panel)测序的个体的详细信息，其中预测的性别与自我声明的性别不匹配。所有使用的样品的中位覆盖度为2116x。所有分类错误的样品的平均覆盖度接近或高于该中位数，表明错误的预测似乎与低于中位数的覆盖度无关。

由于一般人口中性别发展障碍的发生频率为1％[Ainsworth，2015]，因此，某些错误分类的例子实际上可能是正确分类的，但却具有不正确的自我声明的性别。

在TCGA测试数据上，性别分类器的精度达到100％。所有300名个体均被正确分类。在运行时间和内存使用方面，性别预测在所有情况下用时不到一分钟，平均内存使用为526MB(图6)。

3.2血统分类器

在来自1000基因组计划的1041个数据集上训练血统分类器。如第2.2节所描述，个体基因型被用作特征。在两组测试数据集上预测了最佳表现的模型：来自1000基因组计划的剩下的694名个体，以及具有被测序的全外显子的300名TCGA个体。

在1000个基因组测试数据上，血统分类器的平均精度达到了99％，在亚洲血统中表现最佳(南亚和东亚均为100％精度)，其次是非洲和南美血统(99％精度)和欧洲血统(98％)。总共694人中只有5个人被错误分类。

在300个TCGA外显子组测试数据集上，血统分类器的准确度略低，为96.33％，共有11个人被错误分类。将这些结果与EthSEQ进行比较[Romanel et al.,2017]，EthSEQ是目前已知的唯一另一种血统预测方法，可以提供合适的预先计算模型，并且可以在单个全外显子BAM文件中即时使用。但是，两种方法的结果高度一致，EthSEQ达到的精度略低(94％)，共有18个人被错误分类。此外，EthSEQ所需的运行时间和内存要高得多：血统分类器的平均运行时间为28秒，平均内存使用量为540MB，而EthSEQ即使使用多线程(4核)，平均仍需要4.8分钟和14.7GB(图6)。

重要的观察结果是，对于错误分类的数据集，两种算法之间的一致性也很高：11名个体中有10名个体的血统预测与TCGA提供的血统不匹配，这也被EthSEQ进行了不同的分类，且10名中有8名被两种方法均预测了同样的血统。这表明，这些人中至少有一部分可能被TCGA错误分类，TCGA的种族信息基于自我声明。在10个一致的案例中有6个被预测为AFR或AMR，这与Mersha等人的观点一致，他们声称，自我声明错误在非裔美国人和拉美裔人口中最普遍。表2显示了错误分类的个体。

表2：有关TCGA个体的详细信息，其中(通过Ansextry、EthSEQ或两者)预测的血统与TCGA自我声明的种族不匹配。TCGA种族类别包括“黑人或非裔美国人”、“白人”和“亚裔人”。白色行对应于Ansextry或EthSEQ与TCGA都不匹配的种族的样品。浅灰色行是仅EthSEQ的预测与TCGA不匹配的样品；并且深灰色行是只有Ansextry的预测与TCGA不匹配的样品。当相应的基因位点没有足够的覆盖度时，从参照中推断出基因型以用于Ansextry预测。所有样品的中位数覆盖度为91x，这表明大多数错误预测的样品具有中值或中位数以上的覆盖度，因此，错误预测似乎与中位数以下的覆盖度无关。同样，在所有样品中用于Ansextry分类的推测的基因型的中位数为390，接近错误预测的Ansextry样品的中位数(393)。在所有300个TCGA样品中，推断的基因型的数量在227(最小)和690(最大值)之间变化，表明10-15％的推断的基因型似乎对Ansextry预测没有负面影响。

有趣的是，唯一一个被AnSextry错误分类而未被EthSEQ错误分类、且被TCGA分类为白人但被血统分类器预测为AMR的个体实际上被预测为概率为54.7％的AMR和45.1％的EUR的混合血统。

4.结论

AnSextry是一种新颖的方法，可根据来自整个外显子组的比对的配对末端读长可靠轻松地确定一个个体的性别和血统，或者在目标大小允许的情况下，展示靶向测序实验。该工具提供了两个基于逻辑回归的基于Python的分类器，血统预测代表了人口遗传学领域中主要基于PCA的方法的替代方法。AnSextry提供了一个现成的参考模型，并且需要最少的用户输入。它快速、准确且易于使用。

免责声明

在本文档中，因为不同的作者出于同一目的使用不同的术语，术语“血统特定”/“种族特定”/“种群特定”可互换使用。

参考文献

1.Lander,E.S.et al.Initial sequencing and analysis of the humangenome.Nature 409:860–921(2001).[PMID:11237011]

2.Church,D.M.et al.Modernizing reference genome assemblies.PLoSBiol.9:e1001091(2011).[PMID:21750661]

3.Harrow,J.et al.GENCODE:the reference human genome annotation forThe ENCODE Project.Genome Res.22:1760-1774(2012).[PMID 22955987]

4.ENCODE Project Consortium.An integrated encyclopedia of DNAelements in the human genome.Nature 489:57-74(2012).[PMID:22955616]

5.1000 Genomes Project Consortium et al.A global reference for humangenetic variation.Nature 526:68-74(2015).[PMID:26432245]

6.Li H&Durbin R.Fast and accurate short read alignment with Burrows-Wheeler transform.Bioinformatics 25:1754-1760(2009).[PMID:19451168]

7.DePristo,M.A.et al.A framework for variation discovery andgenotyping using next-generation DNA sequencing data.Nat.Genet.43:491-498(2011).[PMID:21478889]

8.Horton,R.et al.Variation analysis and gene annotation of eight MHChaplotypes:the MHC Haplotype Project.Immunogenetics 60:1-18(2008).[PMID:18193213]

9.Pei,B.et al.The GENCODE pseudogene resource.Genome Biol.13:R51(2012).[PMID:22951037]

10.Degner,J.F.et al.Effect of read-mapping biases on detectingallele-specific expression from RNA-sequencing data.Bioinformatics 25:3207-3212(2009).[PMID:19808877]

11.Brandt,D.Y.C.et al.Mapping Bias Overestimates Reference AlleleFrequencies at the HLA Genes in the 1000 Genomes Project Phase I Data.G3 5:931-941(2015).[PMID:25787242]

12.Novak A.；Hickey G.；Garrison E.；Blum S.；Connelly A.；Dilthey A；Eizenga J.；Elmohamed M.；Guthrie S.；Kahles A.；Keenan S.；e Kelleher J.；KuralD.；Li H.；Lin M.；Miga K.；Ouyang N.；Rakocevic G.；Smuga-Otto M.；Zaranek A.；Durbin R.；McVean G.；Haussler D.；(https://www.biorxiv.org/content/biorxiv/ early/2017/01/18/101378.full.pdf)

13.Paten B,Novak AM,Eizenga JM,Garrison E.Genome graphs and theevolution of genome inference.Genome Res.5:665-676(2017)[PMID:28360232]

14.Snyder M.,et al.Personal genome sequencing:current approaches andchallenges.Genes Dev.5,423-431(2010)[PMID:20194435]

15.Young,A.L.et al.A new strategy for genome assembly using shortsequence reads and reduced representation libraries.Genome Res 2:249-256(2010)[PMID:20123915]

16.Flicek,P&Birney,E.Sense from sequence reads:methods for alignmentand assembly.Nat Methods.6:S6-S12(2009)[PMID 19844229]

17.Chen R.&Butte A.J.The reference human genome demonstrates highrisk of type 1 diabetes and other disorders.Pac Symp Biocomput.2011:231-242(2011)[PMID:21121051]

18.International Human Genome Sequencing Consortium.2001.Initialsequencing and analysis of the human genome.Nature 409:860-921(2001)[PMID:11237011]

19.International Human Genome Sequencing Consortium.2004.Finishingthe euchromatic sequence of the human genome.Nature 431:931-945(2004)[PMID:15496913]

20.Schneider V.A.et al.Evaluation of GRCh38 and de novo haploidgenome assemblies demonstrates the enduring quality of the referenceassembly.Genome Res.5:849-864.(2017)[PMID:28396521]

21.[Editorial(October 2010)."E pluribus unum".Nature Methods.5:331.doi:10.1038/nmeth0510-331.(2010)[PMID:20440876]

22.Nielsen R.,Paul J.S.,Albrechtsen A.,Song Y.S.Genotype and SNPcalling from next-generation sequencing data.Nat.Rev.Genet.12:443-45.(2011)[PMID:21587300]

23.Fakhro,K.A.,Staudt M.R.,Ramstetter M.D.,Robay A.,Malek J.A.,BadiiR.,et al.The Qatar genome:a population-specific tool for precision medicinein the Middle East.Hum.Genome Var.3:16016 Human Genome Variation(2016)3,16016doi:10.1038/hgv.2016.16；published online 30 June 2016(2016)[PMID:27408750]

24.Zayed H.The Qatar genome project:translation of whole-genomesequencing into clinical practice.Int J Clin Pract.10:832-834 doi:10.1111/ijcp.12871.Epub 2016 Sep (2016)[PMID:27586018]

25.Sanger F.,et al.DNA sequencing with chain-terminatinginhibitors.Proc Natl Acad Sci USA.74:5463-5467.(1977)[PMID:271968]

26.Venter,J.C.et al.The Sequence of the Human Genome.Science 291:1304-1351.(2001)[PMID:11181995]

27.Petrovski S&Goldstein D.B.Unequal representation of geneticvariation across ancestry groups creates healthcare inequality in theapplication of precision medicine.Genome Biol 2016；17:157.doi:10.1186/s13059-016-1016-y.(2016)[PMID:27418169]

28.Koboldt DC,Ding L,Mardis ER,Wilson RK.Challenges of sequencinghuman genomes.Brief Bioinform.11:484-498.(2010)[PMID:20519329]

29.Dewey F.E.,Chen R.,Cordero S.P.,Ormond K.E.,Caleshu C.,KarczewskiK.J.et al.Phased whole-genome genetic risk in a family quartet using a majorallele reference sequence.PLoS Genet.2011 Sep；7(9):e1002280.doi:10.1371/journal.pgen.1002280.Epub 2011 Sep 15.(2011)[PMID:21935354]

30.Cao H,Wu H,Luo R,Huang S,Sun Y,Tong X et al.De novo assembly of ahaplotype-resolved human genome.Nat Biotechnol 33:617-622.(2015)[PMID:26006006]

31.Wu L.,Yavas G.,Hong H.,et al.Direct comparison of performance ofsingle nucleotide variant calling in human genome with alignment-based andassembly-based approaches.Sci Rep.2017 Sep 8；7(1):10963.doi:10.1038/s41598-017-10826-9.(2017)[PMID:28887485]

32.Meyer,L.R.et al.The UCSC Genome Browser database:extensions andupdates 2013.Nucleic acids research 41:D64-D69(2013).[PMID:23155063]

33.Sudmant,P.H.et al.An integrated map of structural variation in 2,504 human genomes.Nature 526:75-81(2015).[PMID:26432246]

34.Iqbal,Z.,Caccamo,M.,Turner,I.,Flicek,P.&McVean,G.De novo assemblyand genotyping of variants using colored de Bruijn graphs.Nature genetics 44:226-232(2012).[PMID:22231483]

35.Cornish-Bowden A.(1985).Nomenclature for incompletely specifiedbases in nucleic acid sequences:recommendations 1984.Nucleic Acids Res.13:3021-3030.(1985)[PMID:2582368]

36.Mersha T.B.,&Abebe T.Self-reported race/ethnicity in the age ofgenomic research:its potential impact on understanding healthdisparities.Hum.Genomics 9:1.(2015)[PMID:25563503]

37.Baye T.M.Inter-chromosomal variation in the pattern of humanpopulation genetic structure.Hum Genomics 5:220-240.(2011)[PMID:21712187]

38.Fondevila M.et al.Revision of the SNPforID 34-plex forensicancestry test:Assay enhancements,standard reference sample genotypes andextended population studies.Forensic Sci Int Genet 7:63-74.(2013)[PMID:22749789]

39.Ainsworth C.Sex redefined.Nature 518:288-291.doi:10.1038/518288a.(2015)[PMID:25693544]

40.Gall J.G.,Pardue M.L.Formation and detection of RNA-DNA hybridmolecules in cytological preparations.Proc.Natl.Acad.Sci.USA 63,Nr.2,1969,S.378–383,[PMID 4895535].

41.Kallioniemi A.et al.Comparative genomic hybridization for molecular cytogeneticanalysis of solid tumors.Science Band 258,Nr.5083,1992,S.818–821.

42.Goodwin S.,McPherson JD,McCombie WR.Coming of age:ten years ofnext-generation sequencing technologies.Nat.Rev.Genet.2016 May 17；17(6):333351

43.Al-Ali M,Osman W.,Tay G.K.,AlSafar H.S.A 1000 Arab genome projectto study the Emirati population.J.Hum.Genet.63(4):533-536(2018).[PMID:29410509]

44.Cancer Genome Atlas Research Network et al.The Cancer Genome AtlasPan-Cancer analysis project.Nat.Genet.,45(10),1113-1120(2013).

45.Rand,K.A.et al.Whole-exome sequencing of over 4100 men of Africanancestry and prostate cancer risk.Hum.Mol.Genet.,25(2),371-381(2016).

46.Wu,C.et al.A Comparison of Association Methods Correcting forPopula-tion Stratification in Case-Control Studies.Ann.Hum.Genet.,75(3),418-427(2011).

47.Romanel,A.et al.EthSEQ:ethnicity annotation from whole exomesequencing data.Bioinformatics,33(15),2402-2404(2017).

实施例2

用PHREG作为NGS读长映射的参照增加了临床相关的生物标记物的覆盖度

我们使用了来自GDC/TCGA[1]的741个种系样品，这些样品已通过全外显子捕获Illumina测序进行了测序。这组数据有155个非洲人(AFR)血统样品、33个拉丁/混血美洲人(AMR)血统样品、179个东亚人(EAS)血统样品、354个欧洲人(EUR)血统样品和20个南亚人(SAS)血统样品。用Novoalign 4.00.01将每个样品与标准人类参照基因组(HRG)GRCh37[3]、用我们的血统分类器分配的PHREG以及HSA PHREG进行比对。通过汇总所有GnomADv2.1血统[4](包括AFR、AMR、EAS、EUR和SAS)的变体数据产生HSA PHREG。

关于这些读长映射策略，我们比较了15,483个致病性ClinVar生物标记物版本2019-12[5]在涵盖1,288个基因的Gencode v31 CDS外显子[6]中的覆盖度。我们发现，当将ClinVar生物标记物与PHREG而不是HRG进行比对时，对以下的覆盖度有所增加：211(AFR)、147(AMR)、121(EAS)、173(EUR)、105(SAS)和162(HSA)(参见表3)。具有增加的覆盖度的大多数变体位于PHREG植入种群特异性核苷酸的位点附近。当将样品的读长映射到其最接近的PHREG时，比对期间的不匹配数会减少，因此覆盖度会增加，从而消除了与HRG比对时覆盖度下降的情况。

总之，我们的分析表明，使用正确的PHREG可以增加覆盖度，从而改善对临床相关生物标记物的检测。

表3(ClinVar_PHREG_coverage_diff_relative.xlsx)说明：

Gencode CDS外显子中的ClinVar生物标记物列表(基因名称|contig|开始|结束)，显示出与PHREG进行比对相比于与PHG进行比对的覆盖度差异。以每个PHREG(AFR、AMR、EAS、EUR、SAS、HSA)针对HRG的比对计算出的覆盖度为中位数，对血统的案例以及所有741个案例(HSA)得出差异。正数表示覆盖度增加，负数表示覆盖度减小。

表3

实施例2的参考文献

[1]https://portal.gdc.cancer.gov

[2]http://www.novocraft.com/products/novoalign

[3]https://www.ncbi.nlm.nih.gov/grc/human

[4]https://gnomad.broadinstitute.org/faq

[5]https://www.ncbi.nlm.nih.gov/clinvar

[6]https://www.gencodegenes.org/human/release_31lift37.html

附录1

附录2

chr1 36768200rs1573020

chr1 159174683rs2814778

chr1 204790977rs2065160

chr2 7149155rs896788

chr2 109513601rs3827760

chr2 136616754rs182549

chr3 168645035rs1498444

chr4 38803255rs4540055

chr4 159181963rs2026721

chr5 33951693rs16891982

chr7 4457003rs917118

chr10 17064992rs7897550

chr10 34755348rs1978806

chr11 32424389rs5030240

chr12 29369871rs10843344

chr12 56603834rs773658

chr13 20901724rs1335873

chr13 22374700rs1886510

chr13 34864240rs2065982

chr14 36170607rs10141763

chr14 101142890rs730570

chr15 28365618rs12913832

chr15 48426484rs1426654

chr16 31079371rs881929

chr16 90105333rs3785181

chr17 75551667rs2304925

chr18 75432386rs1024116

chr19 42410331rs2303798

chr20 38849642rs1321333

chr21 16685598 rs722098

chr21 17710424 rs239031

chr21 25672460 rs2572307

chr22 26350103 rs5997008

chr22 47836412 rs2040411

Claims (16)

1.一种对人类核酸样品进行基因组和/或基因分析的方法，包括以下步骤：

a)提供一组人类参照基因组；

b)测试所述人类核酸样品的性别和/或血统；

c)基于步骤b)中的所述性别和/或血统测试的结果，从所述一组人类参照基因组中选择一个或多个种群特异性的人类参照基因组(PHREG)；以及

d)将所述人类核酸样品与选择的PHREG进行比对。

2.根据权利要求1所述的方法，其中在多数等位基因水平或非罕见等位基因水平上进行所述比对。

3.根据权利要求1或2所述的方法，包括额外的步骤：

e)对与选择的PHREG比对的人类核酸样品进行变体识别。

4.根据权利要求3所述的方法，其中在多数等位基因水平或非罕见等位基因水平上进行所述变体识别。

5.根据前述权利要求中任一项所述的方法，其中步骤a)中提供的所述人类参照基因组是公开的人类参照基因组或衍生自公开的人类参照基因组。

6.根据前述权利要求中任一项所述的方法，其中步骤a)包括将所述人类参照基因组调整至编码水平，所述编码水平包括独特的核苷酸代码或模糊的核苷酸代码。

7.根据前述权利要求中任一项所述的方法，其中步骤a)中提供的所述人类参照基因组是PHREG。

8.根据前述权利要求中任一项所述的方法，其中所述性别测试包括以下一个或多个步骤：

测试X染色体和/或Y染色体上的性别特异性基因中的至少一个位置；利用X染色体和/或Y染色体上的人类基因组样品的比对差异；细胞遗传学测试；FISH分析；CGH分析。

9.根据前述权利要求中任一项所述的方法，其中所述血统测试基于在人类核酸样品上使用的机器学习算法，或基于利用血统特异性变体的另一个分类方案。

10.根据前述权利要求中任一项所述的方法，其中所述血统测试包括使用至少一个基因组位置的基因型，和/或测试SNP阵列或SNP芯片，和/或测试来自Sanger测序或质谱的标记物。

11.根据前述权利要求中任一项所述的方法，其中所述血统测试包括测试选自下组基因的至少一个基因：ABL2、ATP1A3、CIC、CYP2C8、CYP2C9、EPHA3、EPHA7、ERBB3、ERG、ETV1、F2、FAS、HFE、IL11RA、IL2RA、ITGB6、KIF11、KIT、KLK3、LRP6、MDM4、NAT2、NTRK2、PDGFB、PIK3R1、PLA2G3、PLAU、PRKCB、RICTOR、SLC7A11、STAT3、T、TSC1、VCAM1、VDR、VEGFB、ACVRL1、AXL、CA9、CALCR、CASP9、ENG、EPHB1、ERBB4、ESR1、FGFR2、HPSE、HSP90AA1、ITK、MRE11A、PLK1、PTPRC、SERPINE1、SMC4、TERT、TLR3、WISP3、WT1、XRCC1、ANGPT2、ARID2、BARD1、CBR3、CDH2、CYP1B1、DDR2、DNMT3A、EPCAM、ERCC2、FANCG、FANCL、GSTP1、IRS2、ITGB1、JAK3、LHCGR、MSH6、NCF2、RNF43、SLC5A5、TMPRSS2、TNFRSF8、AKT1、CD248、CD4、ESR2、EZH2、IGF1R、ITGAV、ITGB2、KLHL6、MAP3K1、MET、MLL、MTHFR、NFKB1、NUP93、PARP8、RB1、RPE65、TSHR、ABL1、BLM、CYP19A1、DPP4、EPHA6、ERBB2、EWSR1、FOXP4、ITGAM、KDM5A、LPA、LTK、MLH1、PBRM1、PHLPP2、SF3B1、TNFRSF10A、ABCG2、ACPP、ADAM15、DPYD、EPHA5、EPHB6、FOLH1、KDR、MSH3、MST1R、NTRK1、ROCK2、SLC6A2、TET2、TGM2、TH、ABCB1、CD22、CD40、CD44、CDH20、CYP11B2、ERCC5、GPR124、IL7R、ITGB3、ITGB5、NCL、NOD2、NR4A1、PGR、PLCG1、PPP2R1A、PRAME、PTCH2、RET、SETD2、XPC、ASXL1、EPHB4、PLA2G6、SYK、TET1、EP300、FLT1、ITGA1、LOXL2、PDGFRB、PIK3CD、SSTR5、TEC、APC、ATR、CLU、CREBBP、CYP2D6、EML4、MMP2、PARP2、PDGFRA、TRPM8、CSF1R、DOT1L、FGFR3、FGFR4、GLP2R、IKBKE、JAK1、NOTCH2、SPEN、SPG7、BRCA1、CYP11B1、GNAS、ITGA5、LTF、NRP2、PTK2B、TNKS、ABCC1、CEACAM5、CYP4B1、EGFR、FLT3、INSR、PTCH1、SMARCA4、ZNF217、BCR、EEF2、SELP、SLCO1B1、ABCC2、FLT4、MTR、IL4R、MTOR、RPTOR、TEK、ATM、CARD11、FANCD2、MEFV、NF1、TP73、BRCA2、CD109、PTPRD、ABCC6、IGF2R、P2RX7、ROS1、ACE、PARP1、PRKDC、CENPE、TSC2、ALK、NOTCH1、TNC、NOTCH3、POLE、MLL2、MYH11、POLD1、GRIN3B、F5、FANCA、LRP1B、LRP2、VWF。

12.根据前述权利要求中任一项所述的方法，其中所述人类核酸样品包含一组来自下一代测序(NGS)程序的读长，并且其中所述比对包括将所述读长映射至选择的PHREG的步骤。

13.一种对人类核酸样品进行基因组和/或基因分析的计算机系统，所述计算机系统包括：

a)第一模块，所述第一模块包含用于提供一组人类参照基因组的计算机指令；

b)第二模块，所述第二模块用于测试人类核酸样品的性别和/或血统；

c)第三模块，所述第三模块包含用于基于所述性别和/或血统测试的结果，从所述一组人类参照基因组中选择一个或多个种群特异性的人类参照基因组(PHREG)的计算机指令；以及

d)第四模块，所述第四模块包含用于将所述人类核酸样品与所述确定的PHREG进行比对的计算机指令。

14.计算机程序，其包含下述指令：当计算机执行所述程序时，所述程序使得所述计算机执行根据权利要求1至12所述的任何一种方法的步骤a)、b)、c)和d)。

15.计算机可读存储介质，其包含下述指令：当计算机执行所述指令时，所述指令使得所述计算机执行根据权利要求1至12所述的任何一种方法的步骤a)、b)、c)和d)。

16.一种治疗患者的方法，包括：

-获取所述患者的疾病指征的识别信息；

-获得所述患者的核酸样品；

-基于权利要求1所述的方法对所述核酸样品进行基因组和/或基因分析；

-获取所述患者的所述疾病指征的可能疗法；

-进行变体识别和解释；

-基于所述变体解释对获取的可能疗法进行分类，其中疗法被分类为对所述患者适用或对所述患者禁用；

-选择一种适用的疗法；

-根据选择的疗法对所述患者进行治疗。

Images