BioUnfold.in
All protocols

Align short reads to a reference with BWA-MEM

Index a reference genome, align paired-end short reads with BWA-MEM, and produce a coordinate-sorted, indexed BAM in one pipe.

SSSudipta SardarUpdated July 21, 2026

Prerequisites

  • bwa installed (`bwa` on your PATH)
  • samtools installed (`samtools` on your PATH)
  • A reference genome in FASTA format
  • Paired-end reads as two gzipped FASTQ files

1Index the reference

BWA needs a set of index files before it can align anything. Run this once per reference genome; you can skip it on later alignments against the same FASTA.

bash
bwa index reference.fasta

Expected output

[bwa_index] Pack FASTA... 0.15 sec
[bwa_index] Construct BWT for the packed sequence...
[bwa_index] 33.42 seconds elapse.
[bwa_index] Update BWT... 0.11 sec
[bwa_index] Pack forward-only FASTA... 0.10 sec
[bwa_index] Construct SA from BWT and Occ... 8.31 sec
[main] Version: 0.7.17-r1188
[main] CMD: bwa index reference.fasta
[main] Real time: 42.156 sec; CPU: 41.980 sec

2Align and sort in one pipe

Stream BWA-MEM's SAM output straight into samtools sort so you never write an intermediate SAM to disk. The `-t 4` and `-@ 4` flags give each tool four threads; the trailing `-` tells samtools to read the alignments from stdin.

bash
bwa mem -t 4 reference.fasta R1.fastq.gz R2.fastq.gz | samtools sort -@ 4 -o aligned.sorted.bam -

Expected output

[M::bwa_idx_load_from_disk] read 0 ALT contigs
[M::process] read 200000 sequences (20200000 bp)...
[M::mem_pestat] # candidate unique pairs for (FF, FR, RF, RR): (5, 89234, 8, 6)
[M::mem_pestat] analyzing insert size distribution for orientation FR...
[M::mem_pestat] (25, 50, 75) percentile: (232, 289, 358)
[M::mem_pestat] mean and std.dev: (297.44, 89.21)
[M::mem_process_seqs] Processed 200000 reads in 15.234 CPU sec, 4.021 real sec
[main] Version: 0.7.17-r1188
[main] CMD: bwa mem -t 4 reference.fasta R1.fastq.gz R2.fastq.gz
[main] Real time: 62.115 sec; CPU: 230.472 sec

3Index the sorted BAM

Create a `.bai` index so downstream tools and genome browsers can jump to any region without scanning the whole file. This prints nothing on success and produces aligned.sorted.bam.bai next to the BAM.

bash
samtools index aligned.sorted.bam

Expected output

(no output; creates aligned.sorted.bam.bai)

Troubleshooting

[E::bwa_idx_load_from_disk] fail to locate the index files

The reference has not been indexed yet. Run `bwa index reference.fasta` first, then re-run the alignment against the same FASTA path.

Downstream GATK or Picard tools complain about missing read groups

Attach a read group during alignment by adding `-R '@RG\tID:sample1\tSM:sample1\tPL:ILLUMINA'` to the bwa mem command, e.g. `bwa mem -t 4 -R '@RG\tID:sample1\tSM:sample1\tPL:ILLUMINA' reference.fasta R1.fastq.gz R2.fastq.gz | samtools sort -@ 4 -o aligned.sorted.bam -`.

References