RNA-seq常用命令(无参)

0.前期准备

先在工作目录下创建以下几个目录:

01.raw_data        #用于存放原始数据
02.fastq           #用于存放fastq格式数据
03.fastqc          #用于存放QC结果
04.trinity_result  #用于存放trinity结果
数据处理的顺序:01.raw_data > 02.fastq > 03.fastqc > 04.trinity_result

1.下载SRA

cd01.raw_data目录下,执行以下命令(二选一):

nohup axel -q ftp://ftp.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByRun/sra/SRR/SRR130/SRR1300434/SRR1300434.sra &

axel是一种多线程下载工具,比wget快;-q参数能让axel进入静默模式(选用)。

nohup wget -c ftp://ftp.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByRun/sra/SRR/SRR130/SRR1300434/SRR1300434.sra &

用wget的时候记得加-c参数,支持断点续传,就会方便很多。

用脚本批量下载数据的代码参考附录里的内容。


2.fastq-dump

  fastq-dump --defline-seq '@$sn[_$rn]/$ri' --split-files file.sra --split-3 -O|--outdir

这个是trinity给的命令,如果后面要用trinity进行无参组装的就需要加上--defline-seq '@$sn[_$rn]/$ri'这串看不懂的“乱码”,就要用这个命令,否则trinity会报错说无法识别这个格式。

经过改良后的命令:

nohup fastq-dump --defline-seq '@$sn[_$rn]/$ri' --split-3 file.sra -o ../02.fastq &

默认输出在上一层目录的02.fastq这个文件夹里。

-o:重定向输出结果的位置。将结果输出到上一级目录下的02.fastq文件夹内。

可选用的参数:

--gzip, --bzip2: 以压缩文件的方式输出结果。有利于减少文件的占用空间。


3.质控

不管什么数据,质控一下还是很有必要的,多看fastqc结果,即是积累经验,也是及时发现数据可能出现的问题。首先cd../02.fastq,之后执行:

  fastqc -t 8 *.fq -o ../03.fastqc

-t: 使用的线程数。现在用的服务器是20 CPU 40线程,以后这个值可以改大一点(range from 1~40,但是一般取10,如果文件比较多可以取20。)

-o: 重定向输出结果的位置。放在03.fastqc里面。

filezilla登陆服务器找到对应文件夹将生成的html文件下载到本地进行查看。如无问题则可进行下一步Trinity的操作,如有问题则进行质控部分余下的操作。

因为篇幅原因,将单独写一个关于如何查看fastqc结果的教程。敬请期待

需要去除接头的情况:

  • cutadapt

http://mp.weixin.qq.com/s/jZ37itt48slYrMNSvRxjNg

  • trimmomatic

https://www.jianshu.com/p/7b5591673255

需要对双端测序数据进行过滤低质量reads操作的情况:

  • sickle

https://mp.weixin.qq.com/s/ALY2FJCO10x02A7cy9qHyA

当你想改善你的数据的问题但是你又新入门啥软件都不会的时候(推荐使用):

  • fastp

https://www.jianshu.com/p/fae5b3c11e74

  fastp -i name.fastq -o name.fastp.fastq  #单端
  fastp -i name_1.fastq -o name_1.fastp.fastq -I name_2.fastq -O name_2.fastp.fastq  #双端

4.Trinity

cd../fastq文件夹后,执行以下命令:

双端测序命令:

nohup Trinity --seqType fq\fa --left ?  --right ?  --CPU 20 --max_memory 50G -output ../trinity_result/?.trinity &
  • left和right不能随意设置,left必须是_1,right必须是 _2.

  • 一般type是fq

  • ?的位置用具体需要操作的文件名代入即可。

单端测序命令:

nohup Trinity --seqType fq\fa --single ?  --CPU 20 --max_memory 50G -output ../trinity_result/?.trinity &

ps:输出结果的名字里必须含有trinity这个关键词,否则无法正常运行。


5.TransDecoder

cd../trinity_result/filename.trinity文件夹下后,找到一个叫做Trinity.fasta的文件,执行:

#step 1: 提取最长的开放阅读框
  TransDecoder.LongOrfs -t ./Trinity.fasta
  
  #step 2: (可选),BlastP搜索和Pfam搜索
    #BlastP搜索:蛋白库搜索, Swissprot (快) or Uniref90 (慢 but more comprehensive)
  blastp -query transdecoder_dir/longest_orfs.pep -db uniprot_sprot.fasta -max_target_seqs 1 -outfmt 6 -evalue 1e-5 -num_threads 10 > blastp.outfmt6
    #Pfam搜索:肽或蛋白域预测,需要安装hmmer3和Pfam数据库
  hmmscan --cpu 8 --domtblout pfam.domtblout /path/to/Pfam-A.hmm transdecoder_dir/longest_orfs.pep
  ​
  #step 3: 将Blast和Pfam搜索结果整合到编码区域选择这个步骤在尝试中发现比较费时间,建议用nohup命令进行操作。
  nohup TransDecoder.Predict -t target_transcripts.fasta --retain_pfam_hits pfam.domtblout --retain_blastp_hits blastp.outfmt6 &
  ​
  nohup TransDecoder.Predict -t trinity.fasta &  #简化版

通过TransDecoder这一步将能获得cds、pep、gff3、bed文件。用filezilla将以下五个文件拿到本地:

  • Trinity.fasta

  • Trinity.fasta.transdecoder.cds

  • Trinity.fasta.transdecoder.pep

  • Trinity.fasta.transdecoder.gff

  • Trinity.fasta.transdecoder.bed

    将这几个文件以可分辨的名字进行命名(如改成SRR1300215.fasta等)。


6.出个统计报告

还是在本目录下,对Trnity.fasta文件进行操作:

  TrinityStats.pl trinity_out_dir/Trinity.fasta >> stat.txt

会生成一个统计报告,并将内容重定向到stat.txt这个文件中。这个文件也记得保存下来。

结果展示。

以SRR1300215为例:

  ################################
  ## Counts of transcripts, etc.
  ################################
  Total trinity 'genes':  1524
  Total trinity transcripts:  1793
  Percent GC: 54.43
  ​
  ########################################
  Stats based on ALL transcript contigs:
  ########################################
  ​
    Contig N10: 1008
    Contig N20: 757
    Contig N30: 597
    Contig N40: 487
    Contig N50: 397
  ​
    Median contig length: 291
    Average contig: 379.04
    Total assembled bases: 679625
  ​
  ​
  #####################################################
  ## Stats based on ONLY LONGEST ISOFORM per 'GENE':
  #####################################################
  ​
    Contig N10: 965
    Contig N20: 715
    Contig N30: 559
    Contig N40: 461
    Contig N50: 382
  ​
    Median contig length: 287
    Average contig: 369.46
    Total assembled bases: 563055</pre>

7.附录

批量下载sra文件的脚本模板(mac版本,仅供参考)

  # 在终端运行
  #用axel下载会快一些
  brew install axel
  for i in {56..62}
  do
  # 也可以用wget 下载
  axel ftp://ftp-trace.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByStudy/sra/SRP/SRP075/SRP075747/SRR35899$i/SRR35899$i.sra
  done

批量处理fastq-dump的脚本模板:(从人家的教程里摘下来的,仅供参考)

  #!/usr/bin/env bash
  for i in {56..62}
  do
  fastq-dump --gzip --split-3 -O /Users/chengkai/Desktop/zhuanlu_files -A SRR35899${i}.sra
  done

8.参考来源

文章:

微信公众号:

  • 生信菜鸟团

  • 生信者言


9.版权信息

Author:陈俊浩 Hanschen

QQ:1020607557

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 204,293评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,604评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 150,958评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,729评论 1 277
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,719评论 5 366
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,630评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,000评论 3 397
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,665评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,909评论 1 299
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,646评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,726评论 1 330
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,400评论 4 321
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,986评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,959评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,197评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 44,996评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,481评论 2 342

推荐阅读更多精彩内容

  • suse linux 常用命令 (1)命令ls——列出文件 ls -la 给出当前目录下所有文件的一个长列表,包括...
    小白_pk阅读 12,345评论 0 5
  • 第十章 使用序列数据 生物信息学的核心问题之一是处理大量的(通常定义糟糕或模糊)文件格式。久而久之,一些特定的简单...
    yangliunk1987阅读 4,984评论 3 53
  • 个人学习批处理的初衷来源于实际工作;在某个迭代版本有个BS(安卓手游模拟器)大需求,从而在测试过程中就重复涉及到...
    Luckykailiu阅读 4,679评论 0 11
  • Spring Cloud为开发人员提供了快速构建分布式系统中一些常见模式的工具(例如配置管理,服务发现,断路器,智...
    卡卡罗2017阅读 134,594评论 18 139
  • State&生命周期 1.当 <Clock /> 被传递给 ReactDOM.render() 时,React 调...
    iambabewin阅读 491评论 0 0