How to Create CSI Index BAM File

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

You have a BAM alignment file and you need region queries to return in milliseconds rather than minutes. This guide walks through creating a CSI index (Coordinate Sorted Index) with samtools, checking that the index is valid, and fixing the errors that stop indexing. By the end you will have a sorted BAM, a .bam.csi file sitting beside it, and a set of commands you can paste into a production pipeline.

You need three things before you start. A coordinate-sorted BAM. A samtools 1.x installation, or htslib and pysam if you prefer Python. Enough free disk space for the index, which is a small fraction of the BAM size.

What a CSI Index Is and Why BAM Files Need One

A BAM file is a binary, block-compressed representation of aligned sequencing reads. It is written in BGZF blocks, a variant of gzip that lets a reader jump to any block without decompressing everything before it. That jump capability only helps if the reader knows which block holds data for a given genomic interval. The index supplies that map.

A CSI index stores a hierarchy of bins plus a linear offset table for each reference sequence listed in the BAM header. For any query interval, HTSlib computes the set of bins that can contain overlapping records, converts those bins into virtual file offsets, and seeks directly to the relevant BGZF blocks. The principle is the same one used by other compressed binary genomic formats, where only the data needed for the current view is transmitted rather than the entire file [1]. The index is a companion file, and without it a region query degrades into a full scan [2].

The same design pattern appears across modern bioinformatics formats. Graph annotation files are sorted, compressed with BGZIP, and indexed with HTSlib and tabix so that a subgraph query does not read the whole file [3]. Tabix itself is the standard query tool for bgzipped, indexed text files [4]. Quality-score compressors add index structures specifically to enable random access decompression [5]. K-mer dictionaries are indexed to make lookups fast [6], and read compressors build indexes so any read can be retrieved without full decompression [7]. The CSI format is the BAM-specific member of that family.

The two numbers inside a CSI index

Every CSI index carries two parameters written into its header. The minimum shift (min_shift) sets the smallest bin size as 2 raised to that power. The level count sets how many binning levels sit above it. Together they determine the largest coordinate the index can represent and the granularity of the bins.

With the samtools default of min_shift = 14, the smallest bin covers 16,384 bp and the coordinate ceiling is 2^29, or roughly 512 Mbp. Raising the shift to 20 gives 1,048,576 bp bins and a ceiling of 2^35, roughly 34 Gbp. That single parameter is the main lever you have when working with unusually long contigs.

What the index does not do

A CSI index does not compress the BAM further, does not change the alignments, and does not validate them. It adds a lookup structure. Tools that read BAM will behave identically on the underlying records whether the index is present or absent. What changes is speed for interval queries and the ability of tools to call fetch() at all.

CSI vs BAI: Choosing the Right Index Format

Two index formats exist for BAM. BAI is the original, and CSI is the newer and more flexible one.

For a human GRCh38 BAM, both formats work. The longest human chromosome, chr1, is about 249 Mbp, comfortably under the BAI ceiling. For a genome with a chromosome longer than 512 Mbp, or a concatenated reference with very long scaffolds, BAI will refuse to index past the limit and CSI is the correct choice.

A common misconception is that samtools index produces a CSI by default. It does not. Plain samtools index sample.bam writes sample.bam.bai. You must pass -c or -csi to get CSI output. The command is silent on success, which is why many users never notice which format they produced. Check the file extension before moving on.

Prerequisites: What You Need Before Indexing

Install samtools

Any of the usual routes work. A conda environment is the most reproducible:

conda create -n hts -c bioconda -c conda-forge samtools
conda activate hts
samtools --version

Expected output begins with a version line:

samtools 1.19
Using htslib 1.19

On Debian or Ubuntu, apt install samtools also works, though the packaged version can lag behind bioconda by several releases. If your pipeline depends on -write-index, check the version first.

Understand the file formats involved

Bioinformatics file formats fall into three useful groups: sequence formats, alignment formats, and variant or annotation formats. CSI indexes apply to the alignment and annotation groups. Knowing which group you are holding prevents most indexing mistakes.

FASTQ vs FASTA

FASTA and FASTQ both store nucleotide sequences as plain text. They differ in what else they carry.

