一、下载参考基因组和参考基因注释文件构建
可以从UCSC、Ensembl、NCBI下载,但是要注意参考基因组与注释文件得配套。
1、Ensembl
参考基因组下载
先进入Ensembl主页,找到需要的物种。我需要的是人类的,主页就有,点击进入Ensembl主页
找到下载基因组fasta文件
基因组文件fasta
点进去就能看见按照染色体分类的基因组文件,下载所需要的染色体。
基因组fasta文件
#复制文件链接,我下载的是第一条染色体
wget https://ftp.ensembl.org/pub/release-108/fasta/homo_sapiens/dna/Homo_sapiens.GRCh38.dna.chromosome.1.fa.gz
注释文件下载
如图位置
gtf文件
下载gtf文件
2、UCSC
参考基因组
进入UCSC主页,点击downloads->Genome DataUCSC
找到自己需要的物质,点击
找到自己需要的版本的sequence
找到自己需要染色体的.fa.gz文件,复制下载地址,按照上面的方法下载
下载注释文件
点击Tools->Table Browser二、bowtie2构建index和Tophat2 mapping
对数据进行解压
gunzip *.fa.gz
将两条染色体合并到一起
cat Homo_sapiens.GRCh38.dna.chromosome.1.fa Homo_sapiens.GRCh38.dna.chromosome.2.fa > hg38.chr12.fa
1、bowtie2构建index
bowtie2官网
https://bowtie-bio.sourceforge.net/bowtie2/index.shtml
Bowtie2是一种用于将测序读数与长参考序列对齐的超快且高效记忆的工具。它特别擅长对齐大约50到100或1000个字符的读数,特别擅长对齐相对较长的(例如哺乳动物)基因组。Bowtie2使用FM索引(基于Burrows-Wheeler Transform 或 BWT)对基因组进行索引,以此来保持其占用较小内存。对于人类基因组来说,内存占用在3.2G左右。Bowtie2 支持间隔,局部和双端对齐模式。
conda install bowtie2
bowtie2-build hg38.chr12.fa hg38 > bwt_index.log 2>&1 &
两条染色体,使用一个CPU构建索引大概需要10分钟左右。
2、Tophat2 mapping
Tophat官网
https://ccb.jhu.edu/software/tophat/index.shtml
相对速度快并且占用内存小的TopHat在RNA-Seq中常常用于跨内含子比对,这里我们来关注一下TopHat2,TopHat2的安装是依赖Bowtie2的(当然Bowtie1也是可行的),TopHat2适合于75bp以上Read的比对。TopHat2是一个多步骤比对的程序,如果基因组的注释文件存在的话,那么它会先将Read比对到转录组上。
安装
wget https://ccb.jhu.edu/software/tophat/downloads/tophat-2.1.1.Linux_x86_64.tar.gz
#解压
tar -zxvf tophat-2.1.1.Linux_x86_64.tar.gz
#加入PATH
export PATH=/root/software/tophat-2.1.1.Linux_x86_64/:$PATH
#测试
tophat2 -h
出现错误,因为我系统默认的是python3的环境,tophat需要的是python2的环境。#创建python2的环境
conda create -n py2 python=2.7
#进入py2
conda activate py2
tophat2 -h
#可以成功调用
但是由于bowtie2没有安装在这个py2的环境下,就再次安装一次
conda install bowtie2
调用tophat2
tophat2 -o /root/rna_seq/tophat_out -p 1 -G /root/rna_seq/ref/gft/Homo_sapiens.GRCh38.108.chr.gtf /root/rna_seq/ref/hg38 /root/rna_seq/cutadapt_out/test_R1_cutadapt.fq.gz /root/rna_seq/cutadapt_out/test_R2_cutadapt.fq.gz &
#-o 输出文件位置
#-p 几个CPU
#gtf文件
#index文件
#reads1
#reads2
tophat输出结果
bam文件就是比对结果,也就是我们需要的文件