How to Create a BED File from a BAM File: Extracting Read Intervals

By Dr. Zubair Khalid, DVM, MS, PhD · 2026-11-10

By the end of this guide you will have a coordinate-sorted, validated BED file built from a BAM alignment file, and you will know how to produce the three BED variants that most workflows need: BED3 for simple intervals, BED6 for named reads with scores and strands, and BED12 for spliced alignments. You will also know how to filter alignments before conversion, how to catch off-by-one coordinate errors, and how to check that the number of BED lines matches the number of alignments you expected.

You need four things on hand before you start: a BAM file, a matching BAM index, samtools, and BEDTools. A reference FASTA is optional for BAM input but required if you work with CRAM or if you want chromosome sizes for coverage tracks. The commands below use only standard tools that ship with samtools and BEDTools, plus the POSIX utilities sort, awk, cut, and wc. Install samtools and BEDTools from Conda, Homebrew, or your distribution package manager. Confirm both are on your PATH by running each with no arguments, which prints usage and version information.

What BED, BAM, and Read Intervals Actually Mean

A BED file is a plain-text tab-delimited file that describes genomic intervals. The name stands for Browser Extensible Data, and the format is specified by the UCSC Genome Browser team. The first three columns are mandatory: chromosome, start, and end. The start coordinate is 0-based, meaning the first base of a chromosome is position 0. The end coordinate is exclusive, so an interval that covers bases 1 through 100 in 1-based counting is written as start 0 and end 100. This half-open convention is the single largest source of confusion when people create a BED file from a BAM file.

A BAM file is a binary, BGZF-compressed version of the Sequence Alignment/Map format. The SAM/BAM specification is maintained by the Global Alliance for Genomics and Health. BAM stores read sequences, quality strings, alignment coordinates, CIGAR strings, mapping qualities, and flags. Unlike BED, SAM and BAM use 1-based inclusive coordinates. A read whose leftmost aligned base is reference position 1001 has a SAM POS field of 1001. The equivalent BED start is 1000. Conversion tools handle this shift, but only if you use them correctly.

A read interval is the span on the reference genome that a sequencing read covers after alignment. For a simple 100 bp read with a 100M CIGAR string aligned at POS 1001, the interval is start 1000 and end 1100 in BED coordinates. For a spliced RNA-seq read with a CIGAR string that includes N operations, the interval may span a large distance even though the read covers only a small amount of sequence. For paired-end data, each mate is a separate alignment record, so each mate produces its own interval unless you deliberately collapse the pair into a fragment.

Several terms appear repeatedly in this workflow. A CIGAR string is the compact alignment description that records matches (M), insertions (I), deletions (D), skips or introns (N), soft clips (S), and hard clips (H). Mapping quality, or MAPQ, is a phred-scaled confidence that the alignment position is correct. A flag is a bitwise field in SAM that records whether a read is paired, unmapped, reverse-complemented, a secondary alignment, a duplicate, or supplementary. The UCSC BED specification defines the optional columns that appear after the first three, including name, score, strand, thick start, thick end, itemRgb, block count, block sizes, and block starts.

BED columns at a glance

BED3 contains only the first three columns. BED6 adds name, score, and strand. BED12 adds the six block-related columns and is used for spliced alignments and transcript models. BEDTools uses these same conventions, and its bamtobed subcommand is the standard way to create a BED file from a BAM file.

FASTQ versus FASTA in this workflow

New users sometimes confuse the raw read formats with the alignment formats. FASTQ stores sequence and per-base quality in a four-line record per read. FASTA stores sequence only, usually with a header line and wrapped sequence lines, and it is the format used for reference genomes. A BAM file is produced by aligning reads from FASTQ files to a reference FASTA. BED is downstream of that alignment and stores intervals, not sequences or qualities. If you need to recover reads from a BAM, samtools fastq writes FASTQ and samtools fasta writes FASTA. If you need a BED file, you convert the alignment records, not the reads.