FASTA has two parts per record: a header line beginning with > and one or more sequence lines. There are no quality values. It is the standard format for reference genomes, transcript sets, and assembled contigs.

>chr1
NNNNNNNNNNNNNNNNNNNN
ACGTACGTACGTACGTACGT

FASTQ has four lines per record: a header beginning with @, the sequence, a separator line beginning with +, and a quality string of the same length as the sequence. Each quality character encodes a Phred score for the corresponding base. FASTQ is the raw output of Illumina, PacBio, and Oxford Nanopore instruments, and it is the input to aligners.

@A00123:45:H7J2LDSX2:1:1101:1000:2000 1:N:0:ATCG
ACGTACGTACGTACGTACGT
+
FFFFFFFFFFFFFFFFFFFF

The quality line is the hardest part of FASTQ to compress losslessly, which is why specialized quality-score compressors exist [5] and why k-mer dictionaries have been used to improve quality-value compression [6]. When FASTQ is aligned, those quality strings move into the BAM QUAL field, so the same compression questions follow the data downstream.

A practical rule: if a file has quality scores, it is FASTQ and it came from a sequencer. If it has no quality scores, it is FASTA and it is a reference or an assembly.

SAM, BAM, and CRAM

SAM is the plain-text alignment format. It has a header of @HD, @SQ, @RG, @PG, and @CO lines followed by one tab-delimited line per alignment. BAM is the binary, BGZF-compressed equivalent of SAM. CRAM is a reference-based compressed alignment format that stores only differences from a reference.

CSI and BAI indexes work on BAM. CRAM uses its own index format, CRAI. You cannot attach a .csi to a CRAM and expect samtools to use it.

BGZF and virtual offsets

BGZF is gzip with extra metadata that records the compressed size of each block. A reader can therefore compute the byte offset of any block without decompressing the preceding ones. A virtual offset combines a block offset with an offset inside the decompressed block, packed into a 64-bit integer. CSI stores 64-bit virtual offsets, while BAI stores 32-bit ones, which is one reason CSI scales better.

Coordinate sorting is mandatory

Indexing only works on a BAM whose records are ordered by reference sequence and then by position. A queryname-sorted BAM, the kind produced by samtools sort -n, cannot be indexed with BAI or CSI. The header line records the sort order:

samtools view -H sample.bam | grep '^@HD'

A well-formed coordinate-sorted BAM shows:

@HD	VN:1.6	SO:coordinate

If SO:queryname or SO:unsorted appears, sort the file before indexing.

samtools sort -@ 8 -m 2G -o sample.sorted.bam sample.bam

The -m 2G flag sets the per-thread memory budget. The -@ 8 flag uses eight threads. The output header will carry SO:coordinate.

How to Create a CSI Index BAM File Step by Step

The sequence below works for a single BAM and scales to batches with a shell loop or a workflow manager. Each step prints or produces something you can check.

Step 1: Confirm the input is BAM and not SAM

samtools quickcheck -v sample.bam

Silence means the file passes a basic structural check. Any printed line indicates a truncated or malformed file. If you have a .sam text file, convert it first:

samtools view -b -o sample.bam sample.sam

SAM text cannot be indexed directly because there are no BGZF blocks to seek into.

Step 2: Confirm coordinate sort order

samtools view -H sample.sorted.bam | grep '^@HD'

Expected:

@HD	VN:1.6	SO:coordinate

If sorting is missing, return to the samtools sort command above. Sorting a 30x human genome BAM usually takes tens of minutes on eight threads and needs temporary disk space, so plan the location with -T.

Step 3: Confirm the header contains @SQ lines

samtools view -H sample.sorted.bam | grep '^@SQ' | head -3

Expected:

@SQ	SN:chr1	LN:248956422
@SQ	SN:chr2	LN:242193529
@SQ	SN:chr3	LN:198295559

A BAM with no @SQ lines, typically an unmapped-reads-only file, cannot be indexed. Index the mapped alignment, not the unmapped BAM.

Step 4: Create the CSI index

samtools index -c -@ 8 sample.sorted.bam

On success the command prints nothing and writes sample.sorted.bam.csi next to the BAM. To control the name or location:

samtools index -c -@ 8 -o /path/to/index/sample.sorted.bam.csi sample.sorted.bam

