Fastq-dump使用

来源:hoptop评论7,899

做生信的基本上都跟NCBI-SRA打过交道,尤其是fastq-dump大家肯定不陌生.NCBI的fastq-dump软件一直被大家归为目前网上文档做的最差的软件之一",而我用默认参数到现在基本也没有出现过什么问题,感觉好像也没有啥问题, 直到今天看到如下内容, 并且用谷歌搜索的时候,才觉得大家对fastq-dump的评价非常很到位.

Fastq-dump使用-图片1

过滤

我们一般使用fastq-dump的方式为

fastq-dump /path/to/xxx.sra 

但是这个默认使用方法得到结果往往很糟, 比如说他默认会把双端测序结果保存到一个文件里, 但是如果你加上--split-3之后, 他会把原来双端拆分成两个文件,但是原来单端并不会保存成两个文件. 还有你用--gzip就能输出gz格式, 能够节省空间的同时也不会给后续比对软件造成压力, 比对软件都支持,就是时间要多一点。

Fastq-dump使用-图片2

--gzip时间

Fastq-dump使用-图片3

没有gzip时间

但是很不幸运,这些东西在官方文档并没有特别说明,你只有通过不断的踩坑才能学到这些小知识。

我建议尽量去EMBL-EBI去下载原始数据,而不是这种神奇的sra格式,尽管有一些下载的数据其实就是从SRA解压而来。

不过要用fastq-dump,那就介绍几个比较重要的参数吧。我会按照不懂也加,不懂别加,有点意思,没啥意义这三个级别来阐述不同参数的重要级.

太长不看版

如果你不想了解一些细节性内容,用下面的参数就行了

fastq-dump --gzip --split-3 --defline-qual '+' --defline-seq '@\$ac-\$si/\$ri'   SRR_ID
# 建议加别名
alias fd='fastq-dump --split-3 --defline-qual '+' --defline-seq '@\\\$ac-\\\$si/\\\$ri' '

数据格式

不懂也加: reads拆分

默认情况下fastq-dump不对reads进行拆分, 对于很早之前的单端测序没有出现问题.但是对于双端测序而言,就会把原本的两条reads合并成一个,后续分析必然会出错.

常见的参数有三类:

  • --split-spot: 将双端测序分为两份,但是都放在同一个文件中
  • --split-files: 将双端测序分为两份,放在不同的文件,但是对于一方有而一方没有的reads直接丢弃
  • --split-3 : 将双端测序分为两份,放在不同的文件,但是对于一方有而一方没有的reads会单独放在一个文件夹里

关于遇到的Rejected 35403447 READS because of filtering out non-biological READS就是因为原来是SE数据,但是用--split-3当作PE数据处理,出现的问题. 看起来好像有问题,但是对后续结果分析没有太多影响.

因此,对于一个你不知道到底是单端还是双端的SRA文件,一律用--split-3.

没事别加: read ID

默认双端测序数据拆分后得到两个文件中同一个reads的名字是一样的,但是加上-I | --readids之后同一个reads的ID就会加上.1和.2进行区分.举个例子

是否有-I参数ID 1ID 2
@SRR5829230.1 1 length=36@SRR5829230.1 1 length=36
@SRR5829230.1.1 1 length=36@SRR5829230.1.2 1 length=36

问题来了, 明明已经可以通过ID后面的"1"和"2"来区分ID, 加这个参数干嘛. 加完之后还会让后续的BWA报错.所以,没事千万别加

有点意义: 原始格式

默认情况下输出的文件的ID都是SRR开头,但其实原始数据名字不是这样子,比如说@ST-E00600:143:H3LJWALXX:1:1101:5746:1016 2:N:0:CCTCCTGA,@HWI-ST620:248:HB11HADXX:2:1101:1241:2082#0/1这种. 如果你想看到那种格式,而不是SRR,你需要怎么做呢?

可以通过如下三个选项进行修改

  • F|--origfmt: 仅保留数据名字
  • --defline-seq <fmt>: 定义readsID的显示方式
  • --defline-qual <fmt>: 定义质量的显示方式

其中fmt按照如下要求定义

Fastq-dump使用-图片4

fmt的写法

虽然看起来有点意思,但是对最后的分析其实没啥帮助.

没啥意义: fasta输出

如果下游分析只需要用到fasta文件,那么用--fasta就行. 当然了也有很多方法能够把fastq转换成fasta,比如说samtools.

过滤

我觉得这部分的参数都没有意义, 毕竟完全可以用专门的质控软件处理reads,不过--skip-technical,是唯一比较重要.

  • 根据ID: -N -X
  • 根据长度: -M
  • 多标签序列: --skip-technical, 这个是唯一有点意思的,就是说如果你原来建库测序使用了多个标签来区分序列, 默认不会输出这个标签. 但是如果不输出标签,我们怎么区分呢? 所以一定要显示声明

有点意思: 输出方式

这部分参数也很重要, 选择是否压缩,还是直接输出到标准输出

  • --gzip, --bzip2: 压缩方式
  • -Z | --stdout : 输出到标准输出
  • -O|--outdir : 输出到指定文件夹

发表评论

匿名网友