0.简介

生信分析中经常要根据指定条件查找相似序列,比如构建多个样品间的非冗余基因集、分析样品间的相似程度等等,cd-hit这款软件就可以用较短的时间解决此类问题,可以对单个数据集进行去冗余,包括DNA/RNA序列和蛋白序列,也可以对两个数据集进行比较。其工作原理可概述为:将所有序列按照参数设定进行聚类,并将每一组聚类中的最长序列作为代表序列进行输出,同时给出每组聚类下的每个序列名可供相似度分析使用。下面我们来简单介绍一下它的使用方法。

1.下载与安装

网址:http://cd-hit.org ;http://www.bioinformatics.org/cd-hit/ ;https://github.com/weizhongli/cdhit/archive/V4.6.2.tar.gz;

这是一个在linux系统下使用的工作,我们可以给自己的电脑装一个双系统或者在windows下使用linux的虚拟机。然后我们可以执行下面的命令进行解压(注意我们要将路径先切换到安装包所在的文件目录下,或者在执行命令时使用完整路径)。

gzip -d cdhit-4.2.tar.gz

然后进入到解压后的文件夹(我解压后的文件夹为cdhit-4.2,同样要注意我们的文件路径问题,如果上面使用的是完整路径,最好这里也使用完整路径,比如我使用完整路径是‘cd /home/zpf/cdhit-4.2’)

cd cdhit-4.2

最后编译一下就可以了,执行make

make

然后我们就可以使用这个工具了。

2.输入文件格式

Cd-hit的输入文件仅有一个fasta格式文件 ,一般来说cd-hit是将几个样品的基因或蛋白序列进行聚类,所以需要将这些样品的序列汇总到一起作为输入文件,可在linux系统下通过cat命令实现:

cat a.fasta b.fasta c.fasta > all.fasta

其中a.fasta,b.fasta,c.fasta为fasta格式的三个样品基因或蛋白序列,all.fasta为汇总后的序列,在分析中作为cd-hit的输入序列。值得注意的是,在三个样品序列中不能有序列名相同的序列,否则会出现错误。因此,一般在分析时会在各样品序列名前添加样品名,这样即可避免重复。序列名是fasta文件中以“>”开头的行空格之前的内容,如下图中蓝色线圈出部分。

图1

3.输出文件

Cd-hit有两个输出文件:一个是只含有所有代表序列(即去冗余后的序列)的fasta文件,其格式参看图1;另一个是以.clstr结尾的聚类信息文件,其格式如图2。

图2

以“>”开头的是一个聚类组。每组下面按序号排列,如上图中Cluster 1组有5个聚类序列。每个聚类序列有一个百分比或“*”,百分比代表该序列与代表序列的相似度,“*”代表该序列即为代表序列。

4.去除冗余的基本思路

首先对所有序列按照其长度进行排序,然后从最长的序列开始,形成第一个序列类,然后依次对序列进行处理,如果新的序列与已有的序列类的代表序列的相似性在cutoff以上则把该序列加到该序列类中,否则形成新的序列类。之所以快主要是两个方面的原因:一个是使用了word过滤方法,即如果两条序列之间的相似性在80%(假设序列长度为100),那么它们至少有60个相同的长度为2的word,至少有40个相同的长度为3的word,至少有20个相同的长度为4的word。基于这个原则,在处理新的序列的时候,如果新的序列与已有序列的相同word的长度不能满足这些要求则不需要进行比对了,这极大的降低了时间消耗;另外一个速度快的原因是使用了index table,可以很快的计算序列之间相同word的数目。

#当序列相似性在80%时,有20个位点是有差异的,极端的情况就是这20个位点对应的长度为2的字符串都不一样,因此是40个不一样,当有更多的不一样时,两条序列的相似性不可能在80%;同理,如果这20个位点对应的长度为4的字符串都不一样,则有80个不一样。

5.使用方法和参数介绍

Cd-hit运行时用很多参数可以进行调整设置,其运行命令为(参数仅为示例)在刚才编译的文件路径下执行:

cd-hit -i all.fasta -o new.fa -c 0.8 -aS 0.8 -d 0

下面简单介绍一下重要的几个参数:

-i:输入文件,fasta格式。

-o:输出文件前缀,输出文件有两个,分别为fasta格式序列文件和以.clstr结尾的聚类信息文件。

-c:较短序列比对到长序列的bp与自身bp数的比值超过该数值则聚类为一组,默认为0.9。

-d:聚类信息文件中各个聚类组中序列名的长度,设为0则将取完整序列名。

-aL:控制代表序列比对严格程度的参数,默认为0,若设为0.8则表示比对区间要占到代表(长)序列的80%。

-AL:控制代表序列比对严格程度的参数,默认为99999999,若设为40则表示代表序列的非比对区间要短于40bp。

-aS:控制短序列比对严格程度的参数,默认为0,若设为0.8则表示比对区间要占到短序列的80%。

-AS:控制短序列比对严格程度的参数,默认为99999999,若设为40则表示短序列的非比对区间要短于40bp。

下图详解了-aL,-AL,-aS,-AS四个参数。

图3

aL = Ra / R

AL = R - Ra

aS = Sa / S

AS = S - Sa

6.cdhit的缺点

1 它不能保证同一个序列类中的序列的相似性都在threshold之上,因为每次比对都是用新序列与序列类的代表序列进行,这就有可能使得序列类中除了代表序列外其他序列之间的相似性在threshold之下。比如A是代表序列,B与A的相似性大于0.95,C与A的相似性也大于0.95,但是这并不能保证B与C的相似性也大于0.95.

