DOCSEARCHbyAlgolia
Submit Docs

BCFtools Documentation

repository·develop·Indexed 21 days ago

https://github.com/samtools/bcftools

A suite of command-line tools for manipulating, analyzing, and calling variants from high-throughput sequencing data, specifically focusing on VCF and BCF file formats. BCFtools provides utilities for SNP/indel calling, copy number variation (CNV) calling, VCF annotation, consensus sequence generation, and format conversion between VCF/BCF and other genomic formats.

Tokens
13.6K
Snippets
40
Records
60
Agent score
69%

What's inside BCFtools

  1. Overview of BCFtools

    develop

    BCFtools is a suite of command-line tools for processing and analyzing variant call format (VCF) and binary variant call format (BCF) files. It consolidates various vcf* commands (such as vcfcheck, vcfmerge, and vcfisec) and includes BCF calling capabilities.

    For comprehensive documentation, users should refer to the official bcftools GitHub page.

  2. Overview of BCFtools capabilities

    develop
    BCFtools is a suite of utilities designed to manipulate variant calls in the Variant Call Format (VCF) and its binary counterpart, BCF. It supports uncompressed and BGZF-compressed files and works transparently with both formats. Many commands are designed to work with Unix pipes, treating - as standard input (stdin) and standard output (stdout).

  3. Understand mpileup isec output files

    develop

    After running the mpileup scripts, the output directories contain .isec directories representing the intersection of calls with the truth set. The files are decomposed into separate alleles (e.g., GT 1/2 becomes two entries) to facilitate easier counting of genotypes.

    File Mapping:

    • 0000.vcf: False Negatives (variants only found in the truth set).
    • 0001.vcf: False Positives (variants only found in your calls).
    • 000[23].vcf: True Variants (variants found in both the truth set and your calls).

  4. Use bcftools plugins

    develop

    Plugins are utilities that can be run using the same framework as standard commands. They are prefixed with a +.

    Usage Patterns:

    • List all available plugins: bcftools plugin -l
    • Run a plugin: bcftools plugin <NAME> <FILE> or bcftools +<NAME> <FILE>
    • Run from an explicit path: bcftools +/path/to/plugin.so <FILE>
    • Stream input: cat in.vcf | bcftools +<NAME>
    • View plugin help: bcftools +<NAME> -h

    Plugin Discovery: By default, system directories are searched. You can override this by setting the BCFTOOLS_PLUGINS environment variable to a colon-separated list of directories. If the variable starts or ends with a colon, system directories are also included at that position.

    bcftools +counts in.vcf

  5. Use bcftools plugins

    develop

    If ENABLE_BCF_PLUGINS is enabled, you can use user-defined plugins. Plugins can be invoked in two ways:

    1. Using the plugin command: bcftools plugin <plugin_name> <args>
    2. Using the + prefix: bcftools +<plugin_name> <args>

    To list available plugins, use:

    bcftools plugin -l

    bcftools plugin <plugin_name>
bcftools +<plugin_name>
bcftools plugin -l

  6. Compare false positives between two runs using bcftools isec

    develop

    To identify which false positives were removed or newly acquired between two different runs (e.g., _old and _new), you must first bgzip and tabix the VCF files from the bcftools.isec directories. Then, use bcftools isec to intersect the false positive files (0001.vcf).

    Workflow:

    1. Navigate to the bcftools.isec directory of each run.
    2. Compress and index the VCFs (0000.vcf, 0001.vcf, 0003.vcf).
    3. Run bcftools isec on the false positive files.

    Results of intersection:

    • _fp/0000.vcf: False positives removed in the new run.
    • _fp/0001.vcf: New false positives acquired in the new run.
    # 1. Prepare old run files
d=_old/bcftools.isec
(cd $d; bgzip < 0000.vcf > 0000.vcf.gz; tabix -f 0000.vcf.gz; \
     bgzip < 0001.vcf > 0001.vcf.gz; tabix -f 0001.vcf.gz; \
     bgzip < 0003.vcf > 0003.vcf.gz; tabix -f 0003.vcf.gz)

# 2. Prepare new run files
d=_new/bcftools.isec
(cd $d; bgzip < 0000.vcf > 0000.vcf.gz; tabix -f 0000.vcf.gz; \
     bgzip < 0001.vcf > 0001.vcf.gz; tabix -f 0001.vcf.gz; \
     bgzip < 0003.vcf > 0003.vcf.gz; tabix -f 0003.vcf.gz)

