
Bio Workflows Metagenomics Pipeline
- 4 installs
- 1.1k repo stars
- Updated July 25, 2026
- gptomics/bioskills
Run an end-to-end metagenomics workflow from FASTQ to taxonomic and functional profiles using Kraken2, Bracken, and HUMAnN.
About
Orchestrates host depletion, Kraken2/Bracken taxonomic classification, MetaPhlAn profiling, and HUMAnN3 functional analysis. A developer uses it to profile metagenomic samples from raw reads to taxonomic and functional output.
- Kraken2/Bracken classification and MetaPhlAn profiling
- HUMAnN3 functional profiling with classification-rate QC gates
Bio Workflows Metagenomics Pipeline by the numbers
- 4 all-time installs (skills.sh)
- Ranked #1,625 of 2,064 Data Science & ML skills by installs in the Skillselion catalog
- Data as of Aug 5, 2026 (Skillselion catalog sync)
npx skills add https://github.com/gptomics/bioskills --skill bio-workflows-metagenomics-pipelineAdd your badge
Show developers this skill is listed on Skillselion. Paste this into your README.
| Installs | 4 |
|---|---|
| repo stars | ★ 1.1k |
| Last updated | July 25, 2026 |
| Repository | gptomics/bioskills ↗ |
What it does
Run an end-to-end metagenomics workflow from FASTQ to taxonomic and functional profiles using Kraken2, Bracken, and HUMAnN.
Files
Version Compatibility
Reference examples tested with: Bowtie2 2.5.3+, Bracken 2.9+, HUMAnN 3.8+, Kraken2 2.1+, MetaPhlAn 4.1+, fastp 0.23+, matplotlib 3.8+, pandas 2.2+, seaborn 0.13+
Before using code patterns, verify installed versions match. If versions differ:
- Python:
pip show <package>thenhelp(module.function)to check signatures - CLI:
<tool> --versionthen<tool> --helpto confirm flags
If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.
Metagenomics Pipeline
"Analyze my metagenomic samples from FASTQ to taxonomic and functional profiles" -> Orchestrate host depletion, Kraken2/Bracken taxonomic classification, MetaPhlAn profiling, HUMAnN3 functional analysis, and AMR gene detection.
Complete workflow from metagenomic FASTQ to taxonomic and functional profiles.
Workflow Overview
FASTQ files
|
v
[1. QC & Host Removal] --> fastp + Bowtie2
|
v
[2. Taxonomic Classification]
|
+---> Kraken2 + Bracken (fast, database-dependent)
|
+---> MetaPhlAn (marker-based, standardized)
|
v
[3. Functional Profiling] --> HUMAnN
|
v
Taxonomic profiles + Pathway abundancesPrimary Path: Kraken2 + Bracken + HUMAnN
Step 1: Quality Control and Host Removal
# QC with fastp
for sample in sample1 sample2 sample3; do
fastp -i ${sample}_R1.fastq.gz -I ${sample}_R2.fastq.gz \
-o trimmed/${sample}_R1.fq.gz -O trimmed/${sample}_R2.fq.gz \
--detect_adapter_for_pe \
--qualified_quality_phred 20 \
--length_required 50 \
--html qc/${sample}_fastp.html
done
# Remove host reads (human example)
for sample in sample1 sample2 sample3; do
bowtie2 -p 8 -x human_index \
-1 trimmed/${sample}_R1.fq.gz \
-2 trimmed/${sample}_R2.fq.gz \
--un-conc-gz host_removed/${sample}_R%.fq.gz \
> /dev/null 2> qc/${sample}_host_removal.log
doneStep 2A: Kraken2 Classification
# Classify reads
for sample in sample1 sample2 sample3; do
kraken2 --db kraken2_db \
--threads 8 \
--paired \
--report kraken/${sample}.report \
--output kraken/${sample}.output \
host_removed/${sample}_R1.fq.gz \
host_removed/${sample}_R2.fq.gz
doneStep 2B: Bracken Abundance Estimation
# Estimate species abundance
for sample in sample1 sample2 sample3; do
bracken -d kraken2_db \
-i kraken/${sample}.report \
-o bracken/${sample}.species.txt \
-r 150 \
-l S \
-t 10
done
# Combine samples into abundance matrix
combine_bracken_outputs.py \
--files bracken/*.species.txt \
-o bracken/combined_species.txtStep 2C: Alternative - MetaPhlAn Profiling
# Profile with MetaPhlAn 4
for sample in sample1 sample2 sample3; do
metaphlan host_removed/${sample}_R1.fq.gz,host_removed/${sample}_R2.fq.gz \
--bowtie2out metaphlan/${sample}.bowtie2.bz2 \
--input_type fastq \
--nproc 8 \
-o metaphlan/${sample}_profile.txt
done
# Merge profiles
merge_metaphlan_tables.py metaphlan/*_profile.txt > metaphlan/merged_abundance.txtStep 3: Functional Profiling with HUMAnN
# Run HUMAnN
for sample in sample1 sample2 sample3; do
# Concatenate paired reads
cat host_removed/${sample}_R1.fq.gz host_removed/${sample}_R2.fq.gz > \
host_removed/${sample}_concat.fq.gz
humann --input host_removed/${sample}_concat.fq.gz \
--output humann/${sample} \
--threads 8 \
--metaphlan-options "--bowtie2db metaphlan_db"
done
# Normalize and join tables
humann_renorm_table --input humann/sample1/sample1_pathabundance.tsv \
--output humann/sample1/sample1_pathabundance_cpm.tsv \
--units cpm
humann_join_tables --input humann \
--output humann/merged_pathabundance.tsv \
--file_name pathabundanceVisualization
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
# Load Bracken species table
species = pd.read_csv('bracken/combined_species.txt', sep='\t', index_col=0)
# Top 20 species heatmap
top20 = species.sum(axis=1).nlargest(20).index
plt.figure(figsize=(12, 8))
sns.heatmap(species.loc[top20], cmap='viridis', annot=False)
plt.title('Top 20 Species Abundance')
plt.tight_layout()
plt.savefig('top20_species_heatmap.pdf')
# Stacked bar plot
species_norm = species.div(species.sum()) * 100
top10 = species_norm.sum(axis=1).nlargest(10).index
other = species_norm.loc[~species_norm.index.isin(top10)].sum()
plot_data = species_norm.loc[top10].T
plot_data['Other'] = other
plot_data.plot(kind='bar', stacked=True, figsize=(10, 6))
plt.ylabel('Relative Abundance (%)')
plt.legend(bbox_to_anchor=(1.05, 1))
plt.tight_layout()
plt.savefig('species_barplot.pdf')Parameter Recommendations
| Step | Parameter | Value |
|---|---|---|
| fastp | --length_required | 50 (metagenomic reads) |
| Kraken2 | --confidence | 0.0 (default) or 0.1 |
| Bracken | -r | Read length (e.g., 150) |
| Bracken | -l | S (species) or G (genus) |
| Bracken | -t | 10 (min reads threshold) |
| MetaPhlAn | --min_cu_len | 2000 (default) |
| HUMAnN | --threads | 8+ |
Troubleshooting
| Issue | Likely Cause | Solution |
|---|---|---|
| Low classification rate | Database mismatch, novel organisms | Try different database, check sample type |
| High unclassified | Novel microbes, host contamination | Remove host, use larger database |
| High host reads | Incomplete host removal | Use multiple host reference genomes |
| HUMAnN slow | Large files | Increase threads, pre-filter reads |
Complete Pipeline Script
#!/bin/bash
set -e
THREADS=8
KRAKEN_DB="kraken2_standard_db"
HOST_INDEX="human_bt2_index"
SAMPLES="sample1 sample2 sample3"
OUTDIR="metagenomics_results"
mkdir -p ${OUTDIR}/{trimmed,host_removed,kraken,bracken,metaphlan,humann,qc}
# Step 1: QC
echo "=== QC ==="
for sample in $SAMPLES; do
fastp -i ${sample}_R1.fastq.gz -I ${sample}_R2.fastq.gz \
-o ${OUTDIR}/trimmed/${sample}_R1.fq.gz \
-O ${OUTDIR}/trimmed/${sample}_R2.fq.gz \
--length_required 50 \
--html ${OUTDIR}/qc/${sample}_fastp.html -w ${THREADS}
done
# Host removal
echo "=== Host Removal ==="
for sample in $SAMPLES; do
bowtie2 -p ${THREADS} -x ${HOST_INDEX} \
-1 ${OUTDIR}/trimmed/${sample}_R1.fq.gz \
-2 ${OUTDIR}/trimmed/${sample}_R2.fq.gz \
--un-conc-gz ${OUTDIR}/host_removed/${sample}_R%.fq.gz \
> /dev/null 2> ${OUTDIR}/qc/${sample}_host.log
done
# Step 2: Kraken2
echo "=== Kraken2 ==="
for sample in $SAMPLES; do
kraken2 --db ${KRAKEN_DB} --threads ${THREADS} --paired \
--report ${OUTDIR}/kraken/${sample}.report \
--output ${OUTDIR}/kraken/${sample}.output \
${OUTDIR}/host_removed/${sample}_R1.fq.gz \
${OUTDIR}/host_removed/${sample}_R2.fq.gz
done
# Bracken
echo "=== Bracken ==="
for sample in $SAMPLES; do
bracken -d ${KRAKEN_DB} \
-i ${OUTDIR}/kraken/${sample}.report \
-o ${OUTDIR}/bracken/${sample}.species.txt \
-r 150 -l S -t 10
done
echo "=== Pipeline Complete ==="
echo "Kraken reports: ${OUTDIR}/kraken/"
echo "Bracken abundances: ${OUTDIR}/bracken/"Related Skills
- database-access/sra-data - Pull metagenomic FASTQ from SRA / ENA (16S amplicon or shotgun)
- database-access/ncbi-datasets-cli - Bulk-pull reference genomes for read mapping
- database-access/remote-homology - DIAMOND --ultra-sensitive for predicted-ORF annotation
- metagenomics/kraken-classification - Kraken2 details
- metagenomics/metaphlan-profiling - MetaPhlAn parameters
- metagenomics/abundance-estimation - Bracken options
- metagenomics/functional-profiling - HUMAnN workflow
- metagenomics/metagenome-visualization - Plotting functions
#!/bin/bash
# Reference: Bowtie2 2.5.3+, Bracken 2.9+, HUMAnN 3.8+, Kraken2 2.1+, MetaPhlAn 4.1+, fastp 0.23+, matplotlib 3.8+, pandas 2.2+, seaborn 0.13+ | Verify API if version differs
# Complete metagenomics workflow: Kraken2 + Bracken + HUMAnN
set -e
THREADS=8
KRAKEN_DB="/path/to/kraken2_standard"
HOST_INDEX="/path/to/human_bt2"
SAMPLES="sample1 sample2 sample3"
OUTDIR="metagenomics_results"
mkdir -p ${OUTDIR}/{trimmed,host_removed,kraken,bracken,humann,qc,viz}
echo "=== Metagenomics Pipeline ==="
echo "Samples: ${SAMPLES}"
# === Step 1: Quality Control ===
echo "=== Step 1: QC ==="
for sample in $SAMPLES; do
echo "QC: ${sample}"
fastp \
-i ${sample}_R1.fastq.gz \
-I ${sample}_R2.fastq.gz \
-o ${OUTDIR}/trimmed/${sample}_R1.fq.gz \
-O ${OUTDIR}/trimmed/${sample}_R2.fq.gz \
--detect_adapter_for_pe \
--qualified_quality_phred 20 \
--length_required 50 \
--complexity_threshold 30 \
--html ${OUTDIR}/qc/${sample}_fastp.html \
-w ${THREADS}
done
# === Step 2: Host Removal ===
echo "=== Step 2: Host Removal ==="
for sample in $SAMPLES; do
echo "Removing host: ${sample}"
bowtie2 -p ${THREADS} -x ${HOST_INDEX} \
-1 ${OUTDIR}/trimmed/${sample}_R1.fq.gz \
-2 ${OUTDIR}/trimmed/${sample}_R2.fq.gz \
--un-conc-gz ${OUTDIR}/host_removed/${sample}_R%.fq.gz \
--very-sensitive \
> /dev/null 2> ${OUTDIR}/qc/${sample}_host_removal.log
# Stats
total=$(zcat ${OUTDIR}/trimmed/${sample}_R1.fq.gz | wc -l)
remaining=$(zcat ${OUTDIR}/host_removed/${sample}_R1.fq.gz | wc -l)
host_pct=$(echo "scale=2; 100 - ($remaining / $total * 100)" | bc)
echo "${sample}: ${host_pct}% host reads removed"
done
# === Step 3: Kraken2 Classification ===
echo "=== Step 3: Kraken2 Classification ==="
for sample in $SAMPLES; do
echo "Classifying: ${sample}"
kraken2 \
--db ${KRAKEN_DB} \
--threads ${THREADS} \
--paired \
--report ${OUTDIR}/kraken/${sample}.report \
--report-minimizer-data \
--output ${OUTDIR}/kraken/${sample}.output \
${OUTDIR}/host_removed/${sample}_R1.fq.gz \
${OUTDIR}/host_removed/${sample}_R2.fq.gz
# Classification rate
classified=$(grep -m1 "root" ${OUTDIR}/kraken/${sample}.report | awk '{print $1}')
echo "${sample}: ${classified}% classified"
done
# === Step 4: Bracken Abundance ===
echo "=== Step 4: Bracken Abundance ==="
for sample in $SAMPLES; do
echo "Estimating abundance: ${sample}"
# Species level
bracken -d ${KRAKEN_DB} \
-i ${OUTDIR}/kraken/${sample}.report \
-o ${OUTDIR}/bracken/${sample}.species.txt \
-w ${OUTDIR}/bracken/${sample}.species.report \
-r 150 -l S -t 10
# Genus level
bracken -d ${KRAKEN_DB} \
-i ${OUTDIR}/kraken/${sample}.report \
-o ${OUTDIR}/bracken/${sample}.genus.txt \
-w ${OUTDIR}/bracken/${sample}.genus.report \
-r 150 -l G -t 10
done
# Combine tables
echo "Combining abundance tables..."
combine_bracken_outputs.py \
--files ${OUTDIR}/bracken/*.species.txt \
-o ${OUTDIR}/bracken/combined_species.txt
# === Step 5: HUMAnN Functional Profiling (Optional) ===
echo "=== Step 5: HUMAnN Functional Profiling ==="
for sample in $SAMPLES; do
echo "Functional profiling: ${sample}"
# Concatenate paired reads
cat ${OUTDIR}/host_removed/${sample}_R1.fq.gz \
${OUTDIR}/host_removed/${sample}_R2.fq.gz > \
${OUTDIR}/host_removed/${sample}_concat.fq.gz
humann \
--input ${OUTDIR}/host_removed/${sample}_concat.fq.gz \
--output ${OUTDIR}/humann/${sample} \
--threads ${THREADS}
rm ${OUTDIR}/host_removed/${sample}_concat.fq.gz
done
# Join HUMAnN tables
humann_join_tables -i ${OUTDIR}/humann -o ${OUTDIR}/humann/merged_pathabundance.tsv \
--file_name pathabundance
humann_join_tables -i ${OUTDIR}/humann -o ${OUTDIR}/humann/merged_genefamilies.tsv \
--file_name genefamilies
# Normalize
humann_renorm_table -i ${OUTDIR}/humann/merged_pathabundance.tsv \
-o ${OUTDIR}/humann/merged_pathabundance_cpm.tsv -u cpm
echo "=== Pipeline Complete ==="
echo "Results:"
echo " Kraken2 reports: ${OUTDIR}/kraken/"
echo " Bracken abundances: ${OUTDIR}/bracken/"
echo " HUMAnN pathways: ${OUTDIR}/humann/"
Metagenomics Pipeline - Usage Guide
Overview
This workflow processes metagenomic sequencing data to produce taxonomic and functional profiles. It supports both Kraken2/Bracken and MetaPhlAn approaches.
Prerequisites
# CLI tools
conda install -c bioconda kraken2 bracken metaphlan humann bowtie2 fastp
# Download databases
kraken2-build --download-taxonomy --db kraken2_db
kraken2-build --download-library bacteria --db kraken2_db
kraken2-build --build --db kraken2_dbQuick Start
Tell your AI agent what you want to do:
- "Run the metagenomics pipeline on my shotgun metagenomic data"
- "Classify my metagenomic reads and get species abundances"
- "Profile pathways in my gut microbiome samples"
Example Prompts
Taxonomic profiling
"Classify my metagenomic reads with Kraken2"
"Use MetaPhlAn to profile my samples"
"Get species-level abundances with Bracken"
Functional profiling
"Run HUMAnN on my metagenome"
"Get pathway abundances for my samples"
Visualization
"Create a stacked bar plot of species abundances"
"Make a heatmap of the top taxa"
Input Requirements
| Input | Format | Description |
|---|---|---|
| FASTQ files | .fastq.gz | Paired-end metagenomic reads |
| Kraken2 database | Directory | Pre-built or custom database |
| Host reference | FASTA | For host read removal |
What the Workflow Does
1. Quality Control - Filter low-quality reads 2. Host Removal - Remove host contamination 3. Classification - Assign taxonomic labels 4. Abundance - Estimate species/genus abundances 5. Functional - Profile metabolic pathways
Kraken2/Bracken vs MetaPhlAn
| Feature | Kraken2/Bracken | MetaPhlAn |
|---|---|---|
| Speed | Very fast | Moderate |
| Database | Large, customizable | Marker genes |
| Accuracy | Database-dependent | Standardized |
| Best for | Large-scale screening | Publication-ready profiles |
Tips
- Host removal: Critical for human-associated samples
- Database choice: Standard for general use, custom for specific environments
- Depth: Aim for >5M reads per sample for reliable profiling
- Replicates: Important for differential abundance analysis