Before You Start: Inputs, Index, and a Sanity Check

The BAM must be readable and, for most downstream operations, coordinate-sorted and indexed. The index is a separate file, usually sample.bam.bai or sample.bai, that lets tools jump to a genomic region without scanning the whole file. BEDTools bamtobed can stream an unsorted BAM, but sorting and indexing first makes every later step faster and catches truncated files early.

Run these checks before conversion:

samtools quickcheck -v sample.bam
samtools idxstats sample.bam | head
samtools flagstat sample.bam

samtools quickcheck -v prints nothing and exits with status 0 when the file has a valid header and is not truncated. A truncated or corrupt BAM produces a message such as [E::hts_open_format] Failed to open file "sample.bam" or [main_samview] fail to read the header. samtools idxstats prints one line per reference sequence plus a line for unmapped reads. The columns are reference name, reference length, number of mapped reads, and number of unmapped reads. If the command fails with [E::idx_find_and_load] Could not retrieve index file for 'sample.bam', create the index:

samtools index sample.bam

samtools flagstat summarizes total reads, mapped reads, properly paired reads, duplicates, secondary alignments, supplementary alignments, and singleton reads. Record these numbers. You will compare them with the BED line count later. If the BAM is not coordinate-sorted, sort it first:

samtools sort -@ 4 -o sample.sorted.bam sample.bam
samtools index sample.sorted.bam

The -@ 4 option requests four additional threads for compression. The output file is coordinate-sorted, and the index is written beside it.

How to Create a BED File from a BAM Step by Step

The walkthrough below produces a filtered BED6 file from a coordinate-sorted BAM. Each step includes the command, the expected output, and the check that tells you the step worked.

Step 1: Inspect the header and confirm chromosome names

samtools view -H sample.sorted.bam | head -n 5

Expected output looks like this, with your own reference names and lengths:

@HD	VN:1.6	SO:coordinate
@SQ	SN:chr1	LN:248956422
@SQ	SN:chr2	LN:242193529

The SN fields are the chromosome names that will appear in column 1 of your BED file. They must match the names in any BED, GTF, or genome file you plan to intersect with. A BAM that uses 1, 2, and 3 will not intersect a BED that uses chr1, chr2, and chr3. Check this before you convert, not after.

samtools view -H sample.sorted.bam | awk '$1=="@SQ" {print $2}' | cut -d: -f2 | sort -u

This prints the unique reference names in the BAM header.

Step 2: Decide which intervals you need

The same BAM can produce several different BED files. Choose before you run the conversion because the options differ.

For ChIP-seq, ATAC-seq, and most DNA-seq work, read alignments as BED6 are the right default. For RNA-seq, spliced blocks or BED12 are usually what you want. For fragment-level analysis, collapse paired reads before you count anything.

Step 3: Convert the BAM to BED6

bedtools bamtobed -i sample.sorted.bam > sample.bed
head -n 3 sample.bed

Expected output:

chr1	1000	1100	SRR12345.1	60	+
chr1	1050	1150	SRR12345.2	60	-
chr1	2000	2080	SRR12345.3	42	+

The columns are chromosome, start, end, read name, MAPQ, and strand. bamtobed reads the CIGAR string and converts the alignment span to BED coordinates. Soft-clipped bases are not included in the interval. Deletions and skips are included in the span unless you use -split. The score column contains the mapping quality, which BEDTools copies from the BAM MAPQ field.

A quick check that the file has the expected number of columns:

awk 'NF != 6 {print NR, NF, $0}' sample.bed | head

No output means every line has six fields. Any output identifies the line number and the field count.

Step 4: Filter alignments before conversion

Most BAM files contain reads you do not want in a BED file: unmapped reads, secondary alignments, supplementary alignments, duplicates, and low-MAPQ alignments. Filtering before conversion is faster and cleaner than filtering the BED afterward.

