我的ChIP-Seq(2): cutadapt/fastp/Tr

过滤软件的比较与选择：cutadapt/fastp/trimmomatic

注：还没有完全搞明白，先总结一下特点和使用，之后再慢慢体会、总结经验
本次只针对双端PE
算法都没好好读，因为看不懂==

首先，我们对数据进行过滤，是为了：

去掉接头
去掉低质量reads
去掉污染序列
在尽量去掉上述序列的同时，保留尽可能多的有用数据，减少损失

CutAdapt，2010

1. 基于Python，作者是个德国人，长得还挺帅气(✧◡✧) 不过都9年过去了，嗯。。
2. 不仅支持illumina，还支持SOLID，454等平台产出的数据
3. 支持输入.gz
4. 需要自己先检测接头类型（fastqc等），然后搜索接头序列是啥，手动输入到参数里。但是有一个参数 -n，若是两种接头，也可以指定然后去除：-n 2
5. 一般命令：

cutadapt -a -A #a是read1的3'接头，A是read2的3'接头(5'接头的反向互补序列)
-e 0.1 -0.5 -m 50 #去除接头后read长度大于50才保留
-o -p #生成文件：过滤后的R1 R2
read1.fastq read2.fastq #输入文件

4. 本次分析没用，所以详细参数可以阅读--help

fastp，2018

1. 基于c++这种强大的语言所以算法比较高效，中科院深圳所发的。还没用过，不过身边做RNA-Seq的俩师兄强烈推荐，有空可以test一下。
2. 主题就是ultra-fast，all-in-one，而且是只处理FASTQ也就是illuminate下机数据
3. 特点：

能进行质控，生成比fastqc美观、全面的报告，但是我看了一遍，不如fastqc直观、fresh-friendly
号称去除低质量序列的方法类似于trimmomatic但是更快
自动识别序列并去除
支持illuminate short read，也一定程度支持Pacbio/Nanopore long reads，具体支持到什么程度，需要试验。
参数众多，但是挺有条理的，而且很多都是默认不是必需参数，不会“新手退散”

4. 最简单的命令：
   fastp -i r1.fq -o rr1.fq -i r2.fq -o rr2.fq
5. 这篇介绍写的不错：知乎
   但他说一般下机数据要经过fastqc+cutadapt+trimmomatic，有点不太理解，要这么麻烦吗？

Trimmomatic，2014

1. 也是很好用的，引用量超高，good at去除低质量reads，只针对illuminate数据

2. 最重要的特点：对数据的处理步骤与参数的顺序有关！
   所以建议先去接头，否则接头被剪更无法有效去除。

3. PE数据常用参数：
   ILLMINACLIP: 注意以下参数以：隔开
   <fastaWithAdaptersEtc>: 指定包含接头和引物序列（所有被视为污染的序列）的 fasta 文件
   <seed mismatches>: 第一步seed搜索时允许的mismatch个数，一般2。
   <palindrome clip threshold>: 指定针对 PE的palindrome clip模式下，需要R1和 R2之间至少多少比对分值，才会进行接头切除，例如30。
   <simple clip threshold>: 指定切除接头序列的最低比对分值，一般7-15之间。
   <minAdapterLength>: 只对 PE 测序的 palindrome clip 模式有效，指定 palindrome 模式下可以切除的接头序列最短长度，默认值是 8。但实际上 palindrome 模式可以切除短至 1bp 的接头污染，所以可以设置为 1。
   <keepBothReads> 重要参数！第一次做的时候没加这个参数，结果20%+的数据Unpaired，扔掉不现实，比对处理太麻烦！正确用法：只对 PE 测序的 palindrome clip 模式有效，R1 和 R2 在去除了接头序列之后剩余的部分是完全反向互补的，默认参数 false，意味着整条去除与 R1 完全反向互补的 R2，当做重复去除掉，但在有些情况下，例如需要用到 paired reads 的 bowtie2 流程，就要将这个参数改为 true，否则会损失一部分 paired reads。

4. 本次所用命令：（也是公司报告中所用的）

java -jar trimmomatic-0.38.jar PE -threads 2 #双端模式，两个线程
ILLUMINACLIP: #顾名思义，去illumina接头
TruSeq3-PE.fa: #接头文件，需要指定全路径
2:30:10 # 默认格式为 2:30:10:8:false，可改做：2:30:10:8:true
LEADING:20 #从reads的起始端开始切除质量值低于设定的阈值的碱基，直到有一个碱基其质量值达到阈值。一般用LEADING:3???
TRAILING:20 #一般用3，因为Illumina 平台有些低质量的碱基质量值被标记为2，所以设置为 3 可以过滤掉这部分低质量碱基
SLIDINGWINDOW:4:20 #滑窗剪切，统计滑窗口中所有碱基的平均质量值，如果低于设定的阈值，则切掉窗口。此处设置4bp窗口，阈值20，一般阈值用15。
MINLEN:50 #可被保留的最短 read 长度

trimmomatic PE模式默认处理2个文件，也就是说，sh脚本中使用本办法只能一次列举R1 R2两个文件，不能 In_R1 In_R2 IP_R1 IP_R2这样四个文件都列出来，事实证明会报错，trimmomatic有点傻傻的不知道第三个开始的文件该干嘛。
所以要批量做，需要写循环，或者是认真阅读使用说明的参数。

5. trimmomatic的更多解读可以参考这个，写得很详细。目前我理解的是以上。

最后附一个图：
出自：Chen et al. Source Code for Biology and Medicine 2014, 9:8. Software for pre-processing Illumina nextgeneration sequencing short read sequences.

几种软件比较

以上。可以test一下trimmomatic的true参数，还有fastp试一下到底强大在哪里。