【3.3.1】聚类蛋白质序列(cd-hit,uclust,kclust,usearch,MMSEQ,LINCLUST,LAST)

蛋白质是促进细胞内大多数生物学过程的关键分子。因此,对它们的发现,注释和表征非常重要。在System Biology中,为了检测同源性,直系同源性,家族,共同结构域或功能相似性而大规模地按序列进行蛋白质聚类正成为一个巨大的挑战,尤其是当生活在-Omics时代时,所产生的序列呈指数增长毋庸置疑。尽管如今已有大量具有不同强度的应用程序可以满足此目的,但仍需要陡峭的学习曲线来熟悉此类工具。当用户在进行任何分析之前迷失自述文件时,通常会退出。为了帮助社区克服这种困扰,本文介绍了将蛋白质聚类为组或家族的工具和方法,并重点介绍了可以在简单的Unix终端中执行的基本命令。值得注意的是,描述了基于图的方法和基于序列的方法。

一、前言

随着每天出现新的和假设的蛋白质,蛋白质格局不断变化。今天,IMG拥有55,482个细菌基因组,1,580个古细菌,258个真核生物,1,222个质粒,7,521个病毒,1,196个基因组片段以及14,265个私有和公共鉴定的基因组和元转录组。用非常近似的估计,这相当于约70百万个非冗余蛋白质,在分离物方面具有100%的相似性,而约30亿个非冗余蛋白质对于满足基因组/ metrantranscriptome方面(来自长度约500的支架)。 UniProtKB / TrEMBL2的2017年2月15日发行版包含77,483,538个序列条目。该数字对应于1,465,039(2%)的古细菌蛋白,49,717,238(64%)细菌蛋白,22,299,253(29%)真核蛋白,2,918,867(4%)病毒蛋白和1,083,141(<1%)等。此外,Uniparc包含148,791,725个蛋白质条目。 UniProt存档(UniParc)是一个全面且非冗余的数据库,其中包含世界上大多数可公开获得的蛋白质序列。蛋白质家族的特征可以是具有相似的序列相似性的分子。值得注意的是,这个生物学问题很难解决,而且对于包含大量蛋白质结构域的真核蛋白质,大多数可用的聚类技术都失败了。然而,仍在继续努力检测最佳和更准确的蛋白质聚类仍然是一个非常活跃的研究领域。例如,PFAM版本31.0是一个包含大量蛋白质家族的数据库,它通过相似的域来组织蛋白质家族,并包含16,712个条目。如今,有几种工具遵循各种方法和策略来进行蛋白质聚类。出色的工具(例如CD-HID, UCLUST, kClust和新开发的MMSEQ / LinClust)采用了基于k-mer和基于动态编程的序列比对方法,而工具,例如MCL聚类算法和其他基于网络拓扑的聚类。在第二种情况下,在聚类之前,需要成对相似性矩阵。尽管可以通过多种方式计算出这种相似性,但使用最广泛的是BLAST + 和LAST。在本文中,为了鼓励用户熟悉几种工具并避免进行故障排除,提供了用于执行此类分析的简单命令行。

二、蛋白序列聚类的常用工具

2.1 基于序列的聚类 Sequence-based clustering

2.1.1 CD-HIT

它将蛋白质聚类为满足用户定义的相似性阈值的聚类。 每个簇具有一个代表性序列。 输入是fasta格式的蛋白质数据集,输出产生两个文件:代表序列的fasta文件(质心,centroids)和带有簇列表的文本文件(.clstr)。 在示例中将蛋白质聚类的基本命令。 相似度低至50%的fasta文件为:

cd-hit -example. fasta -o clusters -c 0.5 -n 2

2.1.2 UCLUST

给定用户定义的身份阈值,UCLUST算法将一组序列划分为簇,在将蛋白的相似度降低至50%时更有效。 这些步骤是:

根据长度来排序序列

uclust --sort example. fasta --output seqs_sorted. fasta

Cluster de novo at 50% identity

uclust --input seqs_sorted.fasta --ucresults. uc --id 0.5

Like in CD-HIT, clusters can be stored in a .clstr file after reformatting the result. uc file:

uclust --uc2clstr results.uc --output clusters. clstr.

2.1.3 KCLUST

这是一种在几天内将大型蛋白质序列数据库(例如UniProt)聚类的方法。 它可以将蛋白质成簇降低至最大成对序列同一性的20%-30%。 例如,要对一组蛋白质进行聚类以使同一性降低至50%,基本命令是:

kClust -iexample.fasta -d tmp –s 0.5。

KCLUST将创建一个/ tmp文件夹,其中包含集群结果。 headers.dmp文件将包含每种蛋白质的索引,而clusters.dmp文件将以两列格式显示聚类结果。 第一列包含每个聚类的质心序列索引,第二列包含该聚类的每个成员的索引。

2.1.4 USEARCH

这是一个独特的序列分析工具,在全球拥有成千上万的用户。 USEARCH提供的搜索和聚类算法通常比BLAST快几个数量级。 运行USEARCH的典型命令是:

usearch -cluster_faste xample.fasta -id 0.5 -centroids out.fasta –ucclusters.uc