samtools view -b -q 20 -F 0x904 sample.sorted.bam > sample.filtered.bam
samtools index sample.filtered.bam
bedtools bamtobed -i sample.filtered.bam > sample.filtered.bed

The -q 20 option keeps alignments with MAPQ of at least 20. The -F 0x904 option excludes reads with any of three flags: 0x4 unmapped, 0x100 secondary, and 0x800 supplementary. The sum 0x904 is the combination of those three bits. If you also want to remove duplicates, use -F 0x904 -F 0x400 only if your samtools version accepts repeated -F options. The portable way is to combine the bits yourself: 0x904 plus 0x400 equals 0xD04, so the command becomes:

samtools view -b -q 20 -F 0xD04 sample.sorted.bam > sample.filtered.bam

Check the flag definitions against the SAM specification before you use a combined value. The common flags are 0x1 paired, 0x2 properly paired, 0x4 unmapped, 0x8 mate unmapped, 0x10 reverse strand, 0x20 mate reverse strand, 0x40 first in pair, 0x80 second in pair, 0x100 secondary, 0x200 quality control failure, 0x400 duplicate, and 0x800 supplementary. A flag filter is a bitmask. Any read whose flag shares a set bit with the mask is excluded.

Step 5: Sort the BED file

Many BEDTools operations require sorted input. Sort by chromosome and then by numeric start:

LC_ALL=C sort -k1,1 -k2,2n sample.filtered.bed > sample.sorted.bed

The LC_ALL=C setting forces byte-order sorting, which keeps chromosome order consistent across machines. The -k1,1 key sorts by chromosome. The -k2,2n key sorts numerically by start within each chromosome. If you skip the n, sort treats starts as text, and 1000 will come before 200. That produces an out-of-order file that bedtools merge or bedtools intersect will reject with a message such as Error: Sorted input specified, but the file sample.bed has the following out of order record.

Check that the sort worked:

bedtools sort -i sample.sorted.bed | cmp -s - sample.sorted.bed
echo $?

A zero exit status means the file is already in BEDTools sort order. A nonzero status means the order differs. The cmp -s command compares silently and returns 0 when the files are identical.

Step 6: Compress and index the BED if you need random access

Plain BED is fine for small files. For large files or for repeated region queries, compress with bgzip and index with tabix:

bgzip -c sample.sorted.bed > sample.sorted.bed.gz
tabix -p bed sample.sorted.bed.gz

bgzip is block gzip, which supports random access. tabix -p bed builds a .tbi index for a BED file. The preset knows the BED coordinate convention, so you do not need to specify 0-based or 1-based settings. Do not use ordinary gzip if you plan to index the file.

Step 7: Validate the result

wc -l sample.sorted.bed
awk '$2 < 0 || $3 <= $2 {print NR, $0}' sample.sorted.bed | head
cut -f1 sample.sorted.bed | sort -u

The first command gives the number of intervals. The second should produce no output. Any line where start is negative or end is less than or equal to start is malformed. The third command lists the chromosomes present. Compare that list with the BAM header list from Step 1. They should match, minus any chromosomes with no alignments.

Step 8: Compare counts with the BAM

samtools view -c -F 0x904 sample.sorted.bam
wc -l sample.sorted.bed

The samtools view -c count and the BED line count should be identical if the BED was made from the same filtered BAM and no alignments were dropped. If the BED has fewer lines, check whether you used -split, which produces multiple BED lines per read, or whether you filtered after conversion. If the BED has more lines, -split is the usual cause.

Command Recipes for Common BED Variants

BED3: minimal intervals

bedtools bamtobed -i sample.sorted.bam | cut -f1-3 > sample.bed3

BED3 is the smallest valid BED file. Use it when you only need intervals for intersection or coverage and do not care about read names or strands.

BED6 with custom score

The default score is MAPQ. To use edit distance instead, add -ed:

bedtools bamtobed -ed -i sample.sorted.bam > sample.editdist.bed

The edit distance is read from the NM tag when it is present. If the tag is missing, the behavior depends on the BEDTools version. Check the output before you rely on it.

