常见的质量控制软件有：Trimmomatic seqtk SOAPnuke
1.Trimmomatic
之前这个软件也介绍过，只是当时是为了知道有这个软件而写，但是这回要把这个软件的细节搞明白。
以下是官方给出了一般操作：
Paired End：
java -jar trimmomatic-0.35.jar PE -phred33 input_forward.fq.gz input_reverse.fq.gz output_forward_paired.fq.gz output_forward_unpaired.fq.gz output_reverse_paired.fq.gz output_reverse_unpaired.fq.gz ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36
参数解读：PE 代表数据为双端测序所得。-phred33代表数据质量格式为phred33，默认是Phred64（小心这个坑呦）。input_forward.fq.gz input_reverse.fq.gz output_forward _ pai red.fq.gz 这是两个输入文件，就是你的双端测序的fq1和fq2。
output_forward_paired.fq.gz output_forward_unpaired.fq.gz output_reverse_paired.fq.gz output_reverse_unpaired.fq .gz 是四个输出文件，可以简单认为两个paired文件为过滤好的，两个unpaired文件为要丢弃的文件。ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 代表接头序列切除参数，TruSeq3-PE.fa是软件自带的一个文件，内容如下：
$ more TruSeq3-PE.fa >PrefixPE/1 TACACTCTTTCCCTACACGACGCTCTTCCGATCT >PrefixPE/2 GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCT
而且自带文件不止这一个，还有很多。适用于不同情况的测序接头，但是很奇怪的是我做的项目数据是PE数据，接头在他自带的接头文件中却是SE接头：
$ more TruSeq3-SE.fa >TruSeq3_IndexedAdapter AGATCGGAAGAGCACACGTCTGAACTCCAGTCAC >TruSeq3_UniversalAdapter AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTA
幸好软件支持自定义接头文件，只要你明白规则（接头名要搞好）就可以了：
‘Palindrome’ trimming is specifically designed for the case of ‘reading through’ a short fragment into the adapter sequence on the other end. In this approach, the appropriate adapter sequences are ‘in silico ligated’ onto the start of the reads, and the combined adapter+read sequences, forward and reverse are aligned. If they align in a manner which indicates ‘read-through’, the forward read is clipped and the reverse read dropped (since it contains no new data).
Naming of the sequences indicates how they should be used. For ‘Palindrome’ clipping, the sequence names should both start with ‘Prefix’, and end in ‘/1’ for the forward adapter and ‘/2’ for the reverse adapter. All other sequences are checked using ‘simple’ mode. Sequences with names ending in ‘/1’ or ‘/2’ will be checked only against the forward or reverse read. Sequences not ending in ‘/1’ or ‘/2’ will be checked against both the forward and reverse read. If you want to check for the reverse-complement of a specific sequence, you need to specifically include the reverse-complemented form of the sequence as well, with another name.
然后，2是比对时接头序列时所允许的最大错配数；30指的是要求PE的两条read同时和PE的adapter序列比对，匹配度加起来超30%，那么就认为这对PE的read含有adapter，并在对应的位置需要进行切除。10和前面的30不同，它指的是，我就什么也不管，反正只要这条read的某部分和adpater序列有超过10%的匹配率，那么就代表含有adapter了，需要进行去除。LEADING:3 代表read开头的碱基是否要被切除的质量阈值。SLIDINGWINDOW:4:15 代表滑动窗口长度的参数，SLIDINGWINDOW:5:20代表窗口长度为5，窗口中的平均质量值至少为20，否则会开始切除。MINLEN:36 代表read被切除后至少需要保留的长度，如果低于该长度，会被丢掉。
Single End：
java -jar trimmomatic-0.35.jar SE -phred33 input.fq.gz output.fq.gz ILLUMINACLIP:TruSeq3-SE:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36
单端测序数据只有一点不同就是不用要指定文件来存放被过滤掉的read信息的，后面直接就接Trimmer信息。
存在的问题：
不知道为什么，在使用这个软件进行质控之后，clean数据还是有大量的有接头文件。在使用Soapnuke之后就没有。
2.Soapnuke
基本操作：
SOAPnuke1.5.5 filter -f AGATCGGAAGAGCACACGTCTGAACTCCAGTCAC -r AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGTA -1 IEToleR1ABD16FAAPEI-202_1.fq.gz -2 EToleR1ABD16FAAPEI-202_2.fq.gz -C IEToleR1ABD16FAAPEI-202_1.clean.fq.gz -D IEToleR1ABD16FAAPEI-202_2.clean.fq.gz -o ../cleandata/ -l 5 -q 0.5 -n 0.03 -i -G -Q 2 -L 150 -5 1
参数解读：filter代表过滤普通序列。-f代表fq1的接头。-r代表fq2的接头。-1代表fq1。-2代表fq2。-C 代表clean fq1。-D 代表clean fq2。-o代表输出文件路径。-l代表判定低质量碱基的阈值，默认为5。-q 代表允许序列低质量比例，默认为50%，一旦高于就剔除。-n 代表N占序列的比例，默认0.05，一旦高于就会剔除序列。-i 代表去除index。-G 如果设置了就表示质量值体系选择为phred33 -Q 代表我们过滤后输出的read要选用哪种质量值体系输出，1代表phred64，2代表phred33，默认依然是phred33。-L 代表 fq1序列最长长度。-5 0代表数据名字为旧式的，如：(默认为0)
old fastq name: @FCD1PB1ACXX:4:1101:1799:2201#GAAGCACG/2
-5 1代表数据名字为新式的，如;
new fastq name: @HISEQ:310:C5MH9ANXX:1:1101:3517:2043 2:N:0:TCGGTCAC
这里有一份我使用这个软件进行过滤的结果：可以看出去除的序列基本上全是有接头的，低质量的很少，这里有个疑问就是，是软件过滤的问题，还是数据本来就不错有待后期再验证。
| Total | Percentage | Counts(fq1) | Percentage_2 | Counts(fq2) | Percentage_3 | |
|---|---|---|---|---|---|---|
| Total filtered reads (%) | 760,442 | 100.00 | 380,221 | 100.00 | 380,221 | 100.00 |
| Reads with adapter (%) | 746,650 | 98.19 | 296,724 | 78.04 | 349,062 | 91.81 |
| Reads with low quality (%) | 0 | 0.00 | 0 | 0.00 | 0 | 0.00 |
| Reads with low mean quality (%) | 0 | 0.00 | 0 | 0.00 | 0 | 0.00 |
| Reads with duplications (%) | 0 | 0.00 | 0 | 0.00 | 0 | 0.00 |
| Read with n rate exceed: (%) | 13,792 | 1.81 | 1,349 | 0.35 | 6,465 | 1.70 |
| Read with small insert size: (%) | 0 | 0.00 | 0 | 0.00 | 0 | 0.00 |
| Reads with PolyA (%) | 0 | 0.00 | 0 | 0.00 | 0 | 0.00 |
3. seqtk
seqtk本来是用于处理FASTA或FASTQ格式序列的快速轻量级工具，涉及到转换，切割等一系列的小操作。只因为有个命令为:trimq 切割序列低质量的碱基才进入我们讨论的范围。
命令十分简单，使用Phred算法从两端修剪低质量：（无法切割带有接头的序列）
seqtk trimfq in.fq > out.fq4.后话
在了解seqtk后，我又花了一天时间看了seqkit和csvkit的文档，他们所做的事就是我们日常的小操作，但是却能十分高效地工作。如果有兴趣的话，也可以研究他他们的源码学习一下，下一篇，我将介绍这款瑞士军刀。
Flystar 