基因组相关文件提取和索引创建

提取文件

提取fasta文件

samtools faidx genome.fa chr1 > genome_chr1.fa
samtools faidx genome.fa chr1:1-1000 > genome_chr1.fa

提取bam文件

samtools view -hb sample.bam chr1:1-1000 > sample_chr1.bam

查看bam文件

samtools tview -p chr19:1000 sample.bam genome.fa (若不加-p,进入查看后按g可以到指定位置)

提取vcf文件

vcftools --gzvcf sample.vcf.gz --chr1 --from-bp 1 --to-bp 1000 --recode --recode-INFO-all --out prefix

将基因组根据窗口大小划分为bed文件

bedtools makewindows -g genome.fa.fai -w 5000 > test.bed

创建索引

Bam文件创建索引

samtools index sample.bam

Vcf文件创建索引

  1. 压缩vcf: bgzip sample.vcf (不能使用gzip压缩)
  2. 解压vcf: bgzip -d sample.vcf.gzgunzip sample.vcf.gz
  3. 创建索引: tabix -p vcf sample.vcf.gz (只能对压缩的vcf创建,索引后缀为.tbi)

GATK创建参考基因组索引

  1. 生成.fai: samtools faidx genome.fa
  2. 生成.dict: gatk CreateSequenceDictionary --REFERENCE genome.fa
  3. 生成.img: gatk BwaMemIndexImageCreator -I genome.fa -O genome.fa.img

Bwa创建参考基因组索引

  1. 基因组>2*10^9bp: bwa index -a bwtsw genome.fa
  2. 基因组<2*10^9bp: bwa index -a genome.fa
  • 生成文件: .amb .ann .bwt .pac .sa
  • 计算参考基因组碱基数: grep -v "^>" genome.fa | awk 'BEGIN{sum=0}{sum+=length($0)}END{print sum}'

Hisat2创建参考基因组索引

hisat2-build -p 4 genome.fa prefix (-p线程数)

Author: Giftbear
Link: https://giftbear.github.io/2023/04/06/基因组相关文件提取和索引创建/
Copyright Notice: All articles in this blog are licensed under CC BY-NC-SA 4.0 unless stating additionally.