KeepNotes blog

Stay hungry, Stay Foolish.

0%

Cutadapt对测序数据质控（去接头）

Posted on Edited on In ,
Symbols count in article: 2.6k Reading time ≈ 2 mins.

简单的说，cutadapt是一款去接头并且也能去低质量碱基的软件。至于接头的产生，一般来说是由于在随机打断时，有些片段的长度过于短，导致插入片段不够，使得测序时会读到测序引物。因此理论上来说接头是不应该被测序仪读到的，但是个别情况读到引物，我们需要通过cutadapt进行去除。 网上多为早些时候的cutadapt教程讲解，对于双端测序的实践较少。因此做个笔记记录下cutadapt对双端测序数据的处理。

软件下载及安装

现在这个软件已经更新到1.14版本了，如果下载的是老版本，会导致部分功能无法使用（特别对于双端测序去接头来说），由于cutadapt是python写的，因此下载地址可以是：https://pypi.python.org/pypi/cutadapt#downloads，我选择下载的是cutadapt-1.14.tar.gz

tar zxvf cutadapt-1.14.tar.gz
cd cutadapt-1.14
python setup.py install --user

然后开始看cutadapt的说明书吧，蛮详细的https://cutadapt.readthedocs.io/en/stable/index.html

对双端数据进行质控

一般现在测序以双端测序为止，那么双端测序的接头是如何产生的呢？可以参考下面这篇博文http://blog.sina.com.cn/s/blog_14d1975e90102x5o1.html

先看下cutadapt是如何去接头的

cutadapt -a ADAPTER_FWD -A ADAPTER_REV -o out.1.fastq -p out.2.fastq reads.1.fastq reads.2.fastq

这里我们主要看到-a 接头，-A 接头的反向互补序列，并且-a和-A均是去除3'端的接头序列，那么以一个接头序列(NEBNext® UltraTM RNA Library Prep Kit RNA)为例：

5' Adapter(5’端接头)
5'AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT

3' Adapter(3’端接头) 5’GATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG

按照我的理解，-a去除的是read1的3'的接头序列，因此-a参数后接3' Adapter(3’端接头)；-A 去除的是read2的3'的接头序列，那么-A参数后应该接5' Adapter(5’端接头) 的反向互补序列，至于为什么是这样，可以简单的理解为由于插入片段过短，测序仪不仅测了插入片段，还测到了接头序列，侧穿了（不知道对不对，查了一些别人的解释，测序原理还不是非常清楚）

所以命令如下：

cutadapt -a GATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG \
-A AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTATCATT \
-o out.1.fq -p out.2.fq read.1.fq read.2.fq

下面几个参数跟去接头比较相关

  • -e 最大错配比例，比如cutadapt在某条序列上检测的接头有15bp长，那么允许这个匹配上的15bp接头中有15*0.1约为1个碱基的错配
  • -m --minimum-length 切除接头后的序列长度的最小值
  • -O --overlap 默认必须至少有3个碱基匹配时才会认为是adapter序列，但有时可以适当的调大
  • --discard-trimmed 去除掉有检测到接头的序列（默认cutadapt只是截掉接头序列以及接头序列以后的序列）
  • --untrimmed-output 将没有接头的序列输出到目标文件中(但是必须要跟-o 一起用)
  • --untrimmed-paired-output 将没有接头的paired序列输出到目标文件中（也要跟-p 一起用）
  • --pair-filter=(any|both) 这个参数很好用，对于双端测序而言，read1和read2都有可能检测到接头。如果选择any，则只要两个中其中一个检测到接头，read1和read2均舍弃；如果选择both，则必须两个都检测到接头，read1和read2才舍弃

因此如果我们只想要输出没有接头的序列并且paired也同样舍弃，总结上述参数，那么可以如下：

cutadapt -a GATCGGAAGAGCACACGTCTGAACTCCAGTCACATCACGATCTCGTATGCCGTCTTCTGCTTG \
-A AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTAGATCTCGGTGGTCGCCGTATCATT \
-e 0.1 -O 5 -m 50 --discard-trimmed --paired-filter=any
-o out.1.fq -p out.2.fq read.1.fq read.2.fq

cutadapt除了去接头序列，也能去除低质量碱基。其具体原理可看https://cutadapt.readthedocs.io/en/stable/guide.html#quality-trimming-algorithm，使用-p 20参数即可，默认是对3'端进行去除低质量的碱基；如果想对5'端进行处理的话，使用-p 20,20 也就是用逗号将3'端和5'端的阈值分开表示

Summary

总体上来说，Cutadapt也是一款不错的去接头软件，也有人将其与Trimmomatic的去接头效果进行比较，差别并不大；但是个人觉得Trimmomatic去低质量的算法更加给力，去的更加彻底。单从去接头的来说，确实也差不多，都可以用。Trimmomatic只能去除Illumina平台的接头，而cutadapt则可以去除任意指定的接头

参考文章：
http://www.bio-info-trainee.com/1920.html
http://blog.sina.com.cn/s/blog_14d1975e90102x5o1.html
https://cutadapt.readthedocs.io/en/stable/guide.html#

本文出自于http://www.bioinfo-scrounger.com转载请注明出处

Related Posts