Step 5: Confirm the file exists and is non-empty

ls -lh sample.sorted.bam sample.sorted.bam.csi

Expected:

-rw-r--r-- 1 user group  12G Nov 10 14:02 sample.sorted.bam
-rw-r--r-- 1 user group  1.4M Nov 10 14:02 sample.sorted.bam.csi

Index size scales with the number of distinct bins that contain reads, not with the number of reads. A whole-genome BAM typically gets an index in the low single-digit megabytes.

Step 6: Verify with idxstats

samtools idxstats sample.sorted.bam | head -5

Expected output has four columns per reference: name, length, mapped read count, unmapped read count.

chr1	248956422	2013456	0
chr2	242193529	1954321	0
chr3	198295559	1743221	0
chr4	190214555	1601223	0
chr5	181538259	1550122	0

If idxstats fails with a message about a missing or unreadable index, the index was not written to the location samtools expects. The default lookup path is <bam>.csi beside the BAM.

Step 7: Run a timed region query

time samtools view -c sample.sorted.bam chr1:1000000-1100000

The answer is a read count, and the elapsed time should be well under a second on a local disk. If the query takes as long as reading the whole file, the index was not picked up. Run the same command under strace or add -verbosity 4 to see which index file samtools opened.

Step 8: Cross-check against a known region

samtools view sample.sorted.bam chr7:55019017-55211628 | wc -l

Pick a region you have counted before from an independent run. A mismatch means the BAM was modified after the index was built. The index is a snapshot of record positions, so any re-sort, filter, or reheader invalidates it.

Key Options and Parameters

The samtools index subcommand exposes a small set of options. These are the ones that matter in practice.

The -m flag only has meaning with CSI. Passing -m also turns on CSI output even if -c is omitted in many builds, but relying on that side effect is fragile. Pass both explicitly.

Understanding min-shift in practice

The shift value sets the smallest bin size. A shift of 14 gives 16 kb bins, which is fine for short-read alignments where read density is roughly uniform. A shift of 20 gives 1 Mb bins and a much larger coordinate ceiling, at the cost of coarser bin granularity and slightly more records examined per query.

For a reference with one or more sequences longer than 512 Mbp, use:

samtools index -c -m 20 -@ 8 sample.sorted.bam

For a reference where the longest sequence is 400 Mbp, the default -m 14 is correct and produces smaller indexes.

Writing the index during sorting

Recent samtools releases can write the index in the same pass as sorting, which avoids a second read of a large BAM.

samtools sort -@ 8 -m 2G --write-index=csi -o sample.sorted.bam sample.bam

The -write-index option accepts csi or bai. Check samtools sort, help on your build to confirm the syntax, because older releases lack the option entirely.

How to Verify the Index Works

Four checks catch nearly every indexing problem. Run them once after building, then again after any pipeline change.

Check 1: File pairing

The index must sit beside the BAM or be passed explicitly. Tools do not consult a central registry.

ls sample.sorted.bam*

Expected:

sample.sorted.bam
sample.sorted.bam.csi

If both a .bai and a .csi exist, remove the one you do not want. Two indexes for one BAM invite ambiguous behavior across tool versions.

Check 2: Metadata retrieval

samtools view -H sample.sorted.bam | wc -l

This does not exercise the index, but a failure here means the BAM itself is broken and the index is irrelevant.

Check 3: Per-reference statistics

samtools idxstats sample.sorted.bam > sample.idxstats.tsv
awk '$3 > 0' sample.idxstats.tsv | wc -l

The count is the number of references with at least one mapped read. A value of zero for a chromosome you expect to be covered points to an alignment or reference mismatch, not an index bug.

Check 4: Query round trip

samtools view sample.sorted.bam chr1:1-100000 | samtools view -c -

The two-step form confirms that the fetched records can be piped into the next stage. If the first command errors with a message about an invalid region or unknown reference, the chromosome name in the query does not match the @SQ names in the header. Check the exact spelling, including whether the reference uses chr1 or 1.

Worked Example: From FASTQ to Indexed BAM

A complete short-read path shows where the index fits in a real workflow.

# 1. Index the reference FASTA
bwa index reference.fa

