vcf for ancestry

A step-by-step guide on how to process a sorted BAM file to extract high-quality, biallelic SNPs using bcftools and tabix. Includes commands for variant calling, filtering by quality/depth, and isolating SNPs for downstream genomic analysis.
anil0w0 avatar
anil0w0
markdown
02/12/2026 11:22 AM
3.3 KB
16
Indexable

Perfect πŸ‘ β€” I’ll show you exactly how to extract and store:

βœ… SNP-only, biallelic VCF βœ… PLINK bed/bim/fam set

starting only from:

  • sorted.bam
  • reference: /home/user/Desktop/Syncthing/Syncthing_pc/genomics/references/hs38d1s_nochr.fa

🧬 STEP 1 β€” Create Variant-Only SNP+INDEL VCF

(This is the base file we derive everything from)

bcftools mpileup \ -f /home/user/Desktop/Syncthing/Syncthing_pc/genomics/references/hs38d1s_nochr.fa \ sorted.bam | \ bcftools call -m -v -Oz -o clinical_raw.vcf.gz

Index:

tabix -p vcf clinical_raw.vcf.gz

🧬 STEP 2 β€” Filter to High-Quality Variants

bcftools filter -Oz -o clinical.vcf.gz \ -e 'QUAL<20 || DP<10' clinical_raw.vcf.gz tabix -p vcf clinical.vcf.gz

Now you have:

clinical.vcf.gz (SNP + INDEL, high quality)


🧬 STEP 3 β€” Extract SNP-Only

bcftools view -v snps clinical.vcf.gz -Oz -o snps_only.vcf.gz tabix -p vcf snps_only.vcf.gz

🧬 STEP 4 β€” Keep Only Biallelic SNPs

(critical for PLINK & ancestry)

bcftools view -m2 -M2 snps_only.vcf.gz -Oz -o ancestry.vcf.gz tabix -p vcf ancestry.vcf.gz

βœ… You now have your MASTER ancestry VCF:

ancestry.vcf.gz
ancestry.vcf.gz.tbi

Size: ~200–600 MB This is your long-term archive.


🧬 STEP 5 β€” Convert to PLINK Format

plink --vcf ancestry.vcf.gz --make-bed --out mygenome

This creates:

mygenome.bed
mygenome.bim
mygenome.fam

Size: ~150–400 MB


πŸ“¦ What to Store Long-Term

Keep:

ancestry.vcf.gz
ancestry.vcf.gz.tbi
mygenome.bed
mygenome.bim
mygenome.fam

That’s it.

You can delete:

  • clinical_raw.vcf.gz
  • intermediate files
  • giant all-sites VCF
  • BAM (if ancestry-only forever)

🧠 Why This Is Future-Proof

From ancestry.vcf.gz you can:

Extract 23andMe SNPs

bcftools view -R 23andme_positions.bed ancestry.vcf.gz

Extract AADR 1240K

plink --bfile mygenome --extract 1240k_snps.txt --make-bed --out aadr

Extract any company's panel

Just subset SNP list.


🟒 Final Storage Size Comparison

FileSize
BAM30–90 GB
CRAM10–25 GB
SNP-only VCF200–600 MB
PLINK set150–400 MB

πŸ‘‰ ~100Γ— smaller than BAM.


If you want, I can now:

  • Create a single one-shot script
  • Validate strand orientation vs 1240K
  • Show how to make mock 23andMe file

Tell me your next target πŸ‘

Leave a Comment