有参转录组比对到参考基因组

一、下载参考基因组和参考基因注释文件构建

可以从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 Data
UCSC

找到自己需要的物质,点击

找到自己需要的版本的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文件就是比对结果,也就是我们需要的文件

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
【社区内容提示】社区部分内容疑似由AI辅助生成,浏览时请结合常识与多方信息审慎甄别。
平台声明:文章内容(如有图片或视频亦包括在内)由作者上传并发布,文章内容仅代表作者本人观点,简书系信息发布平台,仅提供信息存储服务。

相关阅读更多精彩内容

友情链接更多精彩内容