# 2. Align paired FASTQ reads and sort in one pipe
bwa mem -t 8 -R '@RG\tID:sample1\tSM:sample1\tPL:ILLUMINA' \
    reference.fa sample_R1.fastq.gz sample_R2.fastq.gz \
  | samtools sort -@ 8 -m 2G -o sample1.sorted.bam -

# 3. Mark duplicates before final indexing
samtools markdup -@ 8 sample1.sorted.bam sample1.md.bam

# 4. Build the CSI index
samtools index -c -@ 8 sample1.md.bam

# 5. Confirm
samtools idxstats sample1.md.bam | head

Two details matter here. samtools markdup preserves coordinate order, so the index can be built on its output without re-sorting. If you use a tool that changes record order, re-sort before indexing.

The same flow in Python with pysam:

import pysam

pysam.index("-c", "-@", "8", "sample1.md.bam")

bam = pysam.AlignmentFile("sample1.md.bam", "rb")
for read in bam.fetch("chr1", 1000000, 1010000):
    print(read.query_name)
bam.close()

The fetch() call only works because the CSI index resolves the region to block offsets. Without it, pysam raises an error rather than scanning.

At the C level, htslib exposes sam_index_build3(fn, fnidx, min_shift, nthreads). Passing a positive min_shift selects CSI, and passing zero selects BAI. That single integer is the difference between the two formats in most codebases.

Indexing VCF, BED, and Other BGZF Files with CSI

The same CSI concepts apply to bgzipped text files, which is where the tabix and bcftools commands come in.

VCF

bgzip -c variants.vcf > variants.vcf.gz
bcftools index --csi variants.vcf.gz

This writes variants.vcf.gz.csi. The tabix route is equivalent:

tabix -C -p vcf variants.vcf.gz

The -C flag selects CSI. The -p vcf preset chooses the column layout for a VCF. Without -C, tabix writes a .tbi index, which has the same coordinate ceiling as BAI.

BED

bgzip -c regions.bed > regions.bed.gz
tabix -C -p bed regions.bed.gz

Why plain gzip fails here

Tabix and CSI require BGZF, not ordinary gzip. A plain .gz file has no block structure, so there is nothing to seek into. If you see an error about a missing BGZF header, recompress with bgzip rather than gzip.

Querying an indexed VCF by region

tabix variants.vcf.gz chr1:1000000-1100000

The output is the matching VCF lines, produced by seeking to the relevant compressed blocks rather than decompressing the whole file. The design mirrors the BAM case, and index-based partial retrieval is what makes hosted, no-download access to large variant files practical [2].

Common Mistakes and Limitations

These are the failures that account for most broken indexes.

Indexing an unsorted BAM

The most common error. samtools reports a message like this:

[E::bam_index_core] the alignment is not sorted (chr2:... vs chr1:...)
samtools index: failed to create index for sample.bam

The fix is to sort first. There is no flag that makes indexing work on an unsorted file. Sorting is a hard requirement because the index assumes records for a reference are contiguous and increasing in position.

Expecting samtools index to make a CSI

Plain samtools index sample.bam creates a BAI. On a genome where BAI is sufficient, nobody notices. On a genome with a 700 Mbp chromosome, the failure appears later during a region query rather than at index time, which makes it easy to misdiagnose.

Hitting the coordinate ceiling

When a contig exceeds what the current parameters allow, indexing fails with a message along these lines:

[E::hts_idx_push] Region 0:... cannot be stored in index
samtools index: failed to create index

The fix is -c -m 20 rather than plain -c. Confirm the longest contig from the @SQ lines before choosing the shift.

Stale index after modifying the BAM

Any operation that rewrites the BAM invalidates the index. Filtering with samtools view -b, reheadering with samtools reheader, or re-sorting all change record positions. Some tools emit a warning:

[W::hts_idx_load3] The index file is older than the data file: sample.bam.csi

Others fail silently and return wrong results. Rebuild the index whenever the BAM changes.

Copying the BAM without the index

An index is not embedded in the BAM. If you move or copy only the .bam to another directory, the index does not follow. Copy with a glob that includes both files, or rebuild on the destination.

Reference name mismatches

samtools view sample.bam chr1:1-1000 fails when the BAM header uses 1 instead of chr1. The error is usually:

[main_samview] region "chr1:1-1000" specifies an invalid region or unknown reference