2 它不能保证一个序列类的病毒与另外一个序列类中的病毒的相似性也在threshold之上,原因还是在于用代表序列代表了整个序列类。

3 基于word filter的方法使得使用每个长度的word能够处理的冗余性水平有限,如使用长度为2的word只能够得到相似性在50%以上的序列,长度为3的word只能够得到相似性在66.7%以上的序列类,类似的,长度为5的word只能够得到相似性在80%以上的序列。在实际应用的时候需要注意选择的word长度与threshold的匹配。

教程 | 如何用cd-hit去除冗余序列?相关推荐

  1. linux cd-hit下载安装,教程 | 如何用cd-hit去除冗余序列?

    原标题:教程 | 如何用cd-hit去除冗余序列? 生信分析中经常要根据指定条件查找相似序列,比如构建多个样品间的非冗余基因集.分析样品间的相似程度等等,cd-hit这款软件就可以用较短的时间解决此类 ...

  2. CD-HIT去除冗余序列

    1.简介 CD-HIT是用于蛋白质序列或核酸序列聚类的工具,根据序列的相似度对序列进行聚类以去除冗余的序列,一般用于构建非冗余的数据集用于后续的实验分析. 基本思路是首先对所有序列按照其长度进行排序, ...

  3. CD-Hit 生信 碱基序列去除冗余的方法

    1.CD-Hit介绍 官方介绍: CD-HIT是一个非常广泛使用的程序,用于蛋白质或核苷酸序列的聚类和比较.最初由李伟忠博士在伯纳姆研究所(现为桑福德伯纳姆医学研究所)亚当·戈兹克博士的实验室开发,C ...

  4. NAR:测序数据鉴别和去除rRNA序列利器RiboDetector

    [编者荐语]: rRNA序列污染是广泛存在于各类高通量测序数据中的.除了在实验建库过程中对文库进行去核糖体的处理,数据分析层面也可通过一些序列比对的策略去除.RiboDetector是邓志罗博士基于深 ...

  5. 去除冗余token的DETR效果怎么样?NUS颜水成团队提出端到端的PnP-DETR结构

    ​作者丨happy 编辑丨极市平台 本文原创首发于极市平台,转载请获得授权并标明出处. 原文链接:https://arXiv.org/abs/2109.10852 语言模型与目标检测这种八竿子打不着的 ...

  6. 求圈地的方块数java_[玩家教程]如何用Residence插件圈地及进行其他操作

    圈地对新手玩家可算是件大事,(领地插件下载)当然对大部分小白玩家来说理解圈地的确不太容易,而领地插件自带的命令显得过于死板,下面假设圈地工具为线,用图解来解释圈地.以下教程以服主的口气编辑,大家可以直 ...

  7. 函数调用关系图如何画_彩铅画入门植物教程 | 如何用彩铅画一株多肉?多肉彩铅画教程步骤图详细...

    画画不难,难的是不拿起手中的笔去画. 彩铅画入门植物教程 | 如何用彩铅画一株多肉?多肉彩铅画教程步骤图详细 多肉的质感如何表达呢?还是那句话:艺术来源于生活,要仔细观察.拿我们今天画的多肉来说,首先 ...

  8. Log4j2的additivity属性(是否去除冗余日志)

    log4j的additivity属性值默认是设置为true的.可参考其api,地址:http://logging.apache.org/log4j/1.2/apidocs/org/apache/log ...

  9. 保姆级教程如何用Xcode搭建python环境

    保姆级教程如何用mac电脑中的Xcode搭建python环境(xcode12) 「mac电脑自带python2.7,你也可以更新你的python版本」 打开Xcode,点击file-new-proje ...

最新文章

  1. 200万年薪,招不到清华姚班毕业生,能上姚班的都是什么人?
  2. 不能往Windows Server 2008 R2 Server中复制文件的解决方法
  3. Micropython教程之TPYBoardv102 DIY蓝牙智能小车实例
  4. Flink 与 TiDB 联合发布实时数仓最佳实践白皮书
  5. Go goroutine
  6. 二分搜索及其扩展(循环递增数组的搜索)
  7. springaop事务逻辑原理_搞懂Spring AOP,这一篇就够了
  8. Asp.NET Core2.0 项目实战入门视频课程_完整版
  9. dev可以运行mysql文件夹_Linux查看mysql 安装路径和运行路径
  10. 计算机未来发展趋势四个字概括,授课教师-世界大学城.doc
  11. java新手笔记18 类比较
  12. (转)量化投资大师采访摘录-詹姆斯·西蒙斯 James Simons
  13. LTE 中的CQI,PMI,RI上报机制
  14. 程序员风格的修真小说之炫小说
  15. 重复字符串的处理问题
  16. 计算机中心云桌面建设方案,多媒体教室云桌面方案.docx
  17. oracle销售退货业务,Oracle EBS OM RMA销售退货异常处理(Datafix)
  18. Java基础之《JVM性能调优(3)—堆》
  19. CSDN如何上传资料
  20. 请用数字灰度传感器和步进电机及A4988控制器和stm32设计一个巡线跑道的范例代码...

热门文章

  1. QT4.8的只能确认不能拒绝的告白小应用
  2. 微信支付的支付流程涉及哪些接口
  3. 学习使用jquery控制select下拉选项的字体样式
  4. Dobot magician机械臂抓取实战---机器人导论(1)
  5. PySCENIC(三):pyscenic单细胞转录因子分析可视化
  6. 几种开源SIP协议栈对比
  7. Java调用Windows系统命令CMD
  8. 2022强网拟态pwn-webheap
  9. Android的WakeLock机制
  10. 大数据存储系统学习笔记(一)