# 3. Intersect false positives
bcftools isec -p _fp _{old,new}/0001.vcf.gz

  7. Run mpileup evaluation scripts

    develop

    The benchmark provides two primary evaluation scripts: run_mpileup_HG002.sh and run_mpileup_HG005.sh.

    • Tuning/Development: Use run_mpileup_HG002.sh for rapid turnaround tuning and parameter setting. Avoid using chromosomes for tuning that you intend to use for final evaluation.
    • Final Evaluation: Use run_mpileup_HG005.sh with a different chromosome (e.g., chr20) to avoid over-fitting to specific regions or samples.

    Usage Examples:

    Training on a small region:

    ./run_mpileup_HG002.sh pb_50x.bam chr1:10000000-20000000 _pb -X pacbio-ccs

    Final evaluation with a specific BCFTOOLS version:

    BCFTOOLS=bcftools.devel ./run_mpileup_HG005.sh illumina_300x.bam chr20:20000000-21000000 _i_dev -L999 -X illumina
    ./run_mpileup_HG005.sh illumina_300x.bam chr20:20000000-21000000 _i_new -L999 -X illumina
    ./run_mpileup_HG002.sh pb_50x.bam chr1:10000000-20000000 _pb -X pacbio-ccs

  8. Generate and visualize VCF statistics with bcftools stats

    develop

    The bcftools stats command parses VCF or BCF files to produce text statistics suitable for machine processing and plotting.

    Usage Modes:

    • Single File: Prints stats for non-reference allele frequency, depth distribution, quality, per-sample counts, singleton stats, etc.
    • Two Files: Generates separate stats for the intersection and the complements. It also prints concordance (Genotype concordance by non-reference allele frequency, Genotype concordance by sample, Non-Reference Discordance) and correlation.

    Key Options:

    • --af-bins <LIST|FILE>: A comma-separated list of allele frequency bins (e.g., 0.1,0.5,1) or a file with one bin per line.
    • --af-tag <TAG>: The INFO tag used for binning (defaults to AC/AN or GT).
    • -1, --1st-allele-only: Considers only the first alternate allele at multiallelic sites.
    • -I, --split-by-ID: Separates stats for "known sites" (ID column set) and "novel sites" (ID column is ".").
    • -E, --exons <file.gz>: A BGZF-compressed, tabix-indexed file with exons (CHR, FROM, TO) for indel frameshift statistics.
    • -F, --fasta-ref <ref.fa>: A faidx-indexed reference file to determine INDEL context.

    Visualization Workflow:

    1. Generate stats: bcftools stats -s - < input.vcf.gz > file.vchk
    2. Plot stats: plot-vcfstats -p <outdir> file.vchk
    # Generate the stats
bcftools stats -s - > file.vchk

# Plot the stats
plot-vcfstats -p outdir file.vchk

  9. Annotate VCF files with bcftools annotate

    develop

    Use bcftools annotate to add or remove annotations from VCF/BCF files. You can use annotation files in VCF, BED, or tab-delimited formats.

    Key Options:

    • -a, --annotations <file>: The annotation source file.
    • -c, --columns <list>: Defines which columns/tags to carry over. Use special prefixes for merge logic:
      • TAG: Overwrite existing values (unless source is .).
      • +TAG: Add only if existing value is ..
      • .TAG: Overwrite even if source is ..
      • .+TAG: Add new tag, never overwrite existing.
      • \-TAG: Overwrite existing, never add if target doesn't exist.
      • =TAG: Append to existing values.
    • -x, --remove <list>: Remove specific annotations (e.g., INFO, FORMAT/GT).
    • -i, --include <EXPRESSION>: Include only sites matching an expression.
    • -e, --exclude <EXPRESSION>: Exclude sites matching an expression.
    • -m, --mark-sites <TAG>: Mark sites present/absent in the annotation file with a new INFO flag.
    # Remove three fields
bcftools annotate -x ID,INFO/DP,FORMAT/DP file.vcf.gz

# Remove all INFO fields and all FORMAT fields except for GT and PL
bcftools annotate -x INFO,^FORMAT/GT,FORMAT/PL file.vcf

# Add ID, QUAL and INFO/TAG, not replacing TAG if already present
bcftools annotate -a src.bcf -c ID,QUAL,+TAG dst.bcf

# Carry over all INFO and FORMAT annotations except FORMAT/GT
bcftools annotate -a src.bcf -c INFO,^FORMAT/GT dst.bcf

  10. Use bcftools CLI commands

    develop

    The bcftools command-line tool is used for variant calling and manipulating VCF (Variant Call Format) and BCF (Binary VCF) files. The general usage pattern is:

    bcftools <command> <argument>

    Most commands accept VCF, bgzipped VCF, and BCF files. File types are detected automatically, including when streaming from a pipe. Indexed files work in all situations, while un-indexed files and streams work in most but not all situations.

    bcftools <command> <argument>

  11. Perform variant calling with bcftools call

    develop

    The bcftools call command is used for SNP/indel calling from the output of samtools mpileup. It supports two main models:

    • Multiallelic calling model (-m/--multiallelic-caller): Recommended for most tasks.
    • Consensus caller model (-c/--consensus-caller): The older calling model.

    Ploidy and Sex Configuration:

    • --ploidy <ASSEMBLY>: Use predefined ploidy assemblies.
    • --ploidy-file <FILE>: Use a custom ploidy definition file (format: CHROM FROM TO SEX PLOIDY).
    • For bcftools call -C <trio>, a PED file is expected where the last column indicates sex (1=male, 2=female).

  12. Download benchmark data using get_data.sh

    develop

    The get_data.sh script downloads necessary files for the mpileup benchmark, including truth sets, BED files, and BAM files. Note that these files can be very large. You can modify the script to use different URLs or download files manually. It is also recommended to obtain a copy of the GRCh38 reference genome, though differences in patches may not affect specific chromosomes like chr1 or chr20.

    ./get_data.sh