Read the @SQ names from the header and match their exact spelling. This bites hardest when merging BAMs produced against different reference builds.

Indexing a SAM text file

samtools index sample.sam

Produces an error stating the file cannot be usefully indexed. Convert to BAM first with samtools view -b.

Using .tbi where CSI is expected

If a VCF was indexed with default tabix settings, the file on disk is variants.vcf.gz.tbi, not .csi. Tools that specifically look for a CSI file will not find one. Use tabix -C or bcftools index, csi when the downstream tool requires CSI.

Confusing CSI with CRAI

CRAM files use .crai indexes. A .csi file next to a CRAM is ignored. The two formats are not interchangeable.

Limitations you should plan around

An index speeds up region queries, but it cannot fix a poor alignment. Whole-genome scans and variant callers that walk every position do not benefit from it. A CSI index also does not reduce BAM size in a meaningful way, and it does not travel inside the BAM, so every destination needs its own index file. Indexing itself streams the entire BAM once, which means building an index on a 100 GB file still reads 100 GB from disk.

Tips for Production Pipelines

A few habits keep indexing from becoming the step that silently breaks.

Validate the sort order in the pipeline rather than on a laptop. Checking @HD SO:coordinate costs one pass over the header and catches upstream changes.

Store the index path explicitly in tools that accept one. Many aligners and viewers accept an index argument, and passing it avoids ambiguity when several indexes coexist.

Rebuild indexes as part of any transformation step, not as an afterthought. If your pipeline is a sequence of samtools view -b, samtools sort, and samtools markdup calls, put samtools index -c at the end of the chain and treat its output as a required artifact.

Use checksums on the pair. Compute md5sum sample.bam sample.bam.csi and store the values. A mismatched index is hard to spot by eye and easy to spot by hash.

Prefer CSI for new projects. It costs nothing extra and avoids a painful migration when a reference update introduces a long contig.

Keep one index format per BAM. Delete the .bai when you switch to .csi so that different tools cannot pick different files.

Frequently Asked Questions

What is a CSI index BAM file?

A CSI index BAM file is a pair: the BAM alignment file and a .bam.csi file that maps genomic coordinates to BGZF block offsets. The CSI file is not a BAM. It is a separate lookup structure that lets tools jump to the right part of the BAM for a given region.

Does samtools index create a CSI by default?

No. Plain samtools index sample.bam writes sample.bam.bai. You must add -c or -csi to get a .bam.csi file.

Can I create a CSI index without sorting the BAM?

No. Both BAI and CSI require coordinate-sorted input. Sorting is the step that makes record positions monotonic, which is the assumption the index is built on.

What is the difference between a .bai and a .csi file?

Both let tools jump to a genomic region. BAI uses 32-bit offsets and a fixed binning scheme capped at about 512 Mbp. CSI uses 64-bit offsets and a configurable binning scheme, so it can handle longer contigs when you raise the minimum shift.

How large can a contig be with a CSI index?

With the default min_shift of 14, the ceiling is 2^29 bp, about 512 Mbp. Raising the shift to 20 moves the ceiling to 2^35 bp, about 34 Gbp, which covers any realistic chromosome.

Why do I get "could not retrieve index file"?

The tool looked for <bam>.csi or <bam>.bai beside the BAM and found neither, or found a file it could not parse. Check the file name, the directory, and whether the index was built before the BAM was last modified.

Can I index a CRAM file with CSI?

No. CRAM uses its own .crai index format. Run samtools index sample.cram and let samtools choose the correct format for the input.

Do I need to delete the BAI when I create a CSI?

You should. Keeping both files lets different tools open different indexes, which produces inconsistent query behavior. Pick one format per BAM and remove the other.

Sources

1. BigWig and BigBed: enabling browsing of large distributed datasets.
2. FAIR privacy-preserving operation of large genomic variant calling format (VCF) data without download or installation.
3. Efficient indexing and querying of annotations in a pangenome graph.
4. STABIX: summary-statistic-based GWAS indexing and compression.
5. LCQS: an efficient lossless compression tool of quality scores with random access functionality.
6. Indexing k-mers in linear space for quality value compression.
7. SCALCE: boosting sequence compression algorithms using locally consistent encoding.