BED12 for spliced alignments

bedtools bamtobed -bed12 -i rna.sorted.bam > rna.bed12

BED12 represents each read as one interval with blocks for the aligned segments. The thickStart and thickEnd fields are set to the start and end of the coding region by convention, and the block count, block sizes, and block starts describe the exons. This is the format to use when you want one line per read and you need to preserve the splice structure.

Split blocks for exon-level intervals

bedtools bamtobed -split -bed12 -i rna.sorted.bam > rna.blocks.bed12

The -split option breaks each alignment at N CIGAR operations. A read that spans three exons produces three BED lines instead of one. This is useful for counting reads per exon, but it inflates the line count relative to the number of reads. Do not compare a split BED line count with a samtools view -c count and expect them to match.

Paired-end fragments

bedtools bamtobed -bedpe -i sample.sorted.bam \
  | awk 'BEGIN {OFS="\t"}
    $1 == $4 {
      s = ($2 < $5 ? $2 : $5)
      e = ($3 > $6 ? $3 : $6)
      print $1, s, e, $7, $8, $9
    }' > sample.fragments.bed

BEDPE columns from bamtobed are chrom1, start1, end1, chrom2, start2, end2, name, score, strand1, and strand2. The awk block keeps pairs on the same chromosome and writes a single BED6 interval from the leftmost start to the rightmost end. The name, score, and strand come from the first mate. For strict proper-pair analysis, filter the BAM with -f 0x2 before conversion.

Unique read start sites

cut -f1,2 sample.sorted.bed | sort -u \
  | awk 'BEGIN {OFS="\t"} {print $1, $2, $2 + 1}' > sample.starts.bed

This produces one 1 bp interval per unique start position. It is common for ATAC-seq and DNase-seq footprinting, where you care about cutsites rather than full read spans. Deduplication by start position is a deliberate choice, so document it.

Coverage track

samtools faidx reference.fa
cut -f1,2 reference.fa.fai > reference.sizes
bedtools genomecov -i sample.sorted.bed -g reference.sizes -bg > sample.coverage.bedgraph

The genome file has two columns: chromosome and length. The -bg output is a BEDGraph file with four columns: chromosome, start, end, and depth. BEDGraph is a BED-derived format, not a BED file with six columns. Keep the distinction clear when you pass it to other tools.

Parameters That Matter

The -split option is not the same as -bed12. -bed12 keeps one line per read and describes blocks within that line. -split emits one line per block. Use -bed12 alone when you want a compact representation. Use -split with or without -bed12 when you want each aligned block to be its own interval.

The -mate1 option is a shortcut that reports only the first mate of each pair. It halves the number of intervals but does not give you the full fragment span. If you need fragment coordinates, use BEDPE and collapse the pair as shown above.

Common Mistakes and Limitations

Off-by-one coordinates

The most common mistake is treating BED start coordinates as 1-based. A SAM record with POS 1001 and a 100M CIGAR string becomes BED start 1000 and end 1100. If you write start 1001 and end 1101, every interval shifts one base to the right. The error is small per read but systematic across the file, and it changes overlap counts near interval boundaries. Use bamtobed rather than hand-written awk converters unless you have a specific reason.

Forgetting to filter secondary and supplementary alignments

A BAM can contain multiple alignment records for the same read. Secondary alignments are alternative mappings. Supplementary alignments are parts of a chimeric or split read. If you convert without filtering, the same read can appear several times in the BED file, and coverage or overlap counts inflate. The mask 0x904 removes unmapped, secondary, and supplementary records in one step. Decide explicitly whether you want duplicates removed as well, because duplicate marking is a separate policy choice.

Counting mates as fragments

In paired-end data, each mate is a separate alignment record. A BED6 file made with default bamtobed has two lines per pair. If you then count overlaps or compute coverage, you are counting read ends, not fragments. For fragment-level analysis, collapse the pair with BEDPE or use a proper-pair filter plus a fragment-building step. State which unit you are counting in your methods.

