一般来讲,我们对测序数据进行QC,就三个大的方向:Quality trimming, Adapter removal, Contaminant filtering,当我们是双端测序数据的时候,去除低质量的reads就容易导致左右两端测序文件不平衡,有一个比较好的软件能解决这个问题,而且软件使用非常简单!
安装代码如下:
## https://github.com/najoshi/sickle
cd ~/biosoft
mkdir sickle && cd sickle
wget https://codeload.github.com/najoshi/sickle/zip/master -O sickle.zip
unzip sickle.zip
cd sickle-master
make
~/biosoft/sickle/sickle-master/sickle -h
这个软件很简单,就是去除低质量的reads,而且还可以保证双端测序的完整性。
本实例的测试数据可以在 http://www.biotrainee.com/jmzeng/reads/test1.fastq 和 http://www.biotrainee.com/jmzeng/reads/test2.fastq
这个软件支持gz压缩格式,我应该压缩好了再上传到我们的云服务器的,这样可以节省流量,这只是一个测试,如果数据传输压力太大了,我们可能会取消链接,改为百度云分享!
~/biosoft/sickle/sickle-master/sickle pe -f test1.fastq -r test2.fastq -t sanger -o trimmed_output_file1.fastq -p trimmed_output_file2.fastq
软件给出的log日志如下:
PE forward file: test1.fastq
PE reverse file: test2.fastqTotal input FastQ records: 200000 (100000 pairs)
FastQ paired records kept: 192262 (96131 pairs)
FastQ single records kept: 3869 (from PE1: 3864, from PE2: 5)
FastQ paired records discarded: 0 (0 pairs)
FastQ single records discarded: 3869 (from PE1: 5, from PE2: 3864)
然后批量查看处理前后的fastqc质量报告:
ls *fastq |xargs -P 5 ~/biosoft/fastqc/FastQC/fastqc
比较所有的fastq文件的结果报告就可以看出它做了什么!
sickle处理前后的结果文件都可以在 http://www.biotrainee.com/jmzeng/sickle/sickle-results.zip 下载。
当然,这只是默认参数的用法,还可以添加很多参数! 比如 -q 30 -l 15 参数解释如下:
-
pe : use paired-end mode
-
-f training/rnaseq/ERR022486_chr22_read1.fastq.gz : the fastq file for read 1
-
-r training/rnaseq/ERR022486_chr22_read2.fastq.gz : the fastq file for read 2
-
-t sanger : the quality encoding. All data downloaded from EBI or NCBI will be "sanger" encoded. For an explanation:http://en.wikipedia.org/wiki/FASTQ_format#Encoding
-
-o ERR022486_chr22_read1_trim.fastq : the output file for trimmed reads from read 1
-
-p ERR022486_chr22_read2_trim.fastq : the output file for trimmed reads from read 2
-
-s ERR022486_chr22_single_trim.fastq : the output file for reads where the mate has failed the quality or length filter
-
-q 30 : the quality value to use. Bases below this will be trimmed, using a sliding window
-
-l 15 : the minimum length allowed after trimming. Here we remove reads with less than 15bp