在这种情况下,示例。 fasta是输入文件,而out是输入文件。 Fasta包含每个群集的质心。 uc最终集群。

2.1.5 MMSEQ / LINCLUST

它声称能够在两天内将数十亿个蛋白质聚类至50%的序列同一性。 它以线性时间运行,并且已经针对UCLUST和CD-HIT进行了测试。 基本命令是:

mmseqs createdb example.fasta DB
mkdir tmp
mmseqs cluster/linclust DB clu tmp --min-seq-id 0.5 --target-cov 0.5
mmseqs createtsv DB DB clu clu.tsv

所有结果都存储在clu.tsv文件中。 第一列代表每个群集的质心,而第二列代表该群集的每个成员。 因此,通过使用命令:

cut -f1 clu.tsv | wc -l

一个可以计算出簇的数量。

2.2 基于图的聚类 Graph-based clustering

将蛋白质聚类成组的另一种方法是利用以二进制形式包含所有成对相似性的构建网络的拓扑特征。 为此,建议使用BLAST +和LAST工具。 BLAST +不能扩展到庞大的数据集,而LAST和mega blast工具可以扩展。 虽然BLAST +更敏感。 对于较小的数据集和更高质量的结果(以较慢的速度为代价),建议使用Smith-Waterman算法的ssearch36或glsearch36实现。

2.2.1 相似之处 Similarities

上一个:步骤1:建立数据库。 假设有一个快速的文件叫做 example. fasta。 首先用户需要使用以下命令构建数据库(DB)

Last db -p -C 2 DB example. fasta

步骤2:必须创建逐对相似性矩阵,并将其用作基于图的聚类的输入。 典型的表示形式是制表符分隔的格式,例如:

Protein A Protein B 0.3
Protein A Protein C 0.2
Protein B Protein D 0.3…..

通过使用last,为了创建这样的矩阵,其中包含有关蛋白质的信息,这些信息的同一性低至50%,该命令必须类似于:

lastal –m200 -pBLOSUM62 -P 0 -f blasttab DB example.fasta | awk '{if($1!~/^#/ && $3>=50) print $1"\t"$2"\t"($3/100)}' >hits.list

hits.list将是一个三列制表符分隔的文件(第一列:源,第二列:目标,第三列:identity)。 创建邻接相似度矩阵的另一种方法是使用BLASTP而不是LAST。

步骤1:用户必须使用以下命令创建blast数据库:

makeblastdb -in example.fasta -dbtype 'prot' -out DB

步骤2:可以通过查询DB数据库来创建二进制邻接表,如下所示:

blastp -query example.fasta -db DB –out fmt '6 std qlenslen' -out output_tmp.list -evalue 1.0e-05

匹配数首先以称为output_tmp.list的BLASTAB文件格式存储。

步骤3:类似于上一个示例,用户可以创建相似度矩阵,该矩阵将用作图聚类输入,例如:

awk '{if($1!~/^#/ && $3>=50) print $1"\t"$2"\t"($3/100)}' output_tmp.list>hits.list

请注意,如果只想考虑最佳匹配,则他/她可能希望首先按照最佳相似性来过滤文件,例如:

sort -nrk 3 hits.list | sort -k 1,2 –u >best_hits.list

2.2.2 网络集群 Network clustering

尽管存在许多图聚类算法,但MCL算法是迄今为止最广泛用于蛋白质聚类的算法。值得注意的是,当涉及到数百万个序列时,MCL对通货膨胀值敏感并且经常对内存贪婪。要运行通胀值为1.8的MCL,典型命令是:

mcl hits.list –I 1.8 --abc -o clusters.txt

群集存储在clusters.txt中。每行都是一个群集,并且每行中的成员都是制表符分隔的。因此,线的数量对应于簇的数量。例如,要过滤单例并保持群集具有两个或多个成员,用户可以过滤为:

awk'{if(NF> = 2)print}'clusters.txt> clusters_after_threshold.txt

与MCL相似,可以像以下那样利用高度可扩展的基于SPICi密度的聚类算法:

spici -ihits.list -o out.spici

密度(Density )默认设置为0.5。要使集群的成员数大于数字,必须使用–s参数。

2.2.3 簇质心 Cluster centroids

对于基于网络的方法,群集的潜在质心可以由群集内某个相似度以上的命中数确定。

三、总结

在本文中,介绍了来自各种工具的非常基本的命令来执行蛋白质序列聚类。 尽管本指南可能对许多用户有很大帮助,但是不应盲目使用这些命令,因为所有工具对参数设置都敏感,结果可能会有很大差异。 USEARCH,UCLUST和CD-HIT应该执行几乎相同的算法。 如果全局比对一致性> 0.5%,则它们会将较短的序列聚类到更长的序列。 但是,kClust,MMSeqs和LAST以本地方式计算序列同一性。 为了直接相互比较聚类,可以使用Rand Index,Variation of Information,F-Score和其他指标。 此类指标在其他地方进行了广泛的解释。

参考资料

药企,独角兽,苏州。团队长期招人,感兴趣的都可以发邮件聊聊:tiehan@sina.cn
个人公众号,比较懒,很少更新,可以在上面提问题,如果回复不及时,可发邮件给我: tiehan@sina.cn