Sorting with the wrong locale or wrong key

sort -k1,1 -k2,2n is the standard BED sort. Omitting the n makes the start column sort as text. Omitting LC_ALL=C can produce different chromosome orders on different machines, especially when chromosome names mix letters and numbers. BEDTools expects chromosome grouping and ascending start positions. A file that looks sorted to a human can still fail the BEDTools order check.

Mixing chromosome naming conventions

A BAM aligned to GRCh38 with Ensembl-style names uses 1, 2, and X. A BED file from UCSC uses chr1, chr2, and chrX. bedtools intersect will return zero overlaps rather than an error when the names do not match. Check the first column of both files before you run an intersection. If the names differ, rename one side deliberately and document the change.

Treating BED score as unrestricted

The UCSC BED specification describes the score column as an integer from 0 to 1000. MAPQ values from 0 to 60 fit comfortably. Edit distances and custom scores may not. If you plan to load the BED in a genome browser, keep the score in range or convert it to a display-friendly value. Downstream command-line tools often ignore the score entirely, so the constraint matters mainly for visualization.

Ignoring spliced reads in RNA-seq

A spliced read with a 50M1000N50M CIGAR string spans 1100 bases on the reference but covers only 100 bases of sequence. If you use the default bamtobed output, the intron is included in the interval. Coverage and overlap calculations then count intronic space as covered. Use -split or -bed12 for RNA-seq. The same issue appears with deletions in DNA-seq, though the spans are usually shorter.

Assuming BED is 1-based because SAM is 1-based

SAM, GFF, and GTF use 1-based coordinates. BED uses 0-based half-open coordinates. VCF uses 1-based coordinates. The conversion between them is not a simple formatting change. When you move data between BED and GFF, add one to the BED start or subtract one from the GFF start, and check the end coordinate convention as well. The UCSC specification is the reference to follow for BED.

Using an unindexed or unsorted BAM for region queries

bamtobed can stream a BAM without an index, but region-restricted samtools commands need an index. If you see [E::idx_find_and_load] Could not retrieve index file, run samtools index on a coordinate-sorted BAM. Indexing an unsorted BAM fails or produces an unusable index. Sort first, then index.

Forgetting the reference for CRAM

CRAM is a reference-based compressed alignment format. Converting CRAM to BED requires the reference FASTA or a reference repository. If the reference is missing, samtools reports a message such as [E::cram_get_ref] Failed to populate reference for id 0. Supply the reference with samtools view -T reference.fa before piping to BEDTools, or set the reference path in the environment according to the samtools documentation. BEDTools bamtobed does not read CRAM directly in all versions, so convert to BAM first when in doubt.

Assuming every BAM line becomes one BED line

The relationship between BAM records and BED lines depends on the options. Default BED6 gives one line per alignment record. -split gives one line per aligned block. -bedpe gives one line per pair. -bed12 gives one line per alignment record with block structure. None of these is wrong, but each answers a different question. Write down the unit before you count.

Overlooking the difference between FASTQ and FASTA inputs

A BAM is built from FASTQ reads aligned to a FASTA reference. If you are troubleshooting a BED file, check both inputs. A reference FASTA with the wrong chromosome names produces a BAM with the wrong names. A FASTQ with adapter contamination produces soft-clipped alignments and intervals that do not reflect the true insert. BED conversion cannot fix upstream problems.

How to Check That the BED File Is Correct

Validation is mostly about counts, coordinates, and chromosome names. Run these checks every time you create a BED file from a BAM.

Line count versus alignment count. Use samtools view -c with the same filter you used for conversion. For default BED6, the counts should match. For -split, the BED count should be greater than or equal to the alignment count.

Coordinate sanity. Check for negative starts and zero-length intervals:

awk '$2 < 0 || $3 <= $2 {print NR, $0}' sample.sorted.bed | head
awk '$3 - $2 > 1000000 {print NR, $0}' sample.sorted.bed | head

The first command should return nothing. The second flags unusually long intervals, which can indicate a spliced read that was not split or a chimeric alignment.

Chromosome list. Compare the unique chromosome names in the BED with the @SQ names in the BAM header. They should be the same set, minus chromosomes with no mapped reads.

Spot check a known read. Pick a read name from the BAM and find it in the BED:

samtools view sample.sorted.bam | grep -m 1 '^SRR12345.1'
grep -m 1 'SRR12345.1' sample.sorted.bed

Compare the SAM POS and CIGAR with the BED start and end. For a simple 100M alignment at POS 1001, the BED line should start at 1000 and end at 1100.

Genome browser check. Load the BED file in IGV or the UCSC Genome Browser with the matching reference genome. Intervals should sit over aligned reads, not offset by one base or shifted to the wrong chromosome. A browser view is a fast visual check for chromosome naming and coordinate errors.

Coverage sanity. Compute the mean depth over a small region and compare it with the expected depth from the library. A sudden tenfold difference usually means duplicate reads, split blocks, or mate double counting.

Practical Tips for Reproducible BED Creation

Record the exact commands and tool versions in a script. BEDTools and samtools behavior can change between major versions, especially around tag handling and split options. A short shell script with the filter mask, MAPQ cutoff, and sort command documented is enough to reproduce the file.

Keep the intermediate BAM. If you filter before conversion, the filtered BAM is the authoritative input for the BED. Saving it lets you regenerate the BED with different options without repeating the alignment filtering.

Use set -euo pipefail at the top of shell scripts so a failed command stops the pipeline. A broken pipe in the middle of a conversion can leave a truncated BED file that looks plausible until you count the lines.

Name files after the unit they contain. sample.reads.bed, sample.fragments.bed, and sample.blocks.bed are clearer than three files called sample.bed. The name prevents the most common downstream analysis error, which is using read-end intervals when fragment intervals were intended.

Keep BED files uncompressed during development and compress only for distribution or archival storage. Plain text is easier to inspect with head, awk, and grep. Once the file is validated, bgzip and tabix add random access without changing the content.

Frequently Asked Questions

What is the difference between a BAM file and a BED file?

A BAM file stores aligned reads with sequence, quality, CIGAR, and flag information. A BED file stores genomic intervals only. BAM is 1-based inclusive, BED is 0-based half-open.

Can I create a BED file from a BAM without BEDTools?

Yes. You can stream samtools view output through awk and convert POS and CIGAR to start and end coordinates. BEDTools bamtobed is safer because it handles CIGAR operations, strand, and BED12 blocks correctly.

Should I use BED6 or BED12?

Use BED6 for simple read alignments and most DNA-seq work. Use BED12 when you need to preserve splice structure or block layout in a single line per read.

How do I convert paired-end reads to fragment intervals?

Use bedtools bamtobed -bedpe, then collapse each pair to the leftmost start and rightmost end with awk. Filter to properly paired reads first if you want only concordant fragments.

Why does my BED start one base before the SAM position?

BED is 0-based and SAM is 1-based. A SAM POS of 1001 becomes a BED start of 1000. The conversion is expected and is handled automatically by bamtobed.

How do I filter low-quality alignments before making a BED file?

Filter the BAM with samtools view -q 20 -F 0x904 before conversion. The -q option sets the minimum MAPQ. The -F mask removes unmapped, secondary, and supplementary alignments.

Can I make a BED file from a CRAM file?

Convert the CRAM to BAM first with samtools view -b -T reference.fa, then run bamtobed. CRAM needs the reference FASTA, and BEDTools does not read every CRAM variant directly.

What is the difference between FASTQ and FASTA when working with BAM files?

FASTQ stores reads with per-base quality scores and is the input to alignment. FASTA stores sequences without quality scores and is used for the reference genome. BAM is the aligned output, and BED is derived from the alignment.