CN102884203A - 用于对查询序列的基因型与亚型进行分类的方法 - Google Patents
用于对查询序列的基因型与亚型进行分类的方法 Download PDFInfo
- Publication number
- CN102884203A CN102884203A CN2010800664360A CN201080066436A CN102884203A CN 102884203 A CN102884203 A CN 102884203A CN 2010800664360 A CN2010800664360 A CN 2010800664360A CN 201080066436 A CN201080066436 A CN 201080066436A CN 102884203 A CN102884203 A CN 102884203A
- Authority
- CN
- China
- Prior art keywords
- hypotype
- sequence
- sequences
- hiv
- genotype
- 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.)
- Granted
Links
Images
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
- G06F18/2132—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on discrimination criteria, e.g. discriminant analysis
-
- G—PHYSICS
- G06—COMPUTING OR CALCULATING; COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/24—Classification techniques
- G06F18/243—Classification techniques relating to the number of classes
- G06F18/2433—Single-class perspective, e.g. one-against-all classification; Novelty detection; Outlier detection
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
- G16B30/10—Sequence alignment; Homology search
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
- G16B40/30—Unsupervised data analysis
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Evolutionary Biology (AREA)
- Bioinformatics & Computational Biology (AREA)
- Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Evolutionary Computation (AREA)
- Artificial Intelligence (AREA)
- General Health & Medical Sciences (AREA)
- Biophysics (AREA)
- Biotechnology (AREA)
- Spectroscopy & Molecular Physics (AREA)
- General Physics & Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- Software Systems (AREA)
- Bioethics (AREA)
- Databases & Information Systems (AREA)
- Epidemiology (AREA)
- Public Health (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
- Peptides Or Proteins (AREA)
Abstract
本发明涉及一种用于对查询序列的基因型与亚型进行分类的方法。更具体地,本发明针对一种用于对查询序列的基因型与亚型进行分类的方法,包括:(i)选择不同病毒的碱基序列作为参考序列,这些病毒的基因型或亚型是已知的,并且通过在所述参考序列的多重比对中计算序列之间的距离而获得一种距离矩阵;以及(ii)开发一种判别方程,该判别方程可以对这些参考序列进行分类,这是通过对通过该距离矩阵的多维定标对所述参考序列成簇而获得的聚簇执行判别分析而实现的,接着根据所述判别方程对一种查询序列的基因型与亚型进行分类。
Description
技术领域
本发明涉及一种用于对查询序列的基因型与亚型进行分类的方法。更具体地,本发明针对一种用于对查询序列的基因型与亚型进行分类的方法,包括:(i)选择不同病毒的碱基序列作为参考序列,这些病毒的基因型或亚型是已知的,并且通过在所述参考序列的多重比对中计算序列之间的距离而获得距离矩阵;以及(ii)开发一种判别方程,该判别方程可以对这些参考序列进行分类,这是通过对通过该距离矩阵的多维定标对所述参考序列成簇而获得的聚簇执行判别分析来实现的,接着根据所述判别方程对查询序列的基因型与亚型进行分类。
背景技术
在理解趋异病毒的进化方面,精确的基因分型(或分亚型)是关键。近来,公共数据库里的病毒序列的数量的迅速增长被注意到。例如,NCBI基因库(NCBI GenBank)拥有的HIV-1与HCV序列条目几乎每三年翻一番。这些病毒还显示出非常好的基因型多样性并且因此已经被分类成组,被称作基因型与亚型(Robertson等人,2000;Simmonds等人,2005)。
因此,基于这些病毒株的序列相似性对它们进行基因分型(或分亚型),在理解它们的进化、流行病学以及研发抗病毒疗法或疫苗方面已经成为最基本的步骤之一。
传统的分亚型方法包括以下:(1)最近邻法,寻找该查询序列与被称作参考的每一亚型的代表的最佳匹配;(2)系统发育方法,寻找该查询序列分支至其上的单系群。由于这些亚型原本已经被定义为单独的聚簇群,所以这些直观上合理的方法已经得以广泛使用并且对于许多案例而言十分成功。
然而,随着序列数目的渐增,观察到不能被确切地分亚型的离群值或对其而言这些方法不适宜的离群值。最近一份将这些不同的自动分亚型方法与HIV-1序列作比较的报告显示,除了亚型B与C之外,它们之中的相符性低于50%(Gifford R、de Oliveira T、Rambaut A、Myers RE、GaleCV、Dunn D、Shafer R、Vandamme AM、Kellam P、Pillay D:UKCollaborative Group on HIV Drug Resistance:Assessment of automatedgenotyping protocols as tools for surveillance of HIV-1genetic diversity.AIDS2006,20:1521-1529)。该不相符性的原因之一要归结于由于重组而引起的增加的趋异性与复杂性。还应注意到,在那些方法中,紧密关联的亚型(B与D)或分享共同起源的亚型(A和CRF01_AE)显示出较差的一致性。
本发明人认为,这一问题的根本是每一亚型的参考序列的数目太少。这些方法使用两至四种手选的参考序列。它们是由各专家在高质量的全基因组序列中仔细挑选的,是要尽量覆盖每一亚型的多样性。然而,利用每一亚型的本质上小数目的参考序列,它们不能解决亚型预测的可信性;低E值的双序列比对或高系统发育树的高引导值(bootstrap value)表明单元操作的可靠性,但是就整体而言并不必然保证一个可信的亚型分类。
对缺少统计置信测度这一问题的认识带来了STAR的引入,这是一种基于特定位点打分矩阵的统计模型的方法,该特定位点打分矩阵是从每一亚型的多重序列比对(MSA)建立的。然而,其当前的实施有一些限制:它仅适用于HIV-1氨基酸序列,以小数目的参考(总共11个亚型的141种)为基础,并且利用少于1000种序列进行了测试。
最近,已经引入了新颖的基于核苷酸组成字符串的基因分型(或分亚型)方法。它的独特在于它绕过了多重序列比对并且仍旧达到高精确度。然而,它也仅使用了42种参考序列并且已经用1156种序列进行了测试。考虑到这些病毒序列数目的爆炸式增长,这些传统方法的测试案例非常少,最多万分之一。
因此,本发明的目的是要提供一种新颖的用于对公知的查询序列的基因型或亚型进行分类的方法。关键是在试图对一种查询序列进行分类之前,评估每一亚型群的聚簇程度如何。考虑这样一个案例,其中这些参考序列大部分都被亚型很好地分开了,除了两种或更多种亚型至少部分地重叠:依赖少数参考的这些方法可能没有注意到这一问题并且可能将高分分配给一种明显的亚型。由于序列范围内的不同突变率,所以每一基因片段的系统发育动力(phylogenetic power)也可能不同。这对于相对短的部分序列来说尤为关键。换言之,如果在基因分型(或分亚型)中仅考虑序列区域的一部分,那么即使这些本应区别成簇的、具有很好特征的参考也不能被分辨出。
这些最近邻法不能评估该背景分类模型的这种有效性,因为它们仅关注查询与参考之间的比对,而不是参考与参考之间。REGA,基于树的方法之一,关注该查询是在由一组参考形成的聚簇的内部还是外部(deOliveira TDeforche K、Cassol S、Salminen M、Paraskevis D、SeebregtsC、Snoeck J、van Rensburg EJ、Wensing AM、van de Vijver DA、BoucherCA、Camacho R、Vandamme AM:An automated genotyping system foranalysis of HIV-1and other microbial sequences.Bioinformatics 2005、21:3797-3800)。然而,就本发明人所知晓的,没有工具定量地报道这样一种测量。
所以,本发明人提出一种方法,该方法基于这些参考序列之间的距离开发了这些背景分类模型,重新评估了它们对于每一查询的有效性,并且就后验概率报告了基因型(或亚型)赋值的统计显著性。
如此,本发明的方法适合于其中许多参考序列可用的案例。本发明通过将主坐标分析(PCoA)与线性判别分析(LDA)(两者是使用生物科学中普遍的应用能很好建立的统计工具)结合起来而实现这些目标。PCoA(也称为经典多维定标(MDS)),将这些序列标绘在高维主坐标空间,同时尽可能地尽力保持它们之间的距离关系。PCoA已经广泛地应用于探索序列集中的全球趋势,在系统发育分析方面对基于树的方法进行了补充。
因为亚型已经被定义为系统发育树中的不同单系类群,所以如果选择一种适当的高维,每一亚型应该在MDS空间里形成良好分离的聚簇。在此类案例中,可以发现一组将这些聚簇分开的超平面并且与这些超平面相关的查询可以得到分类。为了这一目的,本发明将LDA(一种直接的并且强大的分类方法)应用于MDS坐标并且将一种查询分配给显示出最高的关系后验概率的基因型(或亚型)。
这种概率在检测任何需要仔细检验的模糊案例时是有用的。本发明的方法通过留一法交叉验证(LOOCV)来测试这些LDA模型,该验证可以用以通过检测误分类率来估测模型有效性。由于这些序列是由坐标来表示的,因此还可以开发一种简单的措施用以检测基因型(或亚型)离群值。
本发明人实质上已经利用所有来自NCBI基因库(核苷酸)与GenPept(蛋白质)的HIV-1和HCV序列对本发明进行了测试。
披露内容
技术问题
本发明的主要目的是提供一种用于对查询序列的基因型与亚型进行分类的方法,包括:(i)选择不同病毒的碱基序列作为参考序列,这些病毒的基因型或亚型是已知的,并且通过在所述参考序列的多重对比中计算序列之间的距离而获得距离矩阵;以及(ii)开发一种判别方程,该判别方程可以对这些参考序列进行分类,这是通过对通过该距离矩阵的多维定标对所述参考序列成簇而获得的聚簇执行判别分析而实现的,接着根据所述判别方程对一种查询序列的基因型与亚型进行分类。
技术解决方案
本发明的上述主要目的可以通过提供一种用于对查询序列的基因型与亚型进行分类的方法来达到,包括:(i)选择不同病毒的碱基序列作为参考序列,这些病毒的基因型或亚型是已知的,并且通过在所述参考序列的多重对比中计算序列之间的距离而获得距离矩阵;以及(ii)开发一种判别方程,该判别方程可以对这些参考序列进行分类,这是通过对通过该距离矩阵的多维定标对所述参考序列成簇而获得的聚簇执行判别分析而实现的,接着根据所述判别方程对一种查询序列的基因型与亚型进行分类。
本发明的方法的步骤(i)可以进一步包括从所述多重比对中除去插入缺失。
另外,本发明的方法的步骤(ii)的多维定标优选地是一种主坐标分析。
此外,本发明的方法的步骤(ii)的判别分析可以选自不同的方法,比如线性判别分析、二次判别分析、最近邻距离法、支持向量机或线性分类。
有利效果
本发明的方法可以被有效地用于通过分析快速进化的病毒(比如HIV-1与HCV)的序列而对病毒的基因型或亚型进行精确分类。另外,本发明的方法对核苷酸和蛋白质(多肽)序列都适用。
而且,可以应用本发明的方法根据多态性标记(比如SNP)的距离矩阵将个别受试对象分类成群组。
附图简要说明
图1示出了根据本发明的用于对病毒的基因型(或亚型)分析进行分类的方法的示意图。这些球形表示已知被成簇为四种群簇A-D的序列,并且这些组群的分界面由隔离圆圈表示。每一群簇里的实心球形分别地表示参考序列,并且查询序列由星形表示。由于查询序列位于群簇B与D之间的分界面内,所以难以查明该查询序列的基因型(或亚型)。另一方面,可以通过最近邻法来将查询序列分配给最邻近参考序列并且这种情况发生在群簇D中。根据最邻近参考序列的距离,而不考虑已知分类方法的序列的聚簇模式,就该参考序列的选择而言,这些结果可以并不稳健(robust)。
图2示出了沿第一(V1)、第二(V2)以及第三(V3)主坐标轴的HIV-1序列的示例性MDS示意图。这些参考序列被示出为根据其亚型进行了颜色编码的小圆圈。为了清楚起见,没有对亚型F-K进行标记。该查询位于亚型B的中间(‘+’)。
图3示出了对每一基因片段而言通过MDS维数K示出的LOOCV错误率。对于(a)HIV-1核苷酸、(b)HIV-1蛋白质、(c)HCV核苷酸以及(d)HCV蛋白质序列的每一基因片段而言,参考序列的预测基因型(或亚型)的LOOCV错误率是通过使该MDS维数K从1到50进行变化来进行测量。一些显示出与众不同的较高错误率的基因片段被标记。与序列类型无关,这些错误率在k=10后都达到稳定期,这些错误率在随后的分析中被使用。
图4示出了沿基因片段的LOOCV错误率的代表性滑动窗口绘图。这些LOOCV错误率是沿(a)HIV-1env核苷酸与(b)HCV e2蛋白质序列的基因片段在滑动窗口中绘制的。对两种情况而言,该MDS维数是设置在k=10。总类表示出于图8与图9中。
图5示出了用于HIV-1“主要”分析的离群值O的密度分布。在测试的161,440个案例中,根据本发明的方法的159,261个预测与LANL亚型信息(实线)相一致,而剩下的则不一致(虚线)。图5是利用在R统计包中执行的核密度估计函数来产生的。通过O>2滤出的部分标为阴影。在过滤出很大部分的不一致案例的同时,一致性案例的丢失被最小化。
图6示出了HIV-1超变异序列的离群度值的盒形图。离群度(O)参数的盒形图是针对由先前研究(Janini M、Rogers M、Birx DR、McCutchan FE:Human immunodeficiency virus type 1DNA sequencesgenetically damaged by hypermutation are often abundant in patient peripheralbolld mononuclear cells and may be generated during near-simultaneousinfection and activation of CD4(+)T cells.J Virol2001,75(17):7973-7986;Gandhi SK、Siliciano JD、Bailey JR、Siliciano RF、Blankson JN:Role ofAPOBEC3G/F-mediated hypermutation in the control of humanimmunodeficiency virus type 1in elite suppressors.J Virol 2008,82(6):3125-3130;Land AM、Ball TB、Luo M、PilonRm、Sandstrom P、Embree JE、Wachihi C、Kimani J、Plummer FA:Human immunodeficiency virus(HIV)type1proviral hypermutation correlates with CD4count in HIV-infectedwomen from Kenya.J Virol2008,82(16):8172-8182)报道的561种无功能性与1,519种功能性序列绘制的,这些研究明确地标记出每一序列是否为“无功能性的”。
图7示出了本发明对HIV-1进行分亚型的网页服务器屏幕截图。图7(a)示出了输入屏并且图7(b)-(d)分别示出了输出的第一页到最后一页。
图8在滑动窗口中示出了针对HIV-1核苷酸(上图)与蛋白质(下图)序列((a)env、(b)gag、(c)nef、(d)pol、(e)vif、(f)vpu)的LOOCV错误率。
图9在滑动窗口中示出了针对HCV核苷酸(上图)与蛋白质(下图)序列((a)utr、(b)arfp、(c)core、(d)e1、(e)e2、(f)ns2、(g)ns3、(h)ns4a、(i)ns4b、(j)ns5a、(k)ns5b、(l)okamoto、(m)p7)的LOOCV错误率。
图10示出了针对该HIV-1“主要”分析的离群度值的柱状图与LOOCV错误率。对于基于本发明的预测与LANL一致的离群度值的分布示出了以大约1.0为中心的尖峰(a),而那些不一致的则示出了直到10.0的很长的尾巴(b)。在过滤掉低可信度的案例(离群度<2.0)之后,对于不一致性预测(d)仍留下比一致性预测(c)相对更多的具有较高错误率的案例。然而,它们的比例不大并且任何基于这些值的过滤方案都没有被执行过。
最佳模式
在下文中,将参考以下实例与附图更详细地描述本发明。这些实例与附图仅给出用于说明本发明而不在于限制本发明。
总体过程
本发明的方法通过创建该查询与参考序列的多重序列比对(MSA)来开始该过程。不像常规的方法,本发明要求大量的参考,它们应该具有高质量并且具有谨慎指定的基因型(或亚型)。洛斯阿拉莫斯国家实验室(Los Alamos National Laboratory(LANL))数据库分配HIV-1(http://www.hiv.lanl.gov/)and HCV(http://hcv.lanl.gov/)序列的这样的MSA。LANL还提供有关该MSA中每一序列的亚型信息。在2007年发布的HIV-1MSAs中包括总共3,591种核苷酸与3,478种蛋白质序列,而在HCV MSAs中总共有3,093种核苷酸与3,077种蛋白质序列。应该注意,对一些亚型而言,在该MSA中发现超过100种序列,同时有极少亚型仅包括少数参考序列。该样品大小失衡是一个严重问题,但是本发明提出一种基于全局方差(global variance)的相当具有启发性的解决方案。为了与其他方法公平比较,本发明人决定将该查询与已经可从公共数据库中得到的参考序列的MSA进行比对,而不是自己创建MSA,由而对该参考MSA表示尊重。这样做具有节省执行时间的优点,这对网络服务器应用程序很关键。对于这一步骤,使用hmmbuild、hmmcalibrate、andhmmalign(http://hmmer.janelia.org/)这套程序。在使用一种PERL脚本去除该MSA中的插入缺失之后,使用具有Jukes-Cantor修正的EMBOSS程序包(http://emboss.sourceforge.net/)的distmat来计算这些序列间的配对距离矩阵。
下一步骤是所谓的主坐标分析(PCoA),它将该距离矩阵转变为其构成与所搜索的坐标的内积相等的矩阵。通过所得到的矩阵的奇异值分解,获得直到指定的较低维的一组特征向量以及相关特征值。然后将配对欧氏距离近似于这些原始距离的那些序列的多维坐标从包括这些特征向量与特征值的简单矩阵运算中恢复。每一特征值是沿由相应特征向量定义的轴而获得的方差量,也称作主坐标(PC)。为了方便,这些特征值按降序排列并且通过采用最高的少数几个来达到维数降低。如果组内变异是忽略不计的,则最高PCs的数目或该MDS维数k应该最多是N-1,其中N是参考组的数目。然而,根据所考虑的序列区域,一种亚型可能显示出复杂的聚簇模式,分为一个以上的群簇,比如亚-亚型。因此,本发明人采用一种经验性方法,该方法针对从1至50范围的k来调查这些参考序列的交叉验证误差。这一步骤是利用R统计系统(http://www.r-project.org/)中的cmdscale来实现的(图2示出了该MDS结果的一个示意图)。然后,本发明的下一步骤是开发判别模型,这些判别模型根据他们的亚型对参考进行最佳地分类并且根据这些模型给该查询分配亚型成员(membership)。在此,可以想象应用除其他以外的不同分类方法,比如K-最近邻近法(K-NN)、支持向量机(SVM)、线性分类器。如果这种MDS步骤真正有效,则这些参考应该根据其亚型成员而被很好地聚簇,并且因此诸如线性判别分析(LDA)或二次判别分析(QDA)这些最简单的方法应该有效。这两者通过使高斯分布函数适于每一组的中心而起作用,两者之间的不同之处在于是使用全局协方差(LDA)还是使用组协方差(QDA)。由于可以预计组内偏差可能组与组之间不同,因此QDA可能更合适。然而,以上提到的样品大小失衡问题阻碍了应用QDA,因为对于一些基因型(或亚型)而言,在小量参考情况下它变得不稳定。另一方面,LDA通常应用全局协方差至所有这些亚型并且因此针对这一问题可以更稳健(robust)。尽管它不如QDA严谨,但是只要这些组偏差彼此之间不是过于不同,则这种启发式方法运行地相当好。一旦基于这些参考序列计算这些线性判别,则属于特定组的后验概率是作为从该查询至该组中心的所谓的马氏距离函数而给出的。对于该查询,之后分配后验(MAP)估计的最大值,也就是,具有最大可能性的亚型。该后验概率是通过与每一亚型的参考数目成比例的前者来进行衡量的。这一步骤是利用R统计系统(http://www.r-project.org/)中的MASS程序包的lda来实现的。
预测模型的交叉验证
这些线性判别模型的有效性是通过这些参考序列的基因型(或亚型)成员的LOOCV来进行评估。对于这些参考中的每一个而言,其基因型(或亚型)是通过从这些参考中的其余参考产生的模型来进行预测的。误分类错误率(它是误分类参考的数量与参与验证的参考总数量的比)是对该背景分类能力的一种敏感量。公共数据库中的许多病毒序列并不是全基因组,而只覆盖了一些基因或一个基因的一部分,并且因此它们的系统发育信号可以不同。因此,本发明人利用LOOCV重新评估了每一预测的分类能力。如果在针对一种给定查询的MDS空间里这些参考序列得不到很好的辨析,则在LOOCV中会很明显,导致高误分类率。
离群值检测
即使通过亚型使这些参考以低LOOCV失误率得以很好地分开,该查询序列本身还是可能异常:它可以是两种或更多种亚型的复合,位于数种亚型的中间(一种重组体情况);它可以仅接近一种亚型群簇(针对这种亚型具有接近1的P值)但是远在该群簇边界之外(一种趋异情况)。在多变量分析的领域内,习惯是通过计算自样品中心的马氏距离并且通过将其与卡方分布进行比较来检测离群值。由于该马氏距离已经结合在LDA后验概率的计算中,因此本发明人提出一种有些不同的量,即,离群度O,它是从该查询至与属于沿该方向的那种亚型的参考的最大趋异值有关的群簇中心的欧氏距离:
其中XQ、XR以及XC分别是该查询、这些参考之一以及该参考组S的中心的MDS向量。该组S包含所有属于已经将该查询分类给其的基因型(或亚型)的所有参考序列。如果O小于1.0,则该查询是很好地在该群簇内部,否则就在外部。本发明人基于此开发了一种简单的启发式过滤器:例如,可以将阈值设置在2.0以容许一些偏差。REGA还通过检查树形拓扑来执行离群值检测方案以查看该查询是在由参考序列组形成的群簇的内部还是外部。
重组体检测的套合分析(Nested Analysis)
用于表征重组体病毒株的标准过程包括沿该序列的靴扫描(bootscanning)以定位该重组点。它仅适用于长序列并且对于依赖于大样品量的工具(比如本发明的方法)而言,要实用地通过互联网而服务,它花费时间太多,除非采用具有数百CPU的群簇场(cluster farm)。与其执行靴扫描,本发明人通过以下途径解决了该重组问题:(a)对于包括多于一个基因的查询而言,逐基因预测亚型;(b)以一种包括重组参考序列的“套合”方式对该分析进行再迭代。
HIV-1与HCV包含顺序的10个基因并且因此对整个基因组序列进行逐基因分析不会花费比单个基因分析长10倍的时间。如果不同亚型被以高可信度分配给了一种查询的不同基因构成部分,则暗示了一种重组情况。对于一些重组体,断点可以发生在一个基因的中间。在此类情况中,有可能的是,分类的后验概率不是仅受一种亚型支配,但是第二个左右会具有一个不可忽略的P值。本发明人通过对具有大于0.01的P值的亚型以及相关的重组体亚型予以注意,以一种“套合”方式对该预测过程进行再迭代。例如,如果A组或G组的P值大于0.01,则这些参考包括CRF02_AG组。
网页服务器开发
已经研发出接受核苷酸序列作为一种查询并且预测该查询的每一基因片段的基因型(或亚型)的阿帕奇(Apache)网页服务器,每个HIV-1与HCV有一个网络服务器。接受氨基酸序列作为一种查询的相应蛋白质版本也已经得以开发。这些可以在http://www.muldas.org/MuLDAS/上免费取得。以PERL编写的每一CGI程序封装了已从HMMER、EMBOSS以及R的各自发布网站上下载的组件程序。由于距离矩阵的运算耗费许多运行时间,因此本发明人将该任务分割为数个(典型地是四个)计算节点,其中每一个计算节点并行地计算这些行的多个部分,并且这些结果通过主节点进行整合。在英特尔至强CPU Linux盒(Intel Xeon CPU Linux box)上,对一段1000-bp的HIV-1核苷酸序列的典型亚型预测要花费大约20秒。这些网页服务器报告该查询的MAP基因型(或亚型)以及每一亚型的后验概率(posterior P)、这些预测模型的留一法交叉验证结果、以及离群值检测结果(图7的屏幕截图)。该查询的3D示意图与前三个PC中的参考是以PNG格式给出并且描述该查询的所有PC以及这些参考的XML文件可以下载,用于随后利用GGobi(http://www.ggobi.org/)的动态互动可视化(Fig.2)。这对于可视地检查聚簇的质量以及对于确定可以导致识别出潜在的新型或重组体的离群检测结果来说尤其有用。对于HIV-1,以上描述的“套合”分析被进行再迭代并且该结果也被报告。
该网站还运行存储了HIV-1亚型与HCV基因型的预计算结果的数据服务器,这些结果是利用与这些预测服务器完全一样的方法预测的。定期地(典型地是每天)下载NCBI基因库与GenPept中HIV-1或HCV的所有新条目,并且预测它们的基因型(或亚型)并存储在这些数据库中。可以通过NCBI GI编号或主入藏号(primary accessions)检索这些结果。还以利用由诸如后验概率、LOOCV率、离群度、基因型(或亚型)、或基因片段这些系统计算的性质来查询这些条目。该检索的结果包括从LANL数据库里读取的基因型(或亚型)信息,如果有的话。
结果
本发明的方法是利用从NCBI基因库与GenPept下载的HIV-1与HCV的序列数据集进行测试的。针对还没有用作参考序列的158,834种HIV-1序列(包括8,832种重组体)以及48,720种HCV序列,从LANL网站上检索核苷酸序列的亚型信息并且将这些亚型信息用于探寻出源自该核苷酸序列的蛋白质序列的亚型信息。对于一些序列而言,这些基因型/亚型是由最初提交者给出的或由LANL分配。
这些测试数据集的基因型(或亚型)命名法
HIV-1序列被分组为M(主要(main))组、N(非主要(non-main))组、U(未经分类(unclassified))组、O(外类群(outgroup))组。多数可用的序列属于M组。由于N组与O组距离M组非常远,因此M组的亚型在包括这些远离组的MDS示意图中不能得到很好的解析。因此,本发明人集中于将M组序列分类为亚型A-D、F-H、J以及K。在M组的亚型中,有时将A与F进一步分别地分开为亚-亚型A1与A2以及F1与F2。
然而,在LANL数据库中仍有一些新序列在亚型等级上被报道。甚至对于包括在由LANL产生的MSA中的序列也是这种情况。
利用本发明针对相对短的序列解析亚-亚型要求一种仅使用相关亚型序列的“套合”分析。由于这些原因,本发明人没有试图去区别亚-亚型并且在亚型等级上对它们进行分类。M组序列的不同亚型可以重组来形成一种新株。
如果在三个以上流行病学上独立的病人中发现这些株,则称它们为流行重组形式((circulating recombinant forms)CRFs)。在这些CRF中,CRF01_AE由A与现在已灭绝的E株重组形成,并且构成一个与A亚型不同的大家族。
M组与CRF01_AE亚型已被称为“主要”亚型并且本发明的方法针对它们作为“主要”分析来进行。表1列出了按照亚型以及所有已经被LANL分类为“主要”组的测试序列的基因片段统计的分项数据(相应的蛋白质序列参考表2)。该分布远不一致,代表了研究偏差:属于亚型H、J以及K的序列稀少;特别对于诸如vif与vpr的辅助蛋白而言,非B株过于稀少,以至于不能精确评估该分类。
表1.HIV-1M组以及CRF01_AE核苷酸序列的基准测试的总结性统计
(a)过滤之前每一亚型的基因片段的数目
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 2,465 | 8,576 | 2,455 | 672 | 128 | 406 | 35 | 15 | 11 | 596 | 15,359 |
| pol | 2,897 | 50,595 | 4,590 | 1,076 | 753 | 1,035 | 75 | 58 | 8 | 4,398 | 65,485 |
| vif | 51 | 1,105 | 98 | 28 | 3 | 11 | 0 | 1 | 0 | 106 | 1,403 |
| vpr | 41 | 1,504 | 43 | 13 | 4 | 3 | 0 | 1 | 0 | 106 | 1,715 |
| 5tat | 66 | 2,160 | 199 | 45 | 4 | 3 | 0 | 1 | 0 | 108 | 2,586 |
| vpu | 44 | 1,036 | 155 | 14 | 5 | 12 | 0 | 1 | 0 | 146 | 1,413 |
| env | 7,225 | 47,551 | 5,198 | 1,853 | 415 | 912 | 104 | 47 | 19 | 2,861 | 66,185 |
| nef | 284 | 6,055 | 641 | 77 | 15 | 26 | 3 | 1 | 1 | 191 | 7,294 |
| 总计 | 13,073 | 118,582 | 13,379 | 3,778 | 1,327 | 2,408 | 217 | 125 | 39 | 8,512 | 161,440 |
| Dist.Acc. | 12,669 | 113,361 | 12,824 | 3,694 | 1,302 | 2,370 | 217 | 118 | 39 | 7,686 | 154,280 |
(b)离群度过滤(outlierness filtering)(O<2.0)之后每一亚型的基因片段数目
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 2,369 | 8,558 | 2,416 | 657 | 91 | 335 | 5 | 0 | 3 | 574 | 15,008 |
| pol | 2,596 | 49,084 | 4,476 | 1,000 | 509 | 733 | 9 | 2 | 0 | 4,160 | 62,569 |
| vif | 51 | 1,104 | 98 | 24 | 2 | 8 | 0 | 0 | 0 | 106 | 1,393 |
| vpr | 41 | 1,496 | 43 | 13 | 2 | 3 | 0 | 0 | 0 | 106 | 1,704 |
| 5tat | 63 | 2,149 | 199 | 44 | 3 | 3 | 0 | 0 | 0 | 108 | 2,569 |
| vpu | 41 | 1,031 | 155 | 14 | 3 | 12 | 0 | 0 | 0 | 146 | 1,402 |
| env | 7,016 | 47,487 | 5,146 | 1,735 | 320 | 708 | 33 | 12 | 7 | 2,823 | 65,287 |
| nef | 277 | 5,969 | 641 | 71 | 13 | 23 | 3 | 1 | 1 | 187 | 7,186 |
| 总计 | 12,454 | 116,878 | 13,174 | 3,558 | 943 | 1,825 | 50 | 15 | 11 | 8,210 | 157,118 |
| Dist.Acc. | 12,069 | 111,773 | 12,622 | 3,476 | 926 | 1,792 | 50 | 15 | 11 | 7,393 | 150,127 |
(c)过滤前与LANL金标准(LANL Gold Standard)的亚型预测一致性(%)
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 94.73 | 99.07 | 99.27 | 97.77 | 89.06 | 90.64 | 71.43 | 66.67 | 81.82 | 92.95 | 97.70 |
| pol | 88.44 | 99.15 | 99.08 | 96.19 | 98.41 | 96.91 | 96.00 | 56.90 | 75.00 | 9825 | 98.47 |
| vif | 100 | 99.37 | 98.98 | 85.71 | 100 | 27.27 | - | 0 | - | 100 | 98.50 |
| vpr | 100 | 99.73 | 95.35 | 100 | 100 | 100 | - | 100 | - | 100 | 99.65 |
| 5tat | 84.85 | 99.68 | 100 | 95.56 | 100 | 100 | - | 0 | - | 100 | 99.23 |
| vpu | 70.45 | 99.90 | 100 | 100 | 80.00 | 100 | - | 100 | - | 100 | 98.94 |
| env | 98.52 | 99.64 | 98.15 | 95.14 | 94.22 | 92.76 | 50.00 | 51.06 | 57.89 | 99.62 | 99.02 |
| nef | 97.89 | 98.35 | 99.69 | 93.51 | 86.67 | 96.15 | 100.00 | 0 | 100 | 99.48 | 98.38 |
| 总计 | 95.40 | 99.32 | 98.80 | 95.84 | 96.01 | 93.98 | 70.05 | 55.20 | 69.23 | 98.46 | 98.65 |
(d)离群度过滤(O<2.0)后与LANL金标准的亚型预测一致性(%)
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 97.89 | 99.15 | 99.54 | 98.02 | 95.60 | 89.55 | 60.00 | - | 66.67 | 93.03 | 98.47 |
| pol | 93.45 | 99.33 | 99.55 | 97.50 | 99.21 | 99.05 | 100 | 0 | - | 98.41 | 99.00 |
| vif | 100 | 99.37 | 98.98 | 100 | 100 | 37.50 | - | - | - | 100 | 99.07 |
| vpr | 100 | 99.73 | 95.35 | 100 | 100 | 100 | - | - | - | 100 | 99.65 |
| 5tat | 85.71 | 99.77 | 100 | 95.45 | 100 | 100 | - | - | - | 100 | 99.38 |
| vpu | 75.61 | 100 | 100 | 100 | 66.67 | 100 | - | - | - | 100 | 99.22 |
| env | 99.26 | 99.67 | 98.70 | 98.73 | 94.06 | 95.20 | 27.27 | 0 | 42.86 | 99.65 | 99.39 |
| nef | 99.64 | 99.56 | 99.69 | 94.37 | 92.31 | 100 | 100 | 0 | 100 | 100 | 99.51 |
| 总计 | 97.66 | 99.49 | 99.22 | 98.15 | 96.92 | 95.56 | 48.00 | 0 | 54.55 | 98.59 | 99.15 |
(e)混淆表(Confusion table)(左边是LANL,顶部是本发明)
| LANL | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| A | 12,472 | 24 | 38 | 80 | 7 | 97 | 71 | 28 | 8 | 248 | 13,073 |
| B | 75 | 117,780 | 205 | 198 | 96 | 25 | 7 | 11 | 28 | 157 | 118,582 |
| C | 24 | 70 | 13,218 | 7 | 19 | 6 | 14 | 8 | 12 | 1 | 13,379 |
| D | 34 | 25 | 5 | 3,621 | 34 | 3 | 8 | 6 | 39 | 3 | 3,778 |
| F | 4 | 16 | 1 | 2 | 1,274 | 7 | 1 | 2 | 19 | 1 | 1,327 |
| G | 43 | 34 | 4 | 1 | 17 | 2,263 | 25 | 10 | 9 | 2 | 2,408 |
| H | 16 | 0 | 1 | 1 | 13 | 18 | 152 | 15 | 1 | 0 | 217 |
| J | 3 | 0 | 1 | 6 | 17 | 6 | 11 | 69 | 11 | 1 | 125 |
| K | 0 | 0 | 0 | 0 | 9 | 0 | 3 | 0 | 27 | 0 | 39 |
| 01AE | 48 | 56 | 16 | 7 | 2 | NA | 1 | NA | 1 | 8,381 | 8,512 |
| 总计 | 12,719 | 118,005 | 13,489 | 3,923 | 1,488 | 2,425 | 293 | 149 | 155 | 8,794 | 161,440 |
(f)离群度过滤(O<2.0)之后的混淆表(左边是LANL,顶部是本发明)
| LANL | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| A | 12,162 | 20 | 38 | 67 | 1 | 48 | 0 | 0 | 0 | 118 | 12,454 |
| B | 68 | 116,278 | 204 | 147 | 53 | 6 | 1 | 0 | 0 | 121 | 116,878 |
| C | 22 | 63 | 13,071 | 3 | 12 | 1 | 1 | 0 | 0 | 1 | 13,174 |
| D | 30 | 22 | 4 | 3,492 | 5 | 2 | 1 | 0 | 0 | 2 | 3,558 |
| F | 3 | 16 | 1 | 2 | 914 | 5 | 1 | 0 | 1 | 0 | 943 |
| G | 35 | 33 | 4 | 0 | 9 | 1,744 | 0 | 0 | 0 | 0 | 1,825 |
| H | 14 | 0 | 1 | 1 | 4 | 6 | 24 | 0 | 0 | 0 | 50 |
| J | 3 | 0 | 1 | 2 | 6 | 2 | 1 | 0 | 0 | 0 | 15 |
| K | 0 | 0 | 0 | 0 | 5 | 0 | 0 | 0 | 6 | 0 | 11 |
| 01AE | 44 | 52 | 15 | 5 | 0 | 0 | 0 | 0 | 0 | 8,094 | 8,210 |
| 总计 | 12,381 | 116,484 | 13,339 | 3,719 | 1,009 | 1,814 | 29 | 0 | 7 | 8,336 | 157,118 |
表2.HIV-1M组以及CRF01_AE蛋白质序列的基准测试的总结性统计
(a)过滤之前每一亚型的基因片段数目
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 2,297 | 7,616 | 2,271 | 630 | 125 | 402 | 34 | 11 | 12 | 530 | 13,928 |
| pol | 2,740 | 47,937 | 4,377 | 1,028 | 747 | 1,021 | 74 | 56 | 8 | 4,301 | 62,289 |
| vif | 48 | 825 | 88 | 26 | 3 | 3 | - | 1 | - | 114 | 1,108 |
| vpr | 40 | 1,428 | 47 | 11 | 4 | 3 | - | 1 | - | 112 | 1,646 |
| tat | 64 | 1,558 | 209 | 40 | 3 | 2 | - | 1 | - | 112 | 1,989 |
| vpu | 45 | 909 | 155 | 12 | 4 | 12 | - | 1 | - | 142 | 1,280 |
| env | 6,044 | 43,500 | 4,748 | 1,720 | 400 | 867 | 92 | 44 | 21 | 2,684 | 60,120 |
| nef | 162 | 3,981 | 406 | 66 | 14 | 22 | 2 | 1 | 1 | 178 | 4,833 |
| 总计 | 11,440 | 107,754 | 12,301 | 3,533 | 1,300 | 2,332 | 202 | 116 | 42 | 8,173 | 147,193 |
| Dist.Acc. | 11,434 | 107,622 | 12,266 | 3,531 | 1,300 | 2,331 | 202 | 116 | 42 | 8,144 | 146,988 |
(b)离群度过滤(O<2.0)之后每一亚型的基因片段数目
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 2,184 | 7,566 | 2,235 | 604 | 89 | 344 | 24 | 10 | 11 | 513 | 13,580 |
| pol | 2,370 | 44,394 | 4,162 | 942 | 425 | 776 | 31 | 17 | 3 | 3,991 | 57,111 |
| vif | 9 | 726 | - | - | - | 3 | - | - | - | - | 738 |
| vpr | 40 | 1,422 | 47 | 11 | 3 | 3 | - | 1 | - | 111 | 1,638 |
| tat | 63 | 1,552 | 209 | 40 | 3 | 1 | - | - | - | 111 | 1,979 |
| vpu | 45 | 906 | 154 | 12 | 2 | 12 | - | - | - | 140 | 1,271 |
| env | 5,766 | 43,416 | 4,690 | 1,571 | 306 | 641 | 50 | 19 | 11 | 2,597 | 59,067 |
| nef | 155 | 3,970 | 396 | 64 | 12 | 21 | 2 | 1 | 1 | 170 | 4,792 |
| 总计 | 10,632 | 103,952 | 11,893 | 3,244 | 840 | 1,801 | 107 | 48 | 26 | 7,633 | 140,176 |
| Dist.Acc. | 10,628 | 103,825 | 11,858 | 3,242 | 840 | 1,800 | 107 | 48 | 26 | 7,605 | 139,979 |
(c)过滤前与LANL金标准的亚型预测一致性(%)
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 83.54 | 93.83 | 97.4 | 86.35 | 73.6 | 75.12 | 26.47 | 0 | 16.67 | 93.21 | 91.33 |
| pol | 65.44 | 96.49 | 95.98 | 86.38 | 85.54 | 59.94 | 60.81 | 23.21 | 12.5 | 93.14 | 93.84 |
| vif | 39.58 | 96.24 | 100 | 96.15 | 100 | 100 | - | 100 | - | 100 | 94.49 |
| vpr | 95 | 94.54 | 80.85 | 45.45 | 0 | 100 | - | 0 | - | 93.75 | 93.5 |
| tat | 96.88 | 100 | 100 | 87.5 | 100 | 100 | - | 100 | - | 99.11 | 99.6 |
| vpu | 93.33 | 98.68 | 95.48 | 91.67 | 50 | 100 | - | 100 | - | 95.77 | 97.58 |
| env | 95.09 | 99.47 | 98.04 | 92.56 | 75 | 84.2 | 34.78 | 20.45 | 52.38 | 99.66 | 98.17 |
| nef | 90.74 | 99.1 | 97.04 | 83.33 | 57.14 | 77.27 | 50 | 0 | 100 | 97.75 | 98.12 |
| 总计 | 85.38 | 97.64 | 97.11 | 89.3 | 80.54 | 72.08 | 43.07 | 21.55 | 35.71 | 95.62 | 95.62 |
(d)离群度过滤(O<2.0)后与LANL金标准的亚型预测一致性(%)
| 基因 | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| gag | 87.32 | 94.14 | 98.12 | 86.26 | 6629 | 72.67 | 0 | 0 | 9.09 | 93.18 | 92.28 |
| pol | 69.66 | 97.49 | 97.43 | 89.17 | 88.24 | 56.44 | 32.26 | 0 | 0 | 93.89 | 95.25 |
| vif | 0 | 100 | - | - | - | 100 | - | - | - | - | 98.78 |
| vpr | 95 | 94.87 | 80.85 | 45.45 | 0 | 100 | - | 0 | - | 94.59 | 93.89 |
| tat | 98.41 | 100 | 100 | 87.5 | 100 | 100 | - | - | - | 100 | 99.7 |
| vpu | 93.33 | 98.79 | 96.1 | 91.67 | 50 | 100 | - | - | - | 95.71 | 97.8 |
| env | 96.39 | 99.52 | 98.7 | 95.1 | 72.55 | 84.24 | 10 | 0 | 45.45 | 99.73 | 98.62 |
| nef | 93.55 | 99.27 | 99.49 | 82.81 | 50 | 80.95 | 50 | 0 | 100 | 98.24 | 98.6 |
| 总计 | 88.44 | 98.2 | 98.09 | 91.21 | 79.29 | 70.18 | 14.95 | 0 | 26.92 | 96.06 | 96.59 |
(e)过滤之前的混淆表(左边是LANL,顶部是本发明)
| LANL | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| A | 9,767 | 95 | 248 | 105 | 98 | 250 | 104 | 26 | 34 | 713 | 11,440 |
| B | 100 | 105,213 | 603 | 1,056 | 264 | 44 | 77 | 54 | 152 | 191 | 107,754 |
| C | 88 | 83 | 11,945 | 26 | 41 | 40 | 14 | 17 | 9 | 38 | 12,301 |
| D | 36 | 182 | 51 | 3,155 | 17 | 10 | 10 | 45 | 12 | 15 | 3,533 |
| F | 19 | 76 | 51 | 8 | 1,047 | 37 | 5 | 2 | 48 | 7 | 1,300 |
| G | 289 | 54 | 100 | 15 | 74 | 1,681 | 23 | 10 | 16 | 70 | 2,332 |
| H | 57 | 12 | 9 | 6 | 13 | 15 | 87 | 1 | 1 | 1 | 202 |
| J | 14 | 8 | 22 | 5 | 10 | 16 | 12 | 25 | 4 | - | 116 |
| K | 1 | 1 | 7 | - | 10 | 7 | - | - | 15 | 1 | 42 |
| 01AE | 199 | 53 | 39 | 1 | 6 | 39 | 13 | 5 | 3 | 7,815 | 8,173 |
| 总计 | 10,570 | 105,777 | 13,075 | 4,377 | 1,580 | 2,139 | 345 | 185 | 294 | 8,851 | 147,193 |
(f)离群度过滤(O<2.0)之后的混淆表(左边是LANL,顶部是本发明)
| LANL | A | B | C | D | F | G | H | J | K | 01AE | 总计 |
| A | 9,403 | 90 | 228 | 96 | 67 | 142 | 5 | - | 4 | 597 | 10,632 |
| B | 78 | 102,076 | 580 | 883 | 149 | 25 | 11 | - | 4 | 146 | 103,952 |
| C | 79 | 79 | 11,666 | 21 | 14 | 18 | 1 | - | - | 15 | 11,893 |
| D | 30 | 173 | 49 | 2,959 | 9 | 8 | 2 | - | 1 | 13 | 3,244 |
| F | 17 | 72 | 48 | 6 | 666 | 24 | 2 | - | - | 5 | 840 |
| G | 271 | 54 | 95 | 14 | 40 | 1,264 | - | - | - | 63 | 1,801 |
| H | 51 | 12 | 9 | 5 | 3 | 10 | 16 | - | - | 1 | 107 |
| J | 13 | 7 | 12 | 5 | 5 | 6 | - | - | - | - | 48 |
| K | - | 1 | 7 | - | 6 | 4 | - | - | 7 | 1 | 26 |
| 01AE | 183 | 51 | 38 | 1 | 6 | 20 | 1 | - | 1 | 7,332 | 7,633 |
| 总计 | 10,125 | 102,615 | 12,732 | 3,990 | 965 | 1,521 | 38 | - | 17 | 8,173 | 140,176 |
HCV序列现在已经被分类为1到6基因型并且它们的亚型以“a”到“k”为后缀(例如,1a、2k、6h等等)。从LANL网站上下载的多重序列比对仅包括每一亚型的少许参考序列,这使得难以在亚型水平上应用本发明。由于这些基因型彼此大致等距,所以本发明是在基因型水平上应用,将来自一种基因型的所有亚型归并为一个组。针对HCV核苷酸与蛋白质序列的基准测试的结果分别列于表3与表4中。
表3.针对HCV核苷酸序列的基准测试的总结性统计
(a)过滤之前每一亚型的基因片段数目
| 基因 | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| 5utr | 1,533 | 306 | 662 | 406 | 56 | 332 | 3,295 |
| arfp | 2,803 | 313 | 581 | 234 | 21 | 434 | 4,386 |
| core | 4,090 | 469 | 720 | 275 | 30 | 492 | 6,076 |
| e1 | 15,248 | 1,699 | 2,251 | 1,680 | 309 | 341 | 21,528 |
| e2 | 18,271 | 1,827 | 2,348 | 1,209 | 286 | 166 | 24,107 |
| ns2 | 1,909 | 11 | 16 | 20 | 1 | 46 | 2,003 |
| ns3 | 2,764 | 65 | 271 | 22 | 107 | 46 | 3,275 |
| ns4a | 565 | 16 | 200 | 20 | 108 | 41 | 950 |
| ns4b | 746 | 11 | 42 | 20 | 108 | 41 | 968 |
| ns5a | 5,488 | 50 | 259 | 20 | 1 | 53 | 5,871 |
| ns5b | 3,299 | 453 | 504 | 442 | 122 | 334 | 5,154 |
| okamoto | 1,992 | 392 | 330 | 379 | 119 | 174 | 3,386 |
| p7 | 2,142 | 9 | 16 | 20 | 1 | 44 | 2,232 |
| 总计 | 60,850 | 5,621 | 8,200 | 4,747 | 1,269 | 2,544 | 83,231 |
| Dist.Acc. | 35,309 | 3,239 | 5,275 | 2,,698 | 611 | 1,202 | 48,334 |
(b)离群度过滤(O<2.0)之后每一亚型的基因片段数目
| 基因 | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| 5utr | 1,340 | 234 | 392 | 16 | 0 | 265 | 2,247 |
| arfp | 2,790 | 308 | 578 | 142 | 3 | 393 | 4,214 |
| core | 4,082 | 462 | 690 | 258 | 15 | 474 | 5,981 |
| e1 | 15,227 | 1,567 | 2,250 | 1,533 | 250 | 339 | 21,166 |
| e2 | 18,119 | 1,499 | 1,934 | 798 | 222 | 157 | 22,729 |
| ns2 | 1,877 | 11 | 16 | 0 | 1 | 38 | 1,943 |
| ns3 | 2,762 | 65 | 269 | 0 | 4 | 45 | 3,145 |
| ns4a | 565 | 16 | 157 | 1 | 37 | 20 | 796 |
| ns4b | 743 | 11 | 36 | 0 | 3 | 41 | 834 |
| ns5a | 5,471 | 50 | 220 | 0 | 1 | 53 | 5,795 |
| ns5b | 3,249 | 383 | 465 | 4 | 3 | 317 | 4,421 |
| okamoto | 1,992 | 391 | 330 | 379 | 117 | 172 | 3,381 |
| p7 | 2,142 | 9 | 10 | 0 | 1 | 34 | 2,196 |
| 总计 | 60,359 | 5,006 | 7,347 | 3,131 | 657 | 2,348 | 78,848 |
| Dist.Acc. | 35,046 | 3,003 | 4,674 | 2,171 | 437 | 1,166 | 46,497 |
(c)过滤前与LANL金标准的亚型预测一致性(%)
| 基因 | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| 5utr | 94.98 | 96.73 | 90.48 | 86.45 | 98.21 | 48.8 | 88.59 |
| arfp | 99.79 | 98.08 | 99.48 | 100 | 100 | 99.54 | 99.61 |
| core | 99.49 | 97.01 | 99.72 | 98.18 | 100 | 98.98 | 99.23 |
| e1 | 99.13 | 99.12 | 99.91 | 61.37 | 44.34 | 99.41 | 95.48 |
| e2 | 97.15 | 61.9 | 99.06 | 31.84 | 19.58 | 80.72 | 90.36 |
| ns2 | 100 | 100 | 100 | 100 | 100 | 100 | 100 |
| ns3 | 100 | 98.46 | 100 | 100 | 100 | 100 | 99.97 |
| ns4a | 100 | 100 | 99 | 90 | 100 | 100 | 99.58 |
| ns4b | 100 | 100 | 92.86 | 100 | 100 | 100 | 99.69 |
| ns5a | 99.98 | 100 | 100 | 100 | 100 | 100 | 99.98 |
| ns5b | 97.3 | 99.78 | 99.4 | 99.55 | 98.36 | 99.7 | 98.1 |
| okamoto | 95.63 | 100 | 99.09 | 99.74 | 98.32 | 100 | 97.25 |
| p7 | 100 | 100 | 100 | 100 | 100 | 100 | 100 |
| 总计 | 98.47 | 86.78 | 98.74 | 67.6 | 67.93 | 91.67 | 95.27 |
(d)离群度过滤(O<2.0)后的亚型预测一致性(%)
(e)过滤之前的混淆表(左边是LANL,顶部是本发明)
| LANL | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| 1 | 59,916 | 748 | 69 | 28 | 17 | 72 | 60,850 |
| 2 | 592 | 4,878 | 58 | 28 | 13 | 52 | 5,621 |
| 3 | 24 | 33 | 8,097 | 28 | 7 | 11 | 8,200 |
| 4 | 1,260 | 63 | 183 | 3,209 | 14 | 18 | 4,747 |
| 5 | 24 | 11 | 180 | 1 | 862 | 191 | 1,269 |
| 6 | 169 | 6 | 32 | 1 | 4 | 2,332 | 2,544 |
| 总计 | 61,985 | 5,739 | 8,619 | 3,295 | 917 | 2,,676 | 83,231 |
(f)混淆表(左边是LANL,顶部是本发明)(离群度<2.0)
| LANL | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| 1 | 59,543 | 721 | 54 | 9 | - | 32 | 60,359 |
| 2 | 500 | 4,472 | 10 | 2 | 1 | 21 | 5,006 |
| 3 | 15 | 28 | 7,299 | - | - | 5 | 7,347 |
| 4 | 1,238 | 50 | 172 | 1,662 | - | 9 | 3,131 |
| 5 | 18 | 11 | 177 | - | 267 | 184 | 657 |
| 6 | 143 | 3 | 28 | - | - | 2,174 | 2,348 |
| 总计 | 61,457 | 5,285 | 7,740 | 1,673 | 268 | 2,,425 | 78,848 |
(g)利用NCBI基因分型工具对(f)中的错配进行再分析
(h)利用REGA基因分型工具对(f)中的错配进行再分析
表4.针对HCV蛋白质序列的基准测试的总结性统计
(a)过滤之前每一亚型的基因片段数目
| 基因 | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| arfp | 426 | 27 | 114 | 5 | 2 | 19 | 593 |
| core | 3,025 | 326 | 541 | 240 | 27 | 413 | 4,572 |
| e1 | 13,183 | 1,342 | 1,618 | 1,344 | 286 | 330 | 18,103 |
| e2 | 15,091 | 1,016 | 1,580 | 1,116 | 260 | 134 | 19,197 |
| ns2 | 2,100 | 411 | 286 | 468 | 37 | 131 | 3,433 |
| ns3 | 2,700 | 64 | 270 | 22 | 94 | 46 | 3,196 |
| ns4a | 535 | 16 | 199 | 20 | 95 | 41 | 906 |
| ns4b | 720 | 21 | 43 | 22 | 95 | 41 | 942 |
| ns5a | 5,036 | 51 | 261 | 20 | 1 | 53 | 5,422 |
| ns5b | 2,701 | 452 | 489 | 398 | 122 | 318 | 4,480 |
| okamoto | 2,563 | 404 | 397 | 360 | 121 | 230 | 4,075 |
| p7 | 1,398 | 602 | 17 | 42 | 1 | 45 | 2,105 |
| 总计 | 49,478 | 4,732 | 5,815 | 4,057 | 1,141 | 1,801 | 67,024 |
| Dist.Acc. | 28,052 | 2,565 | 3,521 | 1,924 | 520 | 1,054 | 37,636 |
(b)离群度过滤(O<2.0)之后每一亚型的基因片段数目
| 基因 | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| arfp | 153 | 2 | 4 | 1 | NA | 5 | 165 |
| core | 2,892 | 315 | 510 | 227 | 20 | 391 | 4,355 |
| e1 | 13,030 | 1,297 | 1,605 | 1,248 | 28 | 326 | 17,534 |
| e2 | 13,840 | 860 | 889 | 624 | 90 | 107 | 16,410 |
| ns2 | 923 | 10 | 12 | NA | 1 | 38 | 984 |
| ns3 | 2,675 | 64 | 236 | NA | 22 | 35 | 3,032 |
| ns4a | 532 | 14 | 194 | NA | 41 | 30 | 811 |
| ns4b | 691 | 11 | 31 | NA | 22 | 40 | 795 |
| ns5a | 4,994 | 46 | 249 | NA | 1 | 53 | 5,343 |
| ns5b | 2,650 | 366 | 375 | 5 | 6 | 258 | 3,660 |
| okamoto | 2,452 | 401 | 397 | 354 | 119 | 230 | 3,953 |
| p7 | 1,133 | 8 | 13 | NA | 1 | 32 | 1,187 |
| 总计 | 45,965 | 3,394 | 4,515 | 2,459 | 351 | 1,545 | 58,229 |
| Dist.Acc. | 27,922 | 2,397 | 3,354 | 1,802 | 301 | 970 | 36,746 |
(c)过滤前与LANL金标准的亚型预测一致性(%)
(d)离群度过滤(O<2.0)后的亚型预测一致性(%)
| 基因 | 1 | 2 | 3 | 4 | 5 | 6 | 总计 |
| arfp | 100.00 | 100.00 | 100.00 | 100.00 | NA | 100.00 | 100.00 |
| core | 99.65 | 97.78 | 95.88 | 82.38 | 100.00 | 99.23 | 98.14 |
| e1 | 99.28 | 98.84 | 99.88 | 48.24 | 100.00 | 99.69 | 95.68 |
| e2 | 97.51 | 75.35 | 93.25 | 0.00 | 1.11 | 100.00 | 91.90 |
| ns2 | 100.00 | 100.00 | 100.00 | NA | 100.00 | 100.00 | 100.00 |
| ns3 | 100.00 | 98.44 | 100.00 | NA | 95.45 | 100.00 | 99.93 |
| ns4a | 100.00 | 100.00 | 99.48 | NA | 97.56 | 100.00 | 99.75 |
| ns4b | 96.96 | 100.00 | 90.32 | NA | 100.00 | 100.00 | 96.98 |
| ns5a | 99.98 | 100.00 | 100.00 | NA | 100.00 | 100.00 | 99.98 |
| ns5b | 98.11 | 100.00 | 99.47 | 80.00 | 83.33 | 99.61 | 98.50 |
| okamoto | 97.72 | 100.00 | 99.24 | 99.72 | 99.16 | 100.00 | 98.46 |
| p7 | 100.00 | 100.00 | 100.00 | NA | 100.00 | 100.00 | 100.00 |
| 总计 | 98.74 | 93.08 | 97.96 | 46.64 | 73.50 | 99.68 | 96.03 |
(e)过滤之前的混淆表(左边是LANL,顶部是本发明)
| LANL | 1 | 2 | 3 | 4 | 5 | 6 | 小计 |
| 1 | 47,561 | 1,129 | 509 | 16 | 147 | 116 | 49,478 |
| 2 | 680 | 3,837 | 163 | 16 | 22 | 14 | 4,732 |
| 3 | 237 | 176 | 5,306 | 77 | 13 | 6 | 5,815 |
| 4 | 1,618 | 17 | 206 | 2,206 | 8 | 2 | 4,057 |
| 5 | 102 | 143 | 29 | 179 | 673 | 15 | 1,141 |
| 6 | 9 | 62 | 34 | 6 | 4 | 1,686 | 1,801 |
| 总计 | 50,207 | 5,364 | 6,247 | 2,500 | 867 | 1,839 | 67,024 |
(f)混淆表(左边是LANL,顶部是本发明)(离群度<2.0)
| LANL | 1 | 2 | 3 | 4 | 5 | 6 | 小计 |
| 1 | 45,388 | 528 | 19 | - | 24 | 6 | 45,965 |
| 2 | 232 | 3,159 | 2 | 1 | - | - | 3,394 |
| 3 | 20 | 51 | 4,423 | 18 | - | 3 | 4,515 |
| 4 | 1,165 | 14 | 130 | 1,147 | 2 | 1 | 2,459 |
| 5 | 26 | 46 | 17 | - | 258 | 4 | 351 |
| 6 | 2 | 1 | 2 | - | - | 1,540 | 1,545 |
| 总计 | 46,833 | 3,799 | 4,593 | 1,166 | 284 | 1,554 | 58,229 |
MDS维数的确定以及模型有效性的评估
这些判别模型仅从这些参考序列建立并且因此它们的有效性与该查询序列本身大不相关。另一方面,该查询对应的基因与基因组的部分对于该区分能力很关键,这是因为该系统发育信号沿基因组而变化。对于核苷酸序列中的一个给定变异,相应的蛋白质序列可以因为负选择或正选择而显示出大不相同的变异。本发明人通过使用LOOCV错误率解决了这一问题,LOOCV错误率是通过基于参考的剩余部分对分类预测中被误分类的参考序列进行计数来测量的。
首先本发明人通过调查每个全基因片段的错误率来寻找最佳MDS维数k,然后调查具有这种k的每一基因片段的滑动窗口中的错误率。预计我们的判别模型的分类能力将通过在越来越高的k中表示这些序列而增加。通过从1到50改变k而由LOOCV运行来对该误分类错误率进行调查。
如在图3中所示出的,这些错误率迅速下降,对于k>=10达到停滞期。除了HCV 5’-UTR核苷酸、HIV-15’-tat核苷酸以及vpr/vpu蛋白质序列之外,在k>=10观察到极好的表现(错误率<5%)。已知HCV 5’-UTRs是高度保守的,难以把它们分类为各基因型。对于HIV-1而言,短基因片段一般先显示较差的表现。虽然在将k从10增加至50的计算开销(computational overhead)方面没有显著的增大,较高的k可能落入过适。所以,该整个分析使用k=10,尽管该预测网页服务器允许改变这一参数。
然后,通过测量滑动窗口(以10bp步长(step)的100bp窗口或以4aa步长的40aa窗口)中的LOOCV错误率对沿基因组或针对每一基因片段的区分能力的变化进行了测量。在图4中示出了针对HIV-1env与HCV e2的代表性示意图(参见图8与图9获得完整列表)。一般而言,沿该基因片段这些错误率相当低,尽管观察到一些明显的高峰。在HIV-1env中所见的主要高峰对应于V3环并且HCV e2蛋白质的表达谱显著地类似于测量其在1b基因型中的序列变异的熵示意图。如果该查询序列主要由这些区域组成,则该高序列变异性可能引起本发明的方法或任何其他基因分型工具的次优表现。在基于树的方法中,这将创建出由混合基因型(或亚型)组成的分支。在此类情况中,该聚簇质量的评估会不明确。
另一方面,本发明的方法提供了一些质量保证工具:LOOCV错误率,成员的后验概率,以及检查这些序列在多维空间中的分布的能力。即使对于该LOOCV错误率是大约10%的情况而言,如果发现该查询是在受其他基因型(或亚型)污染最小的多维空间的一个区域内,则该基因分型(或分亚型)可以是有效的。
性能测试
只要该查询的主要部分包含良好的系统发育信号,以上提出的问题不会是一个严重问题。这可以通过针对所有真实世界的序列运行根据本发明的基因分型(或分亚型)程序并对这些LOOCV错误率进行列表来最佳地解决。表5示出了利用HIV-1与HCV的所有非重组核苷酸与蛋白质序列的此类程序的总结。这些LOOCV错误率非常低:平均值以及中值小于5%。对于HIV-1核苷酸与HCV蛋白质序列而言,这些案例中超过90%具有小于3%的LOOCV错误率,而相应的HCV核苷酸与HIV-1蛋白质的90%百分位数(percentile)大于10%。
表5.基准测试结果的总结性统计1
1截至2009年1月20日的数据。
2仅主要分析(M组以及CRF01)。
3最大后验概率。
已经证明,根据本发明的线性模型得到了很好地验证,然后调查了分类的后验概率:这些案例中多于92%示出0.90或更高的最大后验概率值,意味着对于多数案例而言的明确响应(表5)。对于HIV-1M组以及CRF01_AE序列而言,根据本发明的预测与那些从LANL检索的预测之间的总体一致率高于95%,而对于HCV而言的那些则大于93%(表5)。
对于HIV-1与HCV核苷酸序列而言,针对每一基因与基因型(或亚型)的一致率分别列于表6与表7中(具体见在表1与表3中)。如果对于任何基因-亚型组合而言仅有少数参考序列可用,则针对该类别的亚型分类的统计模型不可靠:例如,从亚型H、J以及K每一个中仅有两至三种参考可用。这些类别中的测试序列也极其稀少(表1(a)与3(a))。除非从这些亚型中发现更多参考,否则它们的分类仍会是一个挑战。
表6.HIV-1核苷酸序列的测试结果§
§仅主要分析(M-组以及CRF01)
*以100X匹配的/总计而给出的%精确度,其中“匹配的”是本发明与LANL之间一致的案例数目。
表7.HCV核苷酸序列的测试结果
*以100X匹配的/总计而给出的%精确度,其中“匹配的”是本发明与LANL之间一致的案例数目。
离群值过滤
由于已经提议将离群度值O(Eq.1)作为该查询与相应参考之间聚簇程度的一种指标,因此检验其分布:对于这些一致性预测而言,O的密度示意图示出了一种以1.0为中心的尖锐峰(实线),而对于不一致的案例而言,观察到一种直到10.0的长尾(虚线)(图5)(相关直方图还见图10)。似乎位于O=2.0处的截止点会去除一半以上的误分类,而对一致的预测牺牲最小。由于这些总体一致率已经非常高,因此利用核苷酸序列仅看到微幅改进。然而,利用蛋白质序列则看到更多可注意到的改进(表5)。对于HIV-1gag与pol蛋白质序列(分别从91.3%与93.8%至92.2%和95.2%)以及HCV e2与ns2蛋白质序列(分别从90.6%与52.1%至91.9%和100%),看到了最显著的一致率改进。
HIV-1套合分析结果的评估
LANL已经将许多HIV-1序列描述为循环重组形式。对于来自总计8,719种此类序列的总共9,145种核苷酸基因片段,亚型是通过根据本发明的“套合”分析而得以分配。在该“主要”分析之后,对具有大于0.01的后验概率的亚型进行了收集并且将源自这些亚型的CRF加入该池中。衍生自这些亚型的CRF参考序列也被加入到该池中。
然后将根据本发明的分类流程应用于该参考池。总计5,068种核苷酸基因片段序列通过了过滤步骤(O≤2.0)并且具有明确响应(后验概率≥0.99),总体精确度93.6%(表8;针对蛋白质序列的统计结果示出在表9中)。应该注意到,目前对于CRF而言,每一基因片段或亚型的参考序列的数目不高,因此应该仔细地解释在此报道的精确度。
表8.“套合”分析HIV-1重组核苷酸序列的基准结果
表9.“套合”分析HIV-1重组蛋白质序列的基准结果
pol序列的高精确度(表8与9)令人振奋,因为在这一片段中的基因是抗病毒治疗的靶标并且近来用以帮助指导治疗方案的抗性筛选产生了日益增长数量的这类序列。
即使有这一成功,但仍然有许多序列未能通过过滤步骤。作为一种分类工具,本发明是开发在已知的亚型集中分配亚型,因此并不是被设计用来检测一种新亚型。
然而,本发明可以在离群度值以及后验概率集方面为这些离群序列的分析提示一些重要的线索。对于最终的重组体分析而言,沿该序列的靴扫描或滑动窗口分析是必要的。
亚型确定的过程
从前述中明显得出,当且仅当所报道的参数(比如后验概率(P)以及离群度(O))是合理的时,必须接受这些预测结果。对于高度可信的基因型(或亚型)分配而言,一种起作用的提议可以是P大于0.99并且O小于2.0。此标准对HCV序列的直接应用达到了大约2.6%的假阳性率,留下大约13.6%未确定。
对HIV-1序列的亚型确定不如对HCV序列的基因分型那样直接,因为前者必须要处理重组体形式的问题。本发明人示出了,本发明达到了对已经分离为无重组体或CRFs的HIV-1序列集的高分类精确性。在实际情况中,在该分析之前,该查询是否为重组体是未知的。本发明运行“主要”分析并且随后运行“套合”分析。然后需要自动化的确定过程以按有序方式对那些统计进行总结。例如,如果从这些“主要”以及“套合”分析而来的结果是不同的,则该用户会感到迷惑。在此,本发明的目的是要将该精确度最大化,而不牺牲太多预测覆盖度。基于以上提到的过滤标准,本发明提议以下策略:(i)仅当其聚簇是紧凑的(tight)(O≤2.0)并且该后验概率是大于或等于0.99时,接受从“套合”分析而来的结果;(ii)否则,仅当该“主要”以及“套合”分析彼此一致并且这些离群度值之一小于或等于2.0时,接受该结果;或者(iii)如果其离群度值小于或等于1.0并且该P值大于或等于0.99,则接受从“主要”分析而来的结果。本发明人对177,198种HIV-1核苷酸序列(基因片段)应用了这一策略,对于这些核苷酸序列而言,亚型信息从LANL可得(表10)。总共138,452种序列以98.9%的精确度通过了该第一步骤,应用于从步骤(i)剩余的38,746种序列的第二步骤产生了27,401种精确度94.6%的序列。这一启发式确定方案结果是,对于序列总数的94.2%而言总体预测精确度为98.0%,留下10,248种序列无亚型分配(5.8%)。尽管可以尝试损失精确度而使预测覆盖度最大化的可替代性策略,但本发明的方法是最小化误分类并且留下“模糊区域”给用户谨慎处理的方法。
表10.对于HIV-1核苷酸基因片段的每一亚型确定步骤的精确度与覆盖度
*在正文中,序列集2、4以及6分别对应确定步骤(i)、(ii)以及(iii)。
与其他方法比较
本发明人已经针对从LANL数据库中可得的事实上的金标准确认了本发明在对HCV序列进行基因分型以及对HIV-1序列进行分亚型方面的表现有效。由于本发明在核苷酸或蛋白质序列方面都显示出极好的表现,因此与其他自动化基因分型(或分亚型)方法比较可以提供信息。
多数常规方法报道了使用LANL的一致率与本发明的一致率相似,即使在这些方法中这些测试之一示出非常不一致的结果。然而,它们的测试案例十分有限,不如本发明的那些测试案例那样规模全。还应该强调的是,所有那些方法是基于序列比对或系统发育领域内建立好的核心算法。按照这样,那些方法的适当执行应该对该查询的分类起良好作用,只要它只使用这些基因型(或亚型)的其中一种良好地聚簇。
因此,理解这些方法在应对一种有疑问的不趋异也不重组的查询序列方面的不同之处会提供更多信息。由于没有此类公共可用的测试组,因此本发明人设计了自己的组:对其而言本发明与LANL不一致的一组序列以及本发明认为为离群值的其他序列。
本发明人利用REGA(最精密的分亚型工具之一)独立运行对该结果进行了验证(详见表11)。当本发明与REGA做出可信预测时,两者之间的一致性极好(96.6%的一致率)。REGA采用一系列检验点来拒绝有疑问的案例:例如,如果解靴带值(bootstrap value)过低或者该查询位于该亚型群簇的外部或基部,则REGA不分配亚型。
表11.对于HIV-1核苷酸序列的本发明与LANL之间的错误匹配的REGA结果
在对其而言REGA结果能够被检索的1,214个此类案例中,大约2/3未能通过REGA的那些置信检验(confidence check)。大约相同数目的案例还在藉由本发明的离群度测试中失败(O>2.0)。然而,这两种过滤结果之间的一致性不高,因为414个案例中有179个通过了本发明的过滤,但未通过REGA的过滤。尽管它们仅组成所有HIV-1基准测试案例的0.1%,但是鉴于此类案例即将出现的增长,彻底理解是重要的。这些想必落入这样的案例中,其中依赖于小数目的参考(每一亚型两至三种)的分亚型方法不能分配确定的亚型给该查询,但是与本发明相似的方法可以做出符合统计显著性的预测。
该混淆表揭示了一些情况,其中误分类案例的数目对该离群度过滤有抗性:例如,最初被LANL分配给基因型4的1,200多种序列被本发明预测为基因型1(表3(e)以及3(f))。利用NCBI基因分型工具以及REGA对这些案例的独立测试产生了使用本发明的结果(85~88%)比使用LANL的结果(11~13%)更高的一致率(表3(g)以及(h))。对这些案例的进一步的更详细分析是必要的。然而,本发明利用HCV进行了成功测试。
超突变产生远距于其祖先并且经常在成分基因(component gene)上的具有功能缺失的准种。有趣的是看看本发明对此类序列的代表性如何。有14例对其序列已经存放在公共档案中的HIV-1超突变的报告。在2,308种此类序列中,从我们的基准测试结果中鉴定出2,279种。由于具有无功能基因成分的序列可能比完整基因更趋异,因此本发明人比较源于相同研究的两组之间的趋异程度。在那14例报告中,有三例具体标记了该序列是否为无功能的。对561种无功能序列以及1,519种功能性序列的离群度O参数进行了测量。如在图6中所示,前者比后者具有明显更高的O值。这证明,本发明在预筛选超突变方面是有用的。
在下文中,描述了利用其他工具对超突变的序列进行分析的实例。贾尼尼(Janini)等人描述了一种由于超突变而编码一个无功能蛋白酶的297bp长的HIV-1pol基因(基因库(GenBank)收藏号AY036374.1以及GI:15192372)(Janini M,Rogers M,Birx DR,McCutchan FE:Humanimmunodeficiency virus type 1DNA sequences genetically damaged byhypermutation are often abundant in patient peripheral blood mononuclear cellsand may be generated during near-simultaneous infection and activation ofCD4(+)T cells.J Virol,2001,75:7973-7986)。最初的基因库记录将它分类为亚型A。尽管NCBI基因分型工具(NCBI Genotyping Tool)也将它分类为明显的亚型A,但是并没有显而易见的超突变迹象。REGA HIV分亚型工具(REGA HIV Subtyping Tool)以相当高的解靴带值(74%)将它分配给亚型A,尽管该系统发育树的拓扑形状是不正常的。本发明以高可信度(P=1.0)也将它分类给亚型A,但是比已知的亚型A参考有10多倍的不同(O组=10.56)。与所有735种参考序列所涵盖的最大半径相比,它几乎是四倍差异(O全部=3.99)。
Claims (4)
1.一种用于对查询序列的基因型与亚型进行分类的方法,包括:
(i)选择不同病毒的碱基序列作为参考序列,这些病毒的基因型或亚型是已知的,并且通过在所述参考序列的多重比对中计算序列之间的距离而获得一种距离矩阵;以及
(ii)开发一种判别方程,该判别方程可以对这些参考序列进行分类,这是通过对通过该距离矩阵的多维定标对所述参考序列成簇而获得的聚簇执行判别分析而实现的,接着根据所述判别方程对一种查询序列的基因型与亚型进行分类。
2.根据权利要求1所述的方法,其中所述步骤(i)进一步包括从所述多重比对中除去插入缺失。
3.根据权利要求1所述的方法,其中所述步骤(ii)的所述多维定标是一种主坐标分析。
4.根据权利要求1所述的方法,其中所述步骤(ii)的所述判别分析是选自包括线性判别分析、二次判别分析、最近邻点距离法、支撑向量机以及线性分类器的组。
Applications Claiming Priority (3)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| KR10-2010-0017999 | 2010-02-26 | ||
| KR1020100017999A KR20110098400A (ko) | 2010-02-26 | 2010-02-26 | 쿼리 서열의 유전형 또는 아형 분류 방법 |
| PCT/KR2010/005337 WO2011105667A1 (ko) | 2010-02-26 | 2010-08-13 | 쿼리 서열의 유전형 또는 아형 분류 방법 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN102884203A true CN102884203A (zh) | 2013-01-16 |
| CN102884203B CN102884203B (zh) | 2014-06-25 |
Family
ID=44507051
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201080066436.0A Expired - Fee Related CN102884203B (zh) | 2010-02-26 | 2010-08-13 | 用于对查询序列的基因型与亚型进行分类的方法 |
Country Status (4)
| Country | Link |
|---|---|
| US (1) | US10169532B2 (zh) |
| KR (1) | KR20110098400A (zh) |
| CN (1) | CN102884203B (zh) |
| WO (1) | WO2011105667A1 (zh) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2019047109A1 (zh) * | 2017-09-07 | 2019-03-14 | 深圳华大基因股份有限公司 | 一种hpv精确分型的生物信息学分析方法及系统 |
| CN113409886A (zh) * | 2021-06-23 | 2021-09-17 | 北京良芯生物科技发展有限公司 | 一种hiv亚型分类系统及分类方法 |
| CN117316294A (zh) * | 2023-12-01 | 2023-12-29 | 中国科学院微生物研究所 | Hiv序列分型方法、设备及存储介质 |
| CN119694392A (zh) * | 2024-12-24 | 2025-03-25 | 杭州无垠科技有限公司 | 一种基于blast比对的序列分型方法 |
Families Citing this family (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2014117296A1 (en) * | 2013-01-31 | 2014-08-07 | Hewlett-Packard Development Company, L.P. | Generating a hint for a query |
| CN104550046B (zh) * | 2013-10-16 | 2017-01-04 | 北大方正集团有限公司 | 电路板的分堆装置和电路板的分堆方法 |
| CN104156633B (zh) * | 2014-08-12 | 2017-03-01 | 上海美吉生物医药科技有限公司 | 基于rad图谱完善ssr图谱的方法 |
| CN107229839B (zh) * | 2017-05-25 | 2020-05-22 | 西安电子科技大学 | 一种基于新一代测序数据的Indel检测方法 |
| US11456057B2 (en) | 2018-03-29 | 2022-09-27 | International Business Machines Corporation | Biological sequence distance explorer system providing user visualization of genomic distance between a set of genomes in a dynamic zoomable fashion |
| CN113420779A (zh) * | 2021-05-23 | 2021-09-21 | 济南浪潮数据技术有限公司 | 一种web请求异常检测的方法、装置、设备及可读介质 |
-
2010
- 2010-02-26 KR KR1020100017999A patent/KR20110098400A/ko not_active Ceased
- 2010-08-13 WO PCT/KR2010/005337 patent/WO2011105667A1/ko not_active Ceased
- 2010-08-13 CN CN201080066436.0A patent/CN102884203B/zh not_active Expired - Fee Related
- 2010-08-13 US US13/581,338 patent/US10169532B2/en not_active Expired - Fee Related
Non-Patent Citations (4)
| Title |
|---|
| GIFFORD ET AL.: "UK Collaborative Group on HIV Drug Resisitance: Assessment of automated genotyping protocols as tools for surveillance of HIV-1 genetic diversity", 《AIDS》 * |
| MYERS ET AL.: "A statistical model for HIV-1 sequence classification using the subtype analyser(STAR)", 《BIOINFORMATICS》 * |
| WU ET AL.: "Nucleotide composition string selection in HIV-1 subtyping using whole genomes", 《BIOINFORMATICS》 * |
| ZHANG ET AL.: "jpHMM at GOBICS: a web server to detect genomic recombinations in HIV-1", 《NUCLEIC ACIDS RESEARCH》 * |
Cited By (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2019047109A1 (zh) * | 2017-09-07 | 2019-03-14 | 深圳华大基因股份有限公司 | 一种hpv精确分型的生物信息学分析方法及系统 |
| CN111032885A (zh) * | 2017-09-07 | 2020-04-17 | 深圳华大基因股份有限公司 | 一种hpv精确分型的生物信息学分析方法及系统 |
| CN111032885B (zh) * | 2017-09-07 | 2024-05-17 | 深圳华大基因股份有限公司 | 一种hpv精确分型的生物信息学分析方法及系统 |
| CN113409886A (zh) * | 2021-06-23 | 2021-09-17 | 北京良芯生物科技发展有限公司 | 一种hiv亚型分类系统及分类方法 |
| CN117316294A (zh) * | 2023-12-01 | 2023-12-29 | 中国科学院微生物研究所 | Hiv序列分型方法、设备及存储介质 |
| CN117316294B (zh) * | 2023-12-01 | 2024-02-20 | 中国科学院微生物研究所 | Hiv序列分型方法、设备及存储介质 |
| CN119694392A (zh) * | 2024-12-24 | 2025-03-25 | 杭州无垠科技有限公司 | 一种基于blast比对的序列分型方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN102884203B (zh) | 2014-06-25 |
| KR20110098400A (ko) | 2011-09-01 |
| WO2011105667A1 (ko) | 2011-09-01 |
| US20130226466A1 (en) | 2013-08-29 |
| US10169532B2 (en) | 2019-01-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN102884203B (zh) | 用于对查询序列的基因型与亚型进行分类的方法 | |
| Lan et al. | A two-phase learning-based swarm optimizer for large-scale optimization | |
| CN110853756B (zh) | 基于som神经网络和svm的食管癌风险预测方法 | |
| Li et al. | Protein contact map prediction based on ResNet and DenseNet | |
| Deng et al. | Identifying stages of kidney renal cell carcinoma by combining gene expression and DNA methylation data | |
| He et al. | Evolutionary graph clustering for protein complex identification | |
| Tanchotsrinon et al. | A high performance prediction of HPV genotypes by Chaos game representation and singular value decomposition | |
| Liu et al. | A comparison of topologically associating domain callers based on Hi-C data | |
| Hickl et al. | Binny: an automated binning algorithm to recover high-quality genomes from complex metagenomic datasets | |
| Suo et al. | Application of clustering analysis in brain gene data based on deep learning | |
| He et al. | A novel alignment-free method for HIV-1 subtype classification | |
| CN115206437A (zh) | 一种线粒体效应分子的智能筛选体系及其构建方法和应用 | |
| Parmar et al. | A review on data balancing techniques and machine learning methods | |
| KR101953651B1 (ko) | 쿼리 서열의 유전형 또는 아형 분류 방법 | |
| Yang et al. | Data-driven identification of SARS-CoV-2 subpopulations using PhenoGraph and binary-coded genomic data | |
| Li et al. | Virtual screening of drug proteins based on imbalance data mining | |
| Wang et al. | A Hardware Trojan Detection and Diagnosis Method for Gate-Level Netlists Based on Machine Learning and Graph Theory | |
| CN112951324A (zh) | 一种基于欠采样的致病同义突变预测方法 | |
| Modak et al. | Application of support vector machines in viral biology | |
| Warrens et al. | Understanding partition comparison indices based on counting object pairs | |
| Xia et al. | High-accuracy identification of incident HIV-1 infections using a sequence clustering based diversity measure | |
| Li et al. | A new ensemble coevolution system for detecting HIV-1 protein coevolution | |
| Kim et al. | A classification approach for genotyping viral sequences based on multidimensional scaling and linear discriminant analysis | |
| CN117831661A (zh) | 对预测分子性质的人工智能模型进行评价的方法 | |
| Liu et al. | CoT: a transformer-based method for inferring tumor clonal copy number substructure from scDNA-seq data |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20140625 Termination date: 20200813 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |