A Real BAM File Example, Field by Field

A BAM file is the BGZF-compressed binary form of SAM, and a single record in it looks like this when you decompress one line with samtools view:

A00627:18:HGV7TDSXX:1:1101:2085:1000	99	chr1	1000210	60	151M	=	1000455	396	GATTACAGG...	FFFFF:FFF...	NM:i:1	MD:Z:75G75	AS:i:146	XS:i:0	RG:Z:sample1.lane1	MC:Z:151M	MQ:i:60

That is one 151 bp read from one cluster on one Illumina flowcell, aligned to GRCh38. A 30x human whole-genome BAM contains roughly 800 million to 1 billion such records. Everything below explains how to read that line and what we check before trusting the file it came from.

What each field in the record means

The first eleven fields are mandatory and tab-delimited. They appear in a fixed order, so you can cut them with awk without parsing anything.

A00627:18:HGV7TDSXX:1:1101:2085:1000 is QNAME, the read name. Both mates of a pair share it, which is how tools collate pairs after sorting by coordinate ¹. The Illumina convention encodes the fields below, so the name itself tells you the read came from lane 1 of flowcell HGV7TDSXX.

• instrument
• run
• flowcell
• lane
• tile
• coordinates

99 is FLAG, a bitfield. 99 = 1 + 2 + 32 + 64. Those bits mean paired, mapped in a proper pair, mate on the reverse strand, and first in pair. Its mate carries 147. Read flags as bits and never as magic numbers, and use samtools flags 99 to decode one you do not recognise.

chr1 and 1000210 are RNAME and POS, the reference contig and the leftmost 1-based coordinate of the alignment. 60 is MAPQ, the phred-scaled probability that the alignment is misplaced. BWA-MEM emits 60 for uniquely placed reads and 0 for reads whose best placement ties another, which is the main signal you have about repeats and mapping ambiguity ².

151M is the CIGAR string: 151 bases matched or mismatched against the reference, with no indels or clipping. Real BAMs are full of strings like 20S131M (soft-clipped) or 75M2D76M (deletion). Sum the M, I, S and = operators to recover read length. Sum M, D, N and = to recover reference span.

=, 1000455 and 396 are RNEXT, PNEXT and TLEN. The = means the mate is on the same contig, the mate starts at 1000455, and the insert spans 396 bp. A TLEN of 0 or an RNEXT naming a different contig is worth noticing, because discordant pairs are the raw material for structural variant calling.

Then come SEQ and QUAL, the read bases and their base qualities as ASCII phred+33. F is Q37 under Illumina’s binned quality scheme. If you see a long run of identical quality characters, the file has been quality-binned, which is normal for NovaSeq output.

Everything after QUAL is an optional tag in TAG:TYPE:VALUE form. NM:i:1 is the edit distance to the reference. MD:Z:75G75 says the mismatch is at read offset 75 against a reference G. AS:i:146 is the alignment score and XS:i:0 the score of the best alternative placement, so the gap between them is what drove MAPQ 60. RG:Z:sample1.lane1 links the read to a read group in the header.

The header tells you which reference and which sample

Before reading alignments, read the header. samtools view -H sample.bam | head -20 gives you something like:

@HD	VN:1.6	SO:coordinate
@SQ	SN:chr1	LN:248956422	M5:2648ae1bacce4ec4b6cf337dcae37816
@RG	ID:sample1.lane1	SM:sample1	LB:lib1	PL:ILLUMINA	PU:HGV7TDSXX.1.CGTACTAG
@PG	ID:bwa	PN:bwa	VN:0.7.17-r1188	CL:bwa mem -t 16 -R ...
@PG	ID:MarkDuplicates	PN:MarkDuplicates	VN:2.27.5

SO:coordinate means the file is sorted by position and can be indexed. LN:248956422 for chr1 identifies GRCh38, and hg19 chr1 is 249250621. That one number resolves the most common source of downstream confusion, because a BED file or VCF built on the wrong assembly will silently produce nonsense. The M5 field is the MD5 of the contig sequence and pins the exact FASTA, including whether alt contigs and decoys were present.

The @RG line carries SM, the sample name that variant callers write into VCF column headers. If SM says sample1 and you expected your own identifier, stop and confirm provenance. The @PG lines record the aligner, its version and its full command line, so you can tell whether duplicates were marked and by which tool.

The first ten minutes with a new BAM

We run the same four commands on every file we receive, in this order.

samtools quickcheck -v sample.bam verifies the header parses and the BGZF end-of-file block is present. It reads only the ends of the file and takes under a second. A truncated download is the single most common defect, and it will otherwise surface hours later as a cryptic error inside a caller.

samtools index -@ 8 sample.bam builds the .bai companion index. BAM is a series of BGZF blocks, each holding at most 64 KB of uncompressed data, and the index maps genomic bins to virtual file offsets so a region query seeks instead of scanning. Use samtools index -c to write a .csi instead when any contig exceeds 512 Mbp, which is the .bai limit.

samtools flagstat -@ 8 sample.bam gives totals in a few minutes. On a clean 30x human library we expect above 99% of reads mapped, above 95% properly paired, and duplicates somewhere between 5% and 20% depending on library input. Duplicates above roughly 30% usually mean low-complexity input rather than a sequencing problem.

samtools idxstats sample.bam reads the index alone and returns per-contig read counts instantly. It is the fastest way to notice that chrY has near-zero coverage or that most reads landed on decoy contigs. For an interactive look at coverage uniformity and insert sizes without downloading anything, bam.iobio streams a sampled subset of the file in a browser ³.

Slicing a region and counting with flags

Random access is the reason to keep BAM coordinate-sorted. To pull the reads over a single gene:

samtools view -b sample.bam chr7:117,480,000-117,670,000 > cftr.bam
samtools view -c -f 2 -F 3840 sample.bam chr7:117,480,000-117,670,000

-f 2 keeps properly paired reads. -F 3840 drops secondary alignments, supplementary alignments, QC-failed reads and duplicates in one mask. The -X flag lets you name the same bits as strings. We prefer that in scripts because -F UNMAP,SECONDARY,DUP is readable a year later.

Rule-based filtering across many regions at once is where dedicated tools help. VariantBam applies region-specific rules in a single pass, which is how you build a small BAM that keeps full detail at your loci of interest and discards the rest ⁴. If you are writing your own tooling rather than shelling out, the BamTools C++ API exposes records as objects and avoids reparsing text ⁵.

For coverage, samtools coverage sample.bam prints mean depth and breadth per contig in one table. Per-base depth over a target uses samtools depth -a -b targets.bed -Q 20 -q 13 sample.bam, where -Q filters mapping quality and -q filters base quality. Report both mean depth and the fraction of bases above 10x, because a 30x mean hides poorly covered regions.

Failure modes we see repeatedly

Reference mismatch is first. A BAM aligned to a chr-prefixed GRCh38 and a VCF built on Ensembl naming will not intersect, and most tools report zero overlaps rather than an error. Compare @SQ names and lengths before anything else.

Sample mixups are second and more consequential. When you hold sequencing from more than one timepoint, confirm that two BAMs are the same person by genotype concordance rather than by filename. BAM-matcher compares a few thousand common SNP sites directly from alignment files and returns an answer in minutes ⁶.

Third is silent truncation and stale indexes. A .bai older than its .bam will return wrong or empty regions, so check modification times when a region query looks suspiciously sparse. Fourth is missing duplicate marking, which inflates apparent depth at PCR-favoured positions and biases allele fractions. The @PG lines tell you whether it happened.

BAM, CRAM, and when to convert

A 30x human BAM runs roughly 60 to 90 GB. CRAM stores the same alignments by reference difference and typically lands near 25 to 30 GB, which matters when you are paying for storage or moving files over a home connection.

samtools view -@ 8 -T GRCh38.fa -C -o sample.cram sample.bam
samtools index sample.cram

We keep CRAM as the archival copy and BAM for active work. The tradeoff is that CRAM needs the exact reference FASTA to decode, so store the FASTA and its .fai alongside the CRAM and record the M5 checksums. Lose the reference and you have lost the reads.

Note that read counts and coverage tell you about data quality. They do not tell you about health. Interpreting a variant you find in these alignments requires confirmation in a clinical laboratory and a conversation with a clinician or genetic counselor.

Questions people also ask

What is the difference between a SAM and a BAM file? SAM is the tab-delimited text form and BAM is the same records encoded in binary and compressed with BGZF. BAM is several times smaller and supports indexed random access, which text SAM cannot.

Can I open a BAM file in a text editor? No. Use samtools view sample.bam | head to see decoded records, or samtools view -h to include the header. Piping raw BAM into less shows compressed bytes.

Do I need the reference FASTA to read a BAM? Not for most operations, because BAM stores read bases explicitly. You need it for CRAM, for samtools mpileup and for anything that compares reads to reference sequence.

Does a BAM contain unmapped reads? Usually yes, with FLAG bit 4 set and RNAME *. Extract them with samtools view -f 4, and check the count before assuming your BAM holds every read that came off the instrument.

What tools read BAM besides samtools? Most callers and viewers do, including desktop viewers that combine variant detection with alignment display ⁷. Pair-collation utilities such as biobambam handle marking and splitting without a full re-sort ¹.

Footnotes

1. German Tischler, Steven Leonard. biobambam: tools for read pair collation based algorithms on BAM files. Source Code for Biology and Medicine, 2014. https://doi.org/10.1186/1751-0473-9-13 ↩ ↩²

2. Jui-Hung Hung, Zhiping Weng. Mapping Short Sequence Reads to a Reference Genome. Cold Spring Harbor Protocols, 2016. https://doi.org/10.1101/pdb.prot093161 ↩

3. Chase A Miller, Yi Qiao, Tonya DiSera, et al. bam.iobio: a web-based, real-time, sequence alignment file inspector. Nature Methods, 2014. https://doi.org/10.1038/nmeth.3174 ↩

4. Jeremiah Wala, Cheng-Zhong Zhang, Matthew Meyerson, et al. VariantBam: filtering and profiling of next-generational sequencing data using region-specific rules. Bioinformatics, 2016. https://doi.org/10.1093/bioinformatics/btw111 ↩

5. Derek W. Barnett, Erik K. Garrison, Aaron R. Quinlan, et al. BamTools: a C++ API and toolkit for analyzing and managing BAM files. Bioinformatics, 2011. https://doi.org/10.1093/bioinformatics/btr174 ↩

6. Paul P.S. Wang, Wendy T. Parker, Susan Branford, et al. BAM-matcher: a tool for rapid NGS sample matching. Bioinformatics, 2016. https://doi.org/10.1093/bioinformatics/btw239 ↩

7. Michael N. Edmonson, Jinghui Zhang, Chunhua Yan, et al. Bambino: a variant detector and alignment viewer for next-generation sequencing data in the SAM/BAM format. Bioinformatics, 2011. https://doi.org/10.1093/bioinformatics/btr032 ↩