Metagenomics Fundamentals

Lesson 7 of 13 · 12 min

Shotgun Taxonomic Profiling with Kraken2 and Bracken

Kraken2 profiles a whole microbial community by breaking every shotgun read into k-mers and looking up, in a prebuilt database, the taxon each k-mer points to; that stored taxon is itself the lowest common ancestor of all genomes carrying the k-mer. Kraken2 then lays those hits on the taxonomy tree and assigns the read to the deepest taxon along the highest-scoring root-to-leaf path, so a read lands only as specifically as its k-mers agree. This classifies the whole community in one fast pass, unlike amplicon methods that target a single marker gene such as 16S. In this lesson you classify reads against the standard database, learn to read the report, and then correct the raw counts into species-level abundances with Bracken. With --paired each fragment is a read pair, and the large per-read --output file is often sent to /dev/null so only the compact --report is kept.

bash
kraken2 \
  --db /refs/k2_standard \
  --threads 8 \
  --confidence 0.1 \
  --paired \
  sample_R1.fastq.gz sample_R2.fastq.gz \
  --report sample.kreport \
  --output sample.kraken
Loading database information... done.
5000000 sequences (1500.00 Mbp) processed in 43.271s (6933.0 Kseq/m, 2079.87 Mbp/m).
  4185000 sequences classified (83.70%)
  815000 sequences unclassified (16.30%)

Reading the report

text
 16.30	815000	815000	U	0	unclassified
 83.70	4185000	1520	R	1	root
 81.20	4060000	2100	D	2	    Bacteria
 42.15	2107500	3400	P	1224	      Proteobacteria
 28.90	1445000	5200	C	1236	        Gammaproteobacteria
 23.00	1150000	6000	O	91347	          Enterobacterales
 18.40	920000	4300	F	543	            Enterobacteriaceae
 12.30	615000	9000	G	561	              Escherichia
 11.05	552500	540000	S	562	                Escherichia coli

The six tab-separated columns are: percent of all fragments under this clade; fragments assigned to the clade rooted here (this taxon plus its descendants); fragments assigned directly to this taxon only; the rank code (U unclassified, R root, D domain, K kingdom, P phylum, C class, O order, F family, G genus, S species, with a trailing number such as D1 for intermediate ranks); the NCBI taxonomy ID; and the indented scientific name. Look at the third column: the Escherichia genus keeps 9000 fragments of its own, reads whose k-mers matched the genus but could not be resolved to one species, so Kraken2 stopped there. Every clade above species holds such reads, which is exactly why raw species counts undercount true abundance and why Bracken is needed next.

The --confidence flag (0 to 1) is the minimum fraction of a read's k-mers that must fall in a taxon's clade, or the call is pushed toward the root or left unclassified. The default 0 is permissive; raising it to about 0.05 to 0.1 trims spurious species at the cost of leaving more reads at higher ranks. False positives also come from the database: a species absent from it has its reads misassigned to the nearest relative, and host or contaminant reads share k-mers with real taxa. MetaPhlAn takes a different route, mapping reads to clade-specific marker genes instead of all k-mers, so it reports fewer spurious species but can miss low-abundance or novel ones.

Correct abundances with Bracken

Kraken2 read counts are not species abundances. Because classification uses the lowest common ancestor, reads whose k-mers are shared across several species stop at the genus or family node, so species counts are under-estimates and do not reflect true proportions. Bracken re-estimates by probabilistically redistributing those higher-rank reads down to species using a per-database k-mer distribution built once with bracken-build, giving both corrected read counts and a relative-abundance fraction per species.

bash
bracken \
  -d /refs/k2_standard \
  -i sample.kreport \
  -o sample.bracken \
  -w sample_bracken.kreport \
  -r 150 \
  -l S \
  -t 10
name	taxonomy_id	taxonomy_lvl	kraken_assigned_reads	added_reads	new_est_reads	fraction_total_reads
Escherichia coli	562	S	552500	318200	870700	0.42180
Klebsiella pneumoniae	573	S	96000	41500	137500	0.06661
Bacteroides uniformis	820	S	78000	22300	100300	0.04859
Faecalibacterium prausnitzii	853	S	61000	18700	79700	0.03861

Try it yourself: Re-run Kraken2 on the same reads with --confidence 0 and again with --confidence 0.15, then diff the species (S) rows of the two reports. Count how many reads shift up to genus or family and whether any low-abundance species vanish. Run Bracken (-r matched to your read length) on each report and watch how fraction_total_reads for Escherichia coli changes between the two settings.