# How to Convert VCF to BED: A Practical Guide By the end of this guide you will be able to take a variant call file, decide deliberately which coordinate convention you need, run the conversion with a command you understand, sort and index the result, and prove that the output is correct before it reaches a downstream tool. You will also recognize the failure modes that produce an empty intersection or a silent one-base shift, because those two problems account for most of the support requests about this conversion. Have these on hand before starting. A bgzipped and indexed VCF (`bcftools`, `tabix`). `bedtools` for intersection testing. The reference FASTA that matches the genome build of your VCF if you plan to normalize indels. Optionally BEDOPS for its `vcf2bed` utility and Python 3 with `pysam` if you prefer a scripted approach. All commands below use short flags so they paste cleanly into a terminal. VCF is described in the literature as the lingua franca of next-generation sequencing, which is a fair summary of its position in genomics workflows [1]. BED is the lingua franca of genomic intervals. The two formats describe overlapping information with different coordinate systems, and the conversion is simple only if you make the coordinate decision explicitly. ## What VCF and BED Actually Store ### VCF: one record per variant A VCF file is tab-delimited text. Metadata lines begin with two hash characters. The column header line begins with a single hash followed by `CHROM`. After the header, each line is one variant record with eight fixed fields, an optional FORMAT field, and one field per sample. | Field | Meaning | Coordinate note | |--|--|--| | CHROM | Chromosome or contig name | Must match the BED naming exactly | | POS | Position of the first base of REF | 1-based | | ID | Identifier, often a dbSNP rs number | A single dot when absent | | REF | Reference allele | Length defines the reference span | | ALT | Alternate allele or comma-separated alleles | Symbolic alleles use angle brackets | | QUAL | Quality score | Not an interval property | | FILTER | Filter status | `PASS` or a filter name | | INFO | Semicolon-separated key=value pairs | `END` matters for structural variants | Three properties of VCF cause most conversion mistakes. POS is 1-based. POS refers to the first base of the reference allele, not to the center of the variant. And a single record can carry several alternate alleles at once, written `G,T` in the ALT field, so one line does not always mean one variant. VCF also carries sample-level genotypes, phasing, and depth. Those are the parts that do not survive a conversion to BED, and they are the reason you should keep the original VCF as the source of truth. VCF is a rich format that has been used as the input side of many downstream conversions, including clinical conversions into HL7 FHIR Genomics reports where each VCF row becomes a structured variant record [1]. The BED format has no equivalent capacity. ### BED: one record per interval BED is also tab-delimited. It requires three columns and permits up to twelve. | Column | Name | Meaning | |--|--|--| | 1 | chrom | Chromosome or contig name | | 2 | chromStart | 0-based start, inclusive | | 3 | chromEnd | 0-based end, exclusive | | 4 | name | Free text label | | 5 | score | 0 to 1000 by convention | | 6 | strand | `+`, `-`, or `.` | | 7 to 12 | thickStart, thickEnd, itemRgb, blockCount, blockSizes, blockStarts | Display and structure fields | The single most important rule in this entire article is that BED counts from zero and excludes its end coordinate. A variant at 1-based position 100 occupies the interval from 0-based 99 up to but not including 100. It is one base wide. Writing `chr1 100 100` produces a zero-width interval that most tools silently ignore, and writing `chr1 100 101` shifts the variant one base to the right. ### Where FASTQ and FASTA fit These four formats sit at different points in a [sequencing pipeline](/blog/guides/sequencing-pipeline), which is why confusing them causes trouble. FASTQ stores raw reads in four lines per record: a header beginning with `@`, the [nucleotide sequence](/knowledge/molecular-biology/nucleotide-sequence), a separator line beginning with `+`, and a per-base quality string in Phred+33 encoding. FASTA stores sequences in two lines per record: a header beginning with `>` and the sequence itself, with no quality information at all. FASTA is the usual container for reference genomes and protein sequences. FASTQ is the usual container for instrument output. VCF and BED come later. A VCF records positions where a sample differs from the reference. A BED records arbitrary intervals with no sequence and no alleles. Neither can be reconstructed from the other. Reading a FASTQ file as if it were FASTA silently folds quality strings into sequence records, which is exactly the class of error that well-built format libraries are designed to prevent by validating structure rather than trusting the extension [2]. ## The Coordinate Decision That Determines Everything ### One base, worked two ways Take this VCF record: ``` chr1 100 rs1234 A G 60 PASS DP=30 ``` Read it as follows. On chromosome 1, at 1-based position 100, the reference base is A and one alternate base is G. The variant covers exactly one reference base. The equivalent BED interval is: ``` chr1 99 100 ``` The arithmetic is `start = POS - 1` and `end = start + length(REF)`, which gives `99 + 1 = 100`. Now take a deletion: ``` chr1 100 . ATG A 50 PASS . ``` The reference allele is three bases long. The variant covers 1-based positions 100, 101, and 102. The equivalent BED interval is: ``` chr1 99 102 ``` The arithmetic is `start = 100 - 1 = 99` and `end = 99 + 3 = 102`. This is the reference span, and it is the convention that `bcftools` reports with its `%END` placeholder. ### Four interval conventions and when each applies There is no universally correct answer to "how wide is a variant in BED." There are four conventions in common use, and each is right for a particular question. | Convention | BED start | BED end | Typical use | |--|--|--|--| | Point interval | POS - 1 | POS | One base per record, SNP-centric overlap counting | | Reference span | POS - 1 | POS - 1 + length(REF) | Gene overlap, region annotation, matching `bcftools %END` | | Alternate span | POS - 1 | POS - 1 + length(ALT) | Insertion-aware overlap, sequence context windows | | INFO/END span | POS - 1 | INFO/END | Structural variants and symbolic alleles | The point interval is the simplest and the most common in tutorials. It is fine when your VCF contains only single nucleotide variants. It becomes misleading for indels, because a three-base deletion represented as a one-base interval will miss annotation features that overlap the deleted bases. The reference span is the safest general-purpose choice. It matches what `bedtools` assumes when it reads a VCF directly, and it matches the `%END` placeholder in `bcftools query`. The alternate span matters when you care about inserted sequence rather than deleted sequence. For an insertion such as `REF=A, ALT=ATTT`, the reference span is one base while the alternate span is four bases. Whichever convention you choose, write it down in the methods section of your analysis. A BED file carries no record of how it was built. ### Multiallelic records A record like this: ``` chr1 100 . A G,T 45 PASS . ``` describes two alternate alleles at one position. Converting it directly yields one interval for the reference allele, which does not distinguish the two alternates. If you need one interval per allele, split first: ```bash bcftools norm -m -any -f reference.fa -Oz -o split.vcf.gz input.vcf.gz bcftools index -t split.vcf.gz ``` The `-m -any` option splits multiallelic sites into biallelic records. After splitting, each record has a single ALT and the interval arithmetic becomes unambiguous. ### Chromosome naming and record order Two BED files can look identical and share zero intervals because one uses `chr1` and the other uses `1`. This produces no error message. `bedtools intersect` returns an empty file and exits successfully, which is the most expensive kind of failure. Check the naming in the VCF header before converting: ```bash bcftools view -h input.vcf.gz | grep '^##contig' | head -5 ``` The output lists the contig names the caller used. Match those names in every BED file you plan to intersect against. If you must rename, do it once, early, and keep the map file. Record order matters too. Most interval tools assume coordinate-sorted input. Sorting with plain `sort` on the whole line is wrong because it sorts the start column as text. Use the numeric key form: ```bash sort -k1,1 -k2,2n variants.bed > variants.sorted.bed ``` ## How to Convert VCF to BED Step by Step ### Step 1: Inventory the input Never convert a file you have not looked at. ```bash bcftools view -h input.vcf.gz | head -30 bcftools view -H input.vcf.gz | head -3 bcftools index -n input.vcf.gz ``` The first command prints the metadata header. The second prints the first three records so you can see whether IDs are present and whether ALT contains commas. The third prints the total number of records, which you will use later as a checksum. A typical header fragment looks like this: ``` ##fileformat=VCFv4.2 ##contig= ##contig= ##INFO= ``` Seeing the `##contig` lines confirms both the naming scheme and the presence of length metadata. Seeing an `##INFO= variants.bed ``` Three placeholders do the work. `%CHROM` prints the contig name. `%POS0` prints the position converted to 0-based, which is exactly the BED start. `%END` prints the end of the reference allele, which is the 0-based exclusive end. The tab escapes inside the format string keep the output column-delimited. The first lines of the output look like this: ``` chr1 99 100 chr1 149 151 chr1 299 300 chr2 499 499 ``` The first record is a single-base variant at 1-based 100. The second is a two-base reference span starting at 1-based 150, so it covers 150 and 151. The last record has identical start and end, which indicates a zero-width or malformed entry and deserves investigation before you proceed. For a point-interval conversion instead, replace `%END` with `%POS`: ```bash bcftools query -f '%CHROM\t%POS0\t%POS\n' normalized.vcf.gz > variants.point.bed ``` ### Step 4: Add name and score columns A three-column BED is often enough, but annotation tools frequently want a label. Add the ID: ```bash bcftools query -f '%CHROM\t%POS0\t%END\t%ID\n' normalized.vcf.gz > variants.bed ``` Where the ID is a dot, substitute a constructed name from the coordinates: ```bash bcftools query -f '%CHROM\t%POS0\t%END\t%ID\t%QUAL\n' normalized.vcf.gz \ | awk 'BEGIN{OFS="\t"} { if ($4 == ".") $4 = $1 ":" ($2 + 1); print }' \ > variants.named.bed ``` The `$2 + 1` rebuilds the 1-based position so the name reads like `chr1:100`, which is what most genome browsers display. The QUAL field lands in the score column, which is legal only if your qualities fall between 0 and 1000. Check that before you rely on it. ### Step 5: Sort and index ```bash sort -k1,1 -k2,2n variants.named.bed > variants.sorted.bed bgzip -c variants.sorted.bed > variants.sorted.bed.gz tabix -p bed variants.sorted.bed.gz ``` The `-p bed` preset tells tabix that the file uses BED coordinates and that the first three columns are chromosome, start, and end. The result is `variants.sorted.bed.gz` and `variants.sorted.bed.gz.tbi`. The index enables region queries later. ### Step 6: Validate before use Four checks catch nearly every error. ```bash wc -l variants.sorted.bed bcftools index -n normalized.vcf.gz ``` These two numbers should match if you emitted one interval per record. A mismatch after splitting means some records were dropped, and a mismatch before splitting means you collapsed multiallelic sites. ```bash awk 'BEGIN{OFS="\t"} $3 <= $2 {print "zero or negative width:", $0}' variants.sorted.bed | head ``` Any output from this command means an interval has no width. Zero-width records are almost always a coordinate bug. ```bash bcftools view -H -r chr1:100-100 normalized.vcf.gz grep -P '^chr1\t99\t100\t' variants.sorted.bed ``` The first command prints the VCF record at 1-based position 100. The second prints the matching BED lines. Confirm they refer to the same variant. Note the shift in the region syntax: `bcftools view -r` takes 1-based inclusive coordinates while BED stores 0-based half-open coordinates. ```bash bedtools intersect -a variants.sorted.bed -b genes.bed -wa | wc -l ``` A non-zero overlap count against a known gene model confirms that chromosome naming and sorting are compatible. Zero overlaps against a gene set you know should overlap is the signal to check the `chr` prefix. ### Step 7: Use the result ```bash bedtools intersect -a variants.sorted.bed -b genes.bed -wa -wb > variants.in.genes.bed ``` The `-wa` flag writes the original variant interval. The `-wb` flag appends the overlapping gene interval. The output is a joined table you can pass to downstream counting or annotation. ## Command Cookbook | Task | Command | |--|--| | VCF to BED, reference span | `bcftools query -f '%CHROM\t%POS0\t%END\n' in.vcf.gz > out.bed` | | VCF to BED, point interval | `bcftools query -f '%CHROM\t%POS0\t%POS\n' in.vcf.gz > out.bed` | | Plain awk on an uncompressed VCF | `awk 'BEGIN{OFS="\t"} !/^#/ {print $1, $2-1, $2}' in.vcf > out.bed` | | BEDOPS conversion | `vcf2bed < in.vcf | sort-bed - > out.bed` | | Rename `chr` prefix | `awk 'BEGIN{OFS="\t"} /^#/ {next} {sub(/^chr/, "", $1); print}' in.vcf` | | Sort a BED | `sort -k1,1 -k2,2n in.bed > out.bed` | | Compress and index | `bgzip -c in.bed > in.bed.gz && tabix -p bed in.bed.gz` | | Count overlaps per variant | `bedtools intersect -a variants.bed -b regions.bed -c` | The awk approach is a useful fallback when `bcftools` is unavailable, but it assumes the file is uncompressed and that POS is the second column, which is true for VCF and false for VCF-like tables exported from other software. It also ignores the reference allele length, so it only produces point intervals. The BEDOPS route is fast and memory-efficient for large files. Piping through `sort-bed` guarantees the coordinate sort that BEDOPS tools expect. The `pysam` route gives you programmatic control: ```python import pysam vcf = pysam.VariantFile("normalized.vcf.gz") with open("variants.bed", "w") as out: for rec in vcf: label = rec.id if rec.id else "." out.write(f"{rec.chrom}\t{rec.start}\t{rec.stop}\t{label}\n") ``` In `pysam`, `rec.start` is the 0-based start of the reference allele and `rec.stop` is the exclusive end. Those are precisely the BED start and end, so no manual arithmetic is needed. This is the advantage of using a library that models the format rather than treating it as text. Libraries built around validated format models exist across bioinformatics, from quantum chemistry to proteomics, and they exist precisely because ad hoc text munging breaks at scale [2][3]. ## Working With Large VCFs `bcftools query` streams records and does not load the file into memory, so a multi-gigabyte VCF converts with modest RAM. For repeated work, convert once to BCF, the binary counterpart of VCF, and query that: ```bash bcftools view -Ob -o input.bcf input.vcf.gz bcftools index input.bcf bcftools query -f '%CHROM\t%POS0\t%END\n' input.bcf > variants.bed ``` Random access requires bgzip compression plus a tabix or CSI index. A plain gzip file cannot be queried by region. Writing a BED without sorting first also wastes time later, because virtually every interval tool either requires sorted input or re-sorts internally. When the conversion is one step in a longer pipeline, wrapping it in a workflow engine pays off. Workflow libraries provide job automation, reproducibility, and audit traces for every output file, which matters when a coordinate convention decision needs to be documented alongside the results [4]. Dedicated toolkits also tend to outperform generic converters as data volume grows, as benchmarks for Hi-C format conversion have shown [5]. The same principle applies here: for a one-off conversion, `bcftools query` is enough, and for a production pipeline, an explicit, tested step is worth the setup. ## Common Mistakes and Limitations **Forgetting the minus one.** The single most common error. VCF POS is 1-based and BED start is 0-based. Every conversion command in this guide handles it, but hand-written awk that prints `$2` in the start column does not. **Writing POS in both start and end columns.** `chr1 100 100` has zero width. Most tools drop it without complaint. The correct point interval is `chr1 99 100`. **Using POS as the end after subtracting one from the start.** `chr1 99 101` is two bases wide and shifts every variant one position to the right. Compute the end as `start + length(REF)`, not as `POS + 1`. **Mismatched chromosome naming.** `chr1` and `1` do not match. `bedtools intersect` returns an empty file with exit status zero. There is no warning. Check the `##contig` header first. **Unsorted input.** `bedtools` with its sorted mode reports that the file has an out of order record and names the offending line. Fix with `sort -k1,1 -k2,2n`, not with a plain `sort`. **Spaces instead of tabs.** A BED file delimited by spaces fails with a malformed BED entry error and an illegal field count. Set `OFS="\t"` in awk or use `bcftools query`, which emits tabs. **Collapsed multiallelic sites.** One record with two ALTs becomes one interval. If your downstream analysis counts alleles rather than sites, split with `bcftools norm -m -any` first. **Symbolic alleles producing one-base intervals.** A record with `ALT=` has a one-character REF allele, so the reference-span calculation yields a single base. Use the INFO/END field instead: ```bash bcftools query -f '%CHROM\t%POS0\t%INFO/END\n' normalized.vcf.gz > sv.bed ``` Records lacking an END value will fail or emit a placeholder. Filter them out first or handle them separately. **Indexing the wrong file.** Running `tabix` on a plain gzip file or on an uncompressed BED produces an htslib open failure. Compress with `bgzip`, not `gzip`. **Confusing BED coordinates with region syntax.** BED is 0-based half-open. `bcftools view -r` and `samtools view` take 1-based inclusive regions. A BED interval `chr1 99 100` is `chr1:100-100` in region syntax, not `chr1:99-100`. **Expecting a round trip.** BED cannot store REF, ALT, genotype, phase, quality, or filter status. Converting BED back to VCF produces a syntactically valid file with no allele information, which is worse than useless for most analyses because it looks like data. Keep the VCF. **Assuming the conversion is lossless for structural variants.** Breakend notation, symbolic alleles, and multi-allelic structural records do not have a clean interval representation. Treat these separately and document what you did. **Formatting errors that survive to the next step.** Conversion tools that skip validation push the failure downstream, where it is harder to diagnose. GUI wrappers around converters reduce typing errors but do not eliminate the need to inspect output, a point made explicit by tools that convert aligned FASTA into population-genetics input formats [6]. Standardized, human-readable formats with explicit schemas are easier to validate than ad hoc text, which is the argument made for newer formats in mass spectrometry and beyond [7]. Applying the same discipline here means checking counts, widths, and a handful of spot coordinates every time. One more limitation worth stating plainly. Interoperability between formats is an active engineering problem, not a solved one, and even bidirectional converters between well-documented formats require careful handling of proprietary extensions [8]. VCF to BED is simpler than most such problems because both formats are open and well specified, but the coordinate mismatch is a real semantic gap. Scripted conversions that move data between formats need explicit tests, and that is true whether the target is a sequence submission checklist [9] or an interval file [10]. ## Frequently Asked Questions ### What is the difference between VCF POS and BED start? VCF POS is 1-based and points at the first base of the reference allele. BED chromStart is 0-based. For a variant at POS 100, the BED start is 99. ### Does converting VCF to BED lose genotype information? Yes, completely. BED has no fields for genotype, allele depth, phase, or filter status. Only coordinates and optional labels survive, so keep the VCF as the authoritative record. ### Should I subtract 1 from POS when making a BED file? Yes, always. Subtract 1 for the start and compute the end as start plus the length of the reference allele. Commands like `bcftools query -f '%POS0'` do the subtraction for you. ### How do I convert only SNPs to BED? Filter first, then convert. Use `bcftools view -v snps input.vcf.gz` to keep only single nucleotide variants, then run the standard query. Filtering by variant type also removes the indel ambiguity from the coordinate decision. ### Why does bedtools intersect return zero overlaps with my VCF-derived BED? The usual cause is a chromosome naming mismatch, such as `chr1` in one file and `1` in the other. The second most common cause is inconsistent sort order. Neither produces an error message. ### Can I convert BED back to VCF? Only syntactically. BED carries no REF or ALT alleles, so a reconstructed VCF has no variant identity. Round-tripping through BED is a data-loss operation, not a format conversion. ### What is the difference between FASTQ and FASTA? FASTQ stores reads in four lines with a per-base quality string. FASTA stores sequences in two lines with no quality data. FASTQ is instrument output and FASTA is typically a reference or an assembled sequence. ### How do I convert a VCF with structural variants to BED? Use the INFO/END field rather than the reference allele length. Query with `%INFO/END` and handle records that lack it separately, since symbolic alleles do not have a reference span that reflects their true size. ## Sources 1. [vcf2fhir: a utility to convert VCF files into HL7 FHIR format for genomics-EHR integration.](https://pubmed.ncbi.nlm.nih.gov/33653260/) 2. [IOData: A python library for reading, writing, and converting computational chemistry file formats and generating input files.](https://pubmed.ncbi.nlm.nih.gov/33368350/) 3. [ThermoRawFileParser: Modular, Scalable, and Cross-Platform RAW File Conversion.](https://pubmed.ncbi.nlm.nih.gov/31755270/) 4. [SciPipe: A workflow library for agile development of complex and dynamic bioinformatics pipelines.](https://pubmed.ncbi.nlm.nih.gov/31029061/) 5. [hictk: blazing fast toolkit to work with .hic and .cool files.](https://pubmed.ncbi.nlm.nih.gov/38913844/) 6. [Fasta2Structure: a user-friendly tool for converting multiple aligned FASTA files to STRUCTURE format.](https://pubmed.ncbi.nlm.nih.gov/38365590/) 7. [HXMS: a standardized file format for HX-MS data.](https://pubmed.ncbi.nlm.nih.gov/41803064/) 8. [cd2sbgnml: bidirectional conversion between CellDesigner and SBGN formats.](https://pubmed.ncbi.nlm.nih.gov/31904823/) 9. [EMBL2checklists: A Python package to facilitate the user-friendly submission of plant and fungal DNA barcoding sequences to ENA.](https://pubmed.ncbi.nlm.nih.gov/30629718/) 10. [anndataR improves interoperability between R and Python in single-cell transcriptomics.](https://pubmed.ncbi.nlm.nih.gov/42108549/)