cutadapt 的安装与使用

cutadapt 是一款质量过滤的软件， 它可以删除adapter, primer. polyA尾等序列；也可以用来去除低质量序列

源代码： https://github.com/marcelm/cutadapt/

官方文档：http://cutadapt.readthedocs.org/en/stable/

安装：

          git clone https://github.com/marcelm/cutadapt/

          cd cutadapt

          python setup.py build_ext -i

          cd bin

          export PATH=$PATH:$PWD

测试：

1） 去除3端引物序列

在建库的过程中，会在DNA序列的3端加上对应的引物序列， 当测序仪的读长read 过长时，比如建库后的DNA fragement 长度是 200bp , 机器读长是250bp, 机器从序列的5端开始测序，因为读长大于序列长度，会读到adapter 序列，测序得到的reads的3端就会包含引物序列，用法：

cutadapt -a AACCGGTT -o output.fastq input.fastq

-a 指定的是3端引物序列

去除3端

如果机器读长大于fragment 的长度，会出现两种情况；

1） 在测序得到的序列的3'末端会出现部分adapter 序列；

2）在3'末端有完整的adapter 序列，在adapter 序列之后还有其他随机测到的序列；

input.fastq的内容如下：

@NB001
ATCATATCGATCTAGCCAACCGGTTACGA
+
/////////////////////////////
@NB002
ATCATATCGATCTAGCCAACCG
+
//////////////////////

3'端adapter 的序列为 AACCGGTT， 运行上述程序, output.fastq 的结果如下：

@NB001
ATCATATCGATCTAGCC
+
/////////////////
@NB002
ATCATATCGATCTAGCC
+
/////////////////

 对上述两种情况，cutadapt 会把adapter 序列都清除干净

2）去除5'端引物序列

在建库过程中，如果在5端添加adapter 序列， 也需要去除5端引物序列

 

通常情况下，5段引物序列位于序列的起始位置，但是当5端引物降解时，可能只有部分引物序列会出现在序列的起始位置；有时5段引物序列还会出现在序列中间，这时引物序列之前的序列都会被切掉

假设引物序列是ADAPTER, 测序得到的序列是 MYSEQUENCE, 5端引物可能出现的情况有以下几种

ADAPTERMYSEQUENCE                       : 正常情况

DAPTERMYSEQUENCE                        ： 部分降解

TERMYSEQUENCE                              ： 降解严重

SOMETHINGADAPTERMYSEQUENCE      : 在5段接头之前混入其他序列，接头序列被污染

 

上面所有的情况，在切除5段引物时，这些污染序列都会被切除，只剩下正常的序列， 用法：

cutadapt -g AACCGGTT -o output.fastq input.fastq

input.fastq 的内容：

@NB001
AACCGGTTATCATATCGATCTAGCC
+
/////////////////////////
@NB002
CGGTTATCATATCGATCTAGCC
+
//////////////////////
@NB003
GTTATCATATCGATCTAGCC
+
////////////////////
@NB004
ACGTAAACCGGTTATCATATCGATCTAGCC
+
//////////////////////////////

output.fastq 的内容：

@NB001
ATCATATCGATCTAGCC
+
/////////////////
@NB002
ATCATATCGATCTAGCC
+
/////////////////
@NB003
ATCATATCGATCTAGCC
+
/////////////////
@NB004
ATCATATCGATCTAGCC
+
/////////////////

 可以看到所有的adapter 序列都被清除干净了

有一种特殊的情况叫做 Anchored adapter , 就是说只有当在reads的3’或者5'端出现长度一致的序列时，才认定这是一段adapter 序列，才会进行相应的去除adapter 序列的操作；

3） 去除barcode 序列或者兼并引物的序列

cutadapt 允许采用IUPAC 的碱基字符表

            

 

字母N代表任意碱基， 在去除barcode 序列是非常有用， 比如 ACGTANNNNNTAC 其中NNNN 代表barcode 序列， 又叫做index 序列

 

cutadapt -a ACGTAANNNNTTAGC -o output.fastq input.fastq

 

兼并引物序列：YACGT； 就是说可能为 CACGT 或者 TACGT

 

cutadapt -a YACGT -o output.fastq input.fastq

 

 

5) 去除双端测序数据的adapter 序列

用法：

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

 

备注：

cutadapt 在搜寻引物时， 默认允许错配和插入缺失， 比如

adapter 序列是ADAPTER, 此时

ADABTER    有一个错配，也会被剪切掉

ADAPTR      有一个缺失，deletion ,也会被剪切掉

ADAPPTER  有一个插入，也会被剪切掉

 

可以采用-e  指定错配的比例， 默认-e 为0.1， 比如adapter  序列长度为9，允许的错配数为 9 * 0.1 = 0.9, 然后直接向下取整后为0， 所以允许的错配数为0；

 

可以采用-no-indels 来禁止插入和缺失，减少错误的剪切情况

因为cutadapt  允许部分匹配，比如 adapter 序列为ADAPTER ， 测序得到的 序列为ATCGATGCTADCGAGCGC，在序列的中间位置AD是adapter 序列的一部分， 此时会把AD以及之后的序列全部剪切掉，这种情况就属于错误的剪切，所以cutadapt 默认必须至少有3个碱基匹配时才会认为是adapter 序列，然后进行切除， 这个参数可以通过 --overlap  参数来指定默认的3适合绝大多数的情况

 

6) 去除低质量序列

使用-q/--trim-qualities 过滤低质量序列， 在去除adapter序列之前就开始过滤低质量序列；默认使用phred quality+33 的方式识别序列质量，如果为phred quality+64的方式，则添加--quality-base=64 参数

用法：

cutadapt -q 10 -o output.fastq input.fastq

 

默认只过滤3端的低质量序列， 如果想要过滤5端低质量序列，需要用逗号隔开

cutadapt -q 15,10 -o output.fastq input.fastq

 5端用15进行过滤，3端用10进行过滤

 

质量过滤的算法

cutadapt 使用的质量过滤的算法和bwa 一致，假设一段序列质量编码为

42, 40, 26, 27, 8, 7, 11, 4, 2, 3

 

质量过滤的阈值为10，则首先减去10

32, 30, 16, 17, -2, -3, 1, -6, -8, -7

 

从末端开始累加，

(70), (38), 8, -8, -25, -23, -20, -21, -15, -7

 

因为-25 最小，所以保留-25 之前的碱基， 即保留前4位碱基

 

7） 去除polyA 尾的序列

当adapter 序列中包含重复序列时，推荐使用A{10}这种写法， 代表10个A，所以可以用来去除序列末端的polyA尾，用法：

cutadapt -a "A{10}" -o output.fastq input.fastq

注意的是，当序列末尾出现10个或者10个以上的A时，都会被去除

input.fastq 的内容：

@NB001
AACCGGTTATCATATCGATCTAGCCAAAAAAAAAA
+
///////////////////////////////////
@NB002
AACCGGTTATCATATCGATCTAGCCAAAAAAAAAAAAAA
+
///////////////////////////////////////

output.fastq的内容：

@NB001
AACCGGTTATCATATCGATCTAGCC
+
/////////////////////////
@NB002
AACCGGTTATCATATCGATCTAGCC
+
/////////////////////////

 可以看到，序列末尾的polyA尾全部被去除了

8） 从序列两端切除固定长度的序列

cutadapt 还可以在序列的3端和5端切除固定长度的序列，使用-u/--cut 参数， 当值为正数时，从序列开头切除， 当值为负数时，从序列末端开始切除

从序列开头切除5bp，用法：

cutadapt -u 5 -o output.fastq  input.fastq

input.fastq 的内容：

 

@NB001
ATGCATGCTAGC
+
////////////

output.fastq的内容：

@NB001
TGCTAGC
+
///////

从序列末端切除5bp，用法：

cutadapt -u -5 -o output.fastq  input.fastq

 

 input.fastq 的内容：

@NB001
ATGCATGCTAGC
+
////////////

output.fastq的内容：

@NB001
ATGCATG
+
///////

9） cutadapt 的其他参数

--length_tag 在输出文件中给出 序列的length 信息

>read1 length=10

ACGTACGTAC

 

 

filter  reads

-m N/ --minimum-length N      去除序列长度低于N的read

--too-short-output  file           将序列长度低于N的read 输出到文件中

-M N/ --maximum-length N      去除长度大于N的read

--too-long-output  file             将序列长度大于N的read 输出到文件中

--discard-trimmed                   去除出现adapter 序列的read

--untrimed-output                  将没有进行出现过adapter 的序列单独输出到一个文件中

--discard-untrimmed               去除没有出现过adapter序列的read

 

cutadapt 的处理顺序：