输入文件有:物种基因组序列文件、物种基因组注释文件、转录本比对后的BAM文件;
以斑胸草雀为例:
-基因组注释文件: GCF_003957565.2_bTaeGut1.4.pri_genomic.gtf
-基因组序列文件: GCF_003957565.2_bTaeGut1.4.pri_genomic.fna
-转录本比对文件: myBAMfile.bam
文件内容一览:
samtools view myBAMfile.bam | less -S
查看BAM文件内容
第一步提取BAM文件
如果要查看所有转录本的比对情况基本是不现实的(一般上百G),受限于电脑性能无法加载这么大的bam文件,所以查看reads的比对情况一般通过提取出map到基因组的其中一个Scaffold的所有reads的bam文件。
示例从BAM文件提取斑胸草雀转录本比对到Scaffold NC_007897.1 的所有reads
samtools sort -@ 24 -o myBAMfile.sort.bam myBAMfile.bam && #bam文件排序
samtools index -@ 24 myBAMfile.sort.bam && #对排序后的bam建立索引
samtools view myBAMfile.sort.bam NC_007897.1 -O BAM -o NC_007897.1.bam && #提取比对到 NC_007897.1 Scaffold 的reads输出到NC_007897.1.bam
echo "FINISHED"
第二步用IGV浏览器浏览提取的reads的比对情况
打开IGV浏览器,首先载入基因组文件,点击Genomes选择加载本地基因组,导入GCF_003957565.2_bTaeGut1.4.pri_genomic.fna文件
然后选择File点击Load from file选择加载本地gtf文件,导入GCF_003957565.2_bTaeGut1.4.pri_genomic.gtf
创建基因组索引
加载BAM文件
工作目录下的 *.bam.bai 文件就是创建好的bam索引,重新加载BAM文件
定位到基因组Scaffold位置,即可浏览比对到该Scaffold的所有reads