搜索
查看: 4282|回复: 3

[mRNA-seq] 转录组入门分析3_fastq进行质控

[复制链接]

9

主题

15

帖子

132

积分

注册会员

Rank: 2

积分
132
发表于 2017-8-18 13:54:52 | 显示全部楼层 |阅读模式
本帖最后由 laofuzi 于 2017-8-18 13:54 编辑

1. fastq-dump程序函数简介
采用sratoolkit软件进行序列格式的转换。Sratoolkit软件先前已经安装,这里直接调用子程序fastq-dump进行分析。fastq-dump主要有以下几个主要参数:
-A|--accession  Replaces accession derived from in filename(s) and deflines (only for single able dump)。主要作用是改名。
-O|--outdir 指定输出路径
--gzip/--bzip2 输出为压缩格式gzip/bzip。fastqc软件可以直接识别gzip压缩的文件。输出文件压缩成bzip2格式(bzip2比传统的gzip或者ZIP的压缩效率更高,但是它的压缩速度较慢(来自Panda姐)。
--split-3多文件选项,Legacy 3-file splitting for mate-pairs: First biological reads satisfying dumping conditions are placed in files *_1.fastq and *_2.fastq If only one biological read is present it is placed in *.fastq Biological reads and above are ignored. 也就是说如果SRA文件中只有一个文件,那么这个参数就会被忽略。如果原文件中有两个文件,那么它就会把成对的文件按*_1.fastq, *_2.fastq这样分开(来自Panda姐)。

2. 安装Emacs
由于需要编辑一个简单的Shell脚本,所以需要编辑器。vim、Emacs或者其他均可。本人倾向于Emacs。在Ubuntu16.04系统下安装Emacs非常的简单,只需要在系统终端中执行以下三条命令即可:
[AppleScript] 纯文本查看 复制代码
sudo add-apt-repository ppa:ubuntu-elisp/ppa
sudo apt-get update
sudo apt-get install emacs-snapshot emacs-snapshot-el


Tips:可能有人会问,我在windows编辑好了脚本,放到Linux中直接运行可否?答:绝对不行。两个系统的换行符是不一样的,请自行百度。

启动Windows Xming,启动Ubuntu终端,输入如下命令,启动Emacs图形界面:
[AppleScript] 纯文本查看 复制代码
DISPLAY=:0 emacs


3. 将sra文件转换为fastq格式
首先,创建一个文件夹存放fastq文件
[AppleScript] 纯文本查看 复制代码
mkdir /mnt/d/AKAP95/fastq/


其次,在emacs中编辑如下脚本:
[AppleScript] 纯文本查看 复制代码
#!/bin/bash
for i in `seq 56 62`
do fastq-dump --gzip --split-3 /mnt/d/AKAP95/data/SRR35899${i}.sra -O /mnt/d/AKAP95/fastq/
done


将文件存为fastq_dump.sh

然后,在终端输入,如下命令:
[AppleScript] 纯文本查看 复制代码
bash /mnt/d/AKAP95/fastq/fastq_dump.sh


QUESTION:这个Shell脚本是参照论坛里其他大神的做法,问题是每次循环一步后,程序就会停止,然后需要按一下回车,程序才会继续进行下一步循环。不知道是不是这个Ubuntu子系统的问题,还是其他原因。单独跑应该不会出现这个问题,回头看看问题出哪。

4. md5sum数据传输完整性验证
[AppleScript] 纯文本查看 复制代码
md5sum /mnt/d/AKAP95/fastq/*.fastq.gz > /mnt/d/AKAP95/fastq/md5sum.txt


5. 了解fastq文件格式
输入如下命令,获取一条fastq文件:
[AppleScript] 纯文本查看 复制代码
zcat /mnt/d/AKAP95/fastq/SRR3589956_1.fastq.gz | head -n 4


Tips_1:压缩的文件不要着急解压,有很多bash命令能够直接用于压缩文件,如zgrep,zcat,zless,zdiff等。
Tips_2:这里采用了管道符将两个命令连接在一起。

Fastq文件格式如下所示:
@SRR3589956.1 D5VG2KN1:224:C4VAYACXX:5:1101:1159:2173 length=51
GGCGAGTGTAGGGCTGGCGCTGCCGGACGCGGTGCTAGTCGCCGGATGAAG
+SRR3589956.1 D5VG2KN1:224:C4VAYACXX:5:1101:1159:2173 length=51
B<BFBFBF0BFFFBFFBBFFIF<FFI<7<<BF<FFFFFFBB<BBBBBBBBB

FASTQ文件每个序列通常为4行,分别为:
@SRR3589956.1 D5VG2KN1:224:C4VAYACXX:5:1101:1159:2173 length=51  #第一行:@字符开头的标题行,分别为:设备名称/run id/flowcell id/flowcell lane/tile number within the flowcell lane/'x'-coordinate of the cluster within the tile/'y'-coordinate of the cluster within the tile/the member of a pair, 1 or 2/Y if the read is filtered, N otherwise/0 when none of the control bits are on, otherwise it is an even number/index sequence
GGCGAGTGTAGGGCTGGCGCTGCCGGACGCGGTGCTAGTCGCCGGATGAAG  #序列
+SRR3589956.1 D5VG2KN1:224:C4VAYACXX:5:1101:1159:2173 length=51   # +
B<BFBFBF0BFFFBFFBBFFIF<FFI<7<<BF<FFFFFFBB<BBBBBBBBB  #碱基质量格式phred+33

6. 质控
采用fastqc程序进行质控。
创建文件夹:
[AppleScript] 纯文本查看 复制代码
mkdir /mnt/d/AKAP95/fastq/fastqc/


启动Xming和Emacs,编辑如下脚本:
[AppleScript] 纯文本查看 复制代码
#!/bin/bash
cd /mnt/d/AKAP95/fastq/
fastqc -o ./fastqc/ -t 4  SRR3589956_1.fastq.gz SRR3589956_2.fastq.gz
fastqc -o ./fastqc/ -t 4  SRR3589957_1.fastq.gz SRR3589957_2.fastq.gz
fastqc -o ./fastqc/ -t 4  SRR3589958_1.fastq.gz SRR3589958_2.fastq.gz
fastqc -o ./fastqc/ -t 4  SRR3589959_1.fastq.gz SRR3589959_2.fastq.gz
fastqc -o ./fastqc/ -t 4  SRR3589960_1.fastq.gz SRR3589960_2.fastq.gz
fastqc -o ./fastqc/ -t 4  SRR3589961_1.fastq.gz SRR3589961_2.fastq.gz
fastqc -o ./fastqc/ -t 4  SRR3589962_1.fastq.gz SRR3589962_2.fastq.gz


将文件另存为fastqc.sh
其中,常用参数有:
-o: 输出路径-
-extract:输出文件是否需要自动解压 默认是--noextract-
-t:线程数,和电脑配置有关,每个线程需要250MB的内存
-c:测序中可能会有污染, 比如说混入其他物种
-a:接头
-q:安静模式

Tips:fastqc可以有多个输入文件

输入如下命令,进行质控分析:

[AppleScript] 纯文本查看 复制代码
bash /mnt/d/AKAP95/fastq/fastqc.sh


Tips:可以启动Xming,输入如下命令,启动fastqc的图形界面,同样可以进行分析:
[AppleScript] 纯文本查看 复制代码
DISPLAY=:0 fastqc


对于每一个文件,将会产生一个html的网页文件和一个zip的压缩文件,其中html文件用浏览器打开就能直观看到数据。对于fastqc文件的解读,论坛里和网上很多,认真读一遍就行。

7. Multiqc查看质控结果
对于多个文件,一个个查看比较费劲,考虑采用Multiqc程序(官网http://multiqc.info/)进行批量查看QC结果。输入如下命令,安装multiqc:
[AppleScript] 纯文本查看 复制代码
conda install -c bioconda multiqc


输入如下命令,获取multiqc的网页结果:
[AppleScript] 纯文本查看 复制代码
cd /mnt/d/AKAP95/fastq/fastqc/
multiqc .


后记:论坛里有大神贴出了Python脚本,可以考虑一试,好好研究。




上一篇:转录组入门分析2_从文章中获取数据
下一篇:转录组入门分析4_参考基因组和注释文件的下载及IGV的...
回复

使用道具 举报

634

主题

1182

帖子

4030

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
4030
发表于 2017-8-19 09:05:06 | 显示全部楼层
看样子, 你还不会用批处理,shell学习要加强
你这个问题很复杂,需要打赏,请点击 http://www.bio-info-trainee.com/donate 进行打赏,谢谢
回复 支持 反对

使用道具 举报

9

主题

15

帖子

132

积分

注册会员

Rank: 2

积分
132
 楼主| 发表于 2017-8-19 10:13:30 | 显示全部楼层
Jimmy 发表于 2017-8-19 09:05
看样子, 你还不会用批处理,shell学习要加强

嗯~ o(* ̄▽ ̄*)o,有您的鞭策,我要更努力
回复 支持 反对

使用道具 举报

0

主题

29

帖子

259

积分

中级会员

Rank: 3Rank: 3

积分
259
发表于 2017-12-19 14:30:19 | 显示全部楼层
。。详细详细。。不错错不错
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|手机版|小黑屋|生信技能树 ( 粤ICP备15016384号  

GMT+8, 2019-10-22 15:15 , Processed in 0.034457 second(s), 25 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.