转录组分析所需软件

Trimmomatic,Hisat2,SAMtools,Stringtie/Htseq,Deseq

首先先下个conda,便于下载某些软件:

1.先将这个链接中的命令跟着走一遍(安装的时候如果失败/出错,重新试一次就可以了,它经常出错):用conda安装RNA-seq所需要的工具

2.安装软件

conda search xxx

conda install xxx

3.使用:刚安装完的时候直接输入软件名就可以直接打开,但重启终端之后需要先启动环境

source activate

conda activate python2

一、数据质控:Trimmomatic

1.下载安装:

wget http://www.usadellab.org/cms/uploads/supplementary/Trimmomatic/Trimmomatic-0.38.zip

unzip Trimmomatic-0.38.zip

2.输入文件:sample_R1.fastq、sample_R2.fastq

(java -jar <path to trimmomatic.jar> )PE [-threads <threads] [-phred33 | -phred64] [-trimlog<logFile>] >] [-basein <inputBase> | <input 1> <input 2>] [-baseout <outputBase> |<paired output 1> <unpaired output 1> <paired output 2> <unpaired output 2> <step 1> <step 2> ...

-threads 设置多线程运行数 ,-phred可不选

产生4个结果文件

参考:高通量测序数据质控神器TrimmomaticTrimmomatic下载安装Trimmomatic基本功能的使用


二、数据比对:HISAT2

需要用到的数据:基因组序列,转录组序列;结果得到SAM文件

基因组比对软件常用bwa,转录组比对软件常用bowtie2、hisat2等,其中有参考基因组的常用hisat2,没有参考基因组的常用bowtie2。

此处是要将转录组比对到基因组上去。

1.可以使用conda下载:

conda install hisat2

PATH=/home/freshman/lbq/hisat2:$PATH

PATH=/home/freshman/lbq/hisat2/*:$PATH

2.建立基因组、转录组索引:生成的文件会默认建立在你现在所在的目录下

(hisat2-build --ss hg19.splice_sites.gtf --exon hg19.exons.gtf genome/hg19/hg19.fa hg19)

hisat2-build –p 4 ~/data/cdna.fasta genome

hisat2-build –p 4 ~/data/dna.fasta gdna

3.对比命令及结果如下:

hisat2 -f -p16 -x ~/data/gdna -U cdna.fasta -S gly.sam


参考:序列比对:Hisat2转录组序列比对


三、将SAM转换为BAM文件:SAMtools

需要HISAT产生的SAM文件,转成BAM文件后对BAM文件进行分析

1.使用conda下载samtools

conda install samtools

2.利用将sam文件转为bam文件(samtools同样可以用conda下载)

samtools view -S gly.sam -b > gly.bam

3.查看bam文件

cd data

samtools flagstat gly.bam


参考:SAMTools


四、重复序列注释:StringTie

需要sort过的BAM文件与GTF文件

参考:Stringtie的使用用法详解

(四五二选一)


五:表达量分析——Htseq

输入文件:SAM和GTF,其中SAM需要samtools重新sort

1.下载安装

wget https://pypi.python.org/packages/46/f7/6105848893b1d280692eac4f4f3c08ed7f424cec636aeda66b50bbcf217e/HTSeq-0.7.2.tar.gz

tar zxf HTSeq-0.7.2.tar.gz

cd HTSeq-0.7.2

python setup.py install --user

vim .bashrc

PATH=/home/freshman/lbq/HTSeq-0.7.2/.local/bin:$PATH

source .bashrc

chmod 777 /home/freshman/lbq/HTSeq-0.7.2/python2/scripts/htseq-count

export PATH=$PATH:/home/freshman/lbq/HTSeq-0.7.2/*

2.计数

htseq-count -f bam -r name -s no -a 10 -t exon -i ID -m union ~/data/gly_sort.bam ~/data/ori/genes.gff3 > ~/data/counts_out.txt

# 命令参数

-f | --format default: sam 设置输入文件的格式,该参数的值可以是sam或bam。

-r | --order default: name 设置sam或bam文件的排序方式,该参数的值可以是name或pos。前者表示按read名进行排序,后者表示按比对的参考基因组位置进行排序。若测序数据是双末端测序,当输入sam/bam文件是按pos方式排序的时候,两端reads的比对结果在sam/bam文件中一般不是紧邻的两行,程序会将reads对的第一个比对结果放入内存,直到读取到另一端read的比对结果。因此,选择pos可能会导致程序使用较多的内存,它也适合于未排序的sam/bam文件。而pos排序则表示程序认为双末端测序的reads比对结果在紧邻的两行上,也适合于单端测序的比对结果。很多其它表达量分析软件要求输入的sam/bam文件是按pos排序的,但HTSeq推荐使用name排序,且一般比对软件的默认输出结果也是按name进行排序的。

-s | --stranded default: yes 设置是否是链特异性测序。该参数的值可以是yes,no或reverse。no表示非链特异性测序;若是单端测序,yes表示read比对到了基因的正义链上;若是双末端测序,yes表示read1比对到了基因正义链上,read2比对到基因负义链上;reverse表示双末端测序情况下与yes值相反的结果。根据说明文件的理解,一般情况下双末端链特异性测序,该参数的值应该选择reverse(本人暂时没有测试该参数)。

-a | --a default: 10 忽略比对质量低于此值的比对结果。在0.5.4版本以前该参数默认值是0。

-t | --type default: exon 程序会对该指定的feature(gtf/gff文件第三列)进行表达量计算,而gtf/gff文件中其它的feature都会被忽略。

-i | --idattr default: gene_id 设置feature ID是由gtf/gff文件第9列那个标签决定的;若gtf/gff文件多行具有相同的feature ID,则它们来自同一个feature,程序会计算这些features的表达量之和赋给相应的feature ID。

-m | --mode default: union 设置表达量计算模式。该参数的值可以有union, intersection-strict and intersection-nonempty。这三种模式的选择请见上面对这3种模式的示意图。从图中可知,对于原核生物,推荐使用intersection-strict模式;对于真核生物,推荐使用union模式。

-o | --samout 输出一个sam文件,该sam文件的比对结果中多了一个XF标签,表示该read比对到了某个feature上。

-q | --quiet 不输出程序运行的状态信息和警告信息。

-h | --help 输出帮助信息。

其中,gene_id这一处如果报错的话,less genes.gff3看看,它那一列的名字是不是叫其他的,改了就行。

结果如下:其中ambiguous表示该read比对到多个gene上;no_feature表示read没有比对到gene上

参考:RNA-seq分析htseq-count的使用使用HTSeq进行有参转录组的表达量计算 | 陈连福的生信博客


六、差异表达分析——Deseq(R包)

输入文件:read count

1.安装:(此处R版本为3.6.1)

R

install.packages("BiocManager")

BiocManager::install(c("DESeq2"))

2.


参考:DEseq2 - 简书DESeq2的基本使用 - 简书用DESeq2包来对RNA-seq数据进行差异分析

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