A Nextflow pipeline for detection and quantification of Clostridioides difficile PaLoc / CdtLoc toxin genes from paired-end shotgun metagenomes. It performs quality control, optional host decontamination, targeted mapping of the toxin-gene panel with a built-in decoy set, whole-genome C. difficile coverage estimation, optional reference-guided assembly for contig-level confirmation, and antimicrobial-resistance gene (ARG) profiling against the CARD database.
This repository includes a custom bioinformatic workflow for our lab (LPM IKEM, Prague, Czechia). All parameters are tuned to work on specific data from our lab. Before any usage, please check that it fits your data as well.
- 🧬 Overview
- 🚀 Quick Usage
- 📦 Requirements
- 🛠️ Installation
- 🔬 Pipeline Description
- 📥 Inputs
- 📤 Outputs
- 📚 References
The pipeline takes paired-end FASTQ files and a samplesheet as input and produces:
- Quality-controlled and (optionally) host-decontaminated reads
- KMA mapping of the PaLoc/CdtLoc toxin-gene panel with per-template identity, breadth, and depth
- A decoy-aware detection call per toxin gene (guarded against LCT-homolog cross-mapping)
- C. difficile genome median depth, used to report a toxigenic fraction only above a minimum coverage
- Merged long- and wide-format PaLoc result tables plus a detection summary
- Optional reference-guided assembly of C. difficile-mapped reads (metaSPAdes) with QUAST quality metrics
- A normalised antimicrobial-resistance gene (ARG) matrix from RGI + CARD
nextflow run main.nf \
--input samplesheet.csv \
--outdir results/User needs conda environment with Nextflow & Singularity installed. Everything else needed for the pipeline execution can be prepared with the setup_pipeline.sh script.
| Tool | Version (container) | Purpose |
|---|---|---|
| Nextflow | ≥ 22.10 | Workflow manager |
| Singularity | any | Container runtime |
| fastp | 0.23.4 / 1.3.3 | Read trimming, filtering, and re-pairing |
| Bowtie2 | 2.5.5 | Host + PhiX decontamination, C. diff coverage alignment |
| KMA | 1.6.11 | Mapping to the PaLoc/CdtLoc + decoy panel |
| SAMtools | 1.23.1 | Coverage, depth, and read extraction |
| BWA-MEM2 | 2.3 | Mapping reads to the C. diff reference (assembly arm) |
| SPAdes | 4.3.0 | Reference-guided assembly |
| QUAST | 5.3.0 | Assembly quality assessment |
| RGI | 6.0.8 | Antimicrobial-resistance gene detection (CARD) |
| pandas | 2.2.1 | Normalisation and merging scripts |
| MultiQC | 1.21 | QC report aggregation |
| Database | Provided via | Used by |
|---|---|---|
PaLoc/CdtLoc panel FASTA (tcd*/cdt* genes + decoy LCT homologs) |
--paloc_db |
KMA |
| C. difficile reference genome FASTA | --cdiff_ref |
Bowtie2 coverage + BWA-MEM2 assembly arm |
CARD JSON (card.json) |
--card_json |
RGI |
Human reference index (Bowtie2, human-t2t-hla-argos985-mycob140) |
--host_genome_index |
Bowtie2 host removal |
PhiX index (Bowtie2, phiX174) |
--bowtie_phix_index |
Bowtie2 PhiX removal |
1. Clone the repository
git clone https://github.com/xpolak37/CdiffScreen.git
cd CdiffScreen2. Install Nextflow & Singularity
conda install -c bioconda nextflow
conda install -c conda-forge singularity3. Point the pipeline at your resources
Edit the paths at the top of nextflow.config (or override them on the command line) so they match your environment:
singularity_cache_dir = '/path/to/singularity_cache'
host_genome_index = '/path/to/bowtie2_human_index'
bowtie_phix_index = '/path/to/bowtie2_phix_index'
paloc_db = '/path/to/paloc_cdtloc_plus_decoys.fasta'
cdiff_ref = '/path/to/cdiff_ref.fasta'
card_json = '/path/to/card.json'flowchart TB
A["Input samplesheet"] --> B["fastp<br/>quality trimming"]
B --> HD{"--host_decontamination?"}
HD -- yes --> H1["Host removal<br/>Bowtie2 (human)"]
H1 --> H2["PhiX removal<br/>Bowtie2 (phiX174)"]
H2 --> H3["Re-pair / sync<br/>fastp"]
HD -- no --> CLEAN
H3 --> CLEAN["Clean reads"]
CLEAN --> K["KMA map<br/>PaLoc + decoy panel"]
CLEAN --> C1["Bowtie2 align<br/>C. diff genome"]
C1 --> C2["Coverage / median depth<br/>SAMtools"]
K --> N["Normalise & report"]
C2 --> N
N --> R["PaLoc result tables<br/>+ detection summary"]
CLEAN --> AS{"--run_assembly?"}
AS -- yes --> M1["Map to C. diff<br/>BWA-MEM2"]
M1 --> M2["Extract mapped pairs<br/>SAMtools"]
M2 --> M3["Reference-guided assembly<br/>metaSPAdes"]
M3 --> M4["Assembly QC<br/>QUAST"]
CLEAN --> G1["RGI bwt<br/>CARD"]
G1 --> G2["Merge & normalise<br/>ARG matrix"]
Quality control (fastp) Raw reads are trimmed and filtered with fastp. Adapter detection, quality/length filtering, polyG/X trimming, deduplication, low-complexity filtering, and overlap correction are all fully configurable (see Inputs).
Host decontamination — optional (--host_decontamination, default true)
Reads mapping to the human reference (T2T-HLA + ARGOS + mycobacterial sequences) are removed with Bowtie2, followed by a second pass to remove PhiX174 spike-in reads. Surviving reads are then re-paired and synchronised with fastp so that read1/read2 stay in register.
PaLoc / CdtLoc mapping (KMA)
The toxin-gene panel FASTA (--paloc_db) is indexed once with kma index, then clean reads are mapped with kma -ipe … -1t1 -ef -mem_mode. This yields per-template .res (identity, template coverage) and .mapstat (fragment counts) files. The panel includes decoy LCT-homolog templates (tpeL, tcsL, tcsH, tcnA by default) so that cross-mapping from related large clostridial toxins can be flagged rather than mis-called.
C. difficile genome coverage (Bowtie2 + SAMtools)
Clean reads are aligned to the C. difficile reference with Bowtie2 (--very-sensitive --no-unal), sorted, and passed through samtools depth to compute the median genome depth per sample.
Normalisation & detection calling
bin/normalise_paloc.py merges the KMA .res/.mapstat outputs with the genome depth values and applies the calling thresholds: a toxin gene is called detected only if it meets --min_breadth and --min_identity, and the toxigenic fraction is reported only when the C. diff median depth clears --min_cdiff_depth. Decoy templates are handled separately. Output is written as a long table, a wide summary, and a detection summary.
Reference-guided assembly — optional (--run_assembly, default true)
Clean reads are mapped to the C. diff reference with BWA-MEM2; properly-paired mapped reads are extracted back to FASTQ with SAMtools and assembled with metaSPAdes (--careful). Samples with fewer than 1,000 mapped read pairs skip assembly gracefully. Resulting contigs/scaffolds are evaluated with QUAST.
Antimicrobial-resistance genes (RGI + CARD)
The CARD JSON is loaded once into a local RGI database (rgi load --local), then each sample is profiled with rgi bwt using the configured aligner (--rgi_aligner, default kma). Per-sample gene- and stats-level outputs are merged and normalised by bin/merge_and_normalize.py into a combined ARG abundance matrix.
MultiQC
A MultiQC process is included to aggregate fastp and QUAST reports into a single HTML report. (Note: the MultiQC aggregation call is currently commented out in main.nf — uncomment it to enable.)
A CSV file with three columns — no spaces, no extra headers:
sample,read1,read2
CD001,/data/CD001_R1.fastq.gz,/data/CD001_R2.fastq.gz
CD002,/data/CD002_R1.fastq.gz,/data/CD002_R2.fastq.gz| Column | Description |
|---|---|
sample |
Unique sample identifier |
read1 |
Absolute path to forward reads (gzipped FASTQ) |
read2 |
Absolute path to reverse reads (gzipped FASTQ) |
| Parameter | Default | Description |
|---|---|---|
--input |
required | Path to samplesheet CSV |
--outdir |
./results |
Output directory |
--paloc_db |
required | FASTA of PaLoc/CdtLoc genes + decoy LCT homologs (nucleotide) |
--cdiff_ref |
required | C. difficile reference genome FASTA |
--card_json |
required | CARD database JSON (for RGI) |
--host_decontamination |
true |
Enable human + PhiX removal |
--host_genome_index |
— | Path to Bowtie2 human genome index directory |
--bowtie_phix_index |
— | Path to Bowtie2 PhiX index directory |
--min_breadth |
80.0 |
Min % of a gene covered to call it "detected" |
--min_identity |
95.0 |
Min % template identity to call "detected" |
--min_cdiff_depth |
5.0 |
Min C. diff median depth to report toxigenic fraction |
--decoys |
tpeL,tcsL,tcsH,tcnA |
Substrings flagging decoy templates |
--run_assembly |
true |
Enable reference-guided assembly + QUAST |
--rgi_aligner |
kma |
Aligner used by rgi bwt |
--singularity_cache_dir |
— | Directory for cached container images |
--max_cpus |
16 |
Global CPU ceiling |
--max_memory |
128.GB |
Global memory ceiling |
--max_time |
48.h |
Global time ceiling |
fastp parameters (key options — see modules/preprocessing.nf for the full list):
| Parameter | Default | Description |
|---|---|---|
--fastp_qualified_quality_phred |
15 |
Minimum base quality score (Q15) |
--fastp_unqualified_percent_limit |
40 |
Max % of bases below quality threshold |
--fastp_n_base_limit |
5 |
Max N bases allowed |
--fastp_length_required |
15 |
Minimum read length to keep |
--fastp_detect_adapter_for_pe |
true |
Auto-detect adapters for paired-end data |
--fastp_dedup |
false |
Remove duplicate reads |
--fastp_low_complexity_filter |
false |
Filter low-complexity reads |
--fastp_correction |
false |
Enable overlap-based base correction |
--fastp_trim_poly_g |
false |
Trim polyG tails (NextSeq/NovaSeq) |
--fastp_extra_args |
"" |
Any additional fastp arguments |
results/
├── fastp/ # Trimmed reads + fastp JSON/HTML reports
├── hostile/ # Host/PhiX-decontaminated, re-synced reads (if --host_decontamination)
├── kma/ # Per-sample KMA .res / .mapstat (PaLoc panel)
├── cdiff_cov/ # C. diff median depth + normalised PaLoc result tables
├── Assembly/ # metaSPAdes contigs & scaffolds (if --run_assembly)
├── assembly/ # QUAST reports (if --run_assembly)
├── rgi_db/ # Local CARD database built by RGI
├── rgi_bwt/ # Per-sample RGI allele/gene/stats outputs
├── arg_merged/ # Merged + normalised ARG matrix
├── multiqc/ # Aggregated MultiQC report (if enabled)
└── pipeline_info/ # Nextflow timeline, report, trace, DAG
| File | Description |
|---|---|
cdiff_cov/paloc_results_long.tsv |
Long-format per-template KMA results merged with genome depth |
cdiff_cov/paloc_summary_wide.tsv |
Wide summary table (templates × samples) |
cdiff_cov/paloc_detection_summary.tsv |
Decoy-aware detection calls per toxin gene |
arg_merged/arg_matrix_normalized.tsv |
Normalised ARG abundance matrix |
arg_merged/arg_matrix_raw.tsv |
Raw ARG count matrix |
arg_merged/arg_gene_class_lookup.tsv |
ARG → drug-class lookup table |
assembly/*.quast_report.txt |
Per-sample QUAST assembly metrics |
pipeline_info/report.html |
Nextflow execution report |
- fastp: Chen S et al. (2018). fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics, 34(17):i884–i890. https://doi.org/10.1093/bioinformatics/bty560
- Bowtie2: Langmead B, Salzberg SL. (2012). Fast gapped-read alignment with Bowtie 2. Nature Methods, 9:357–359. https://doi.org/10.1038/nmeth.1923
- KMA: Clausen PTLC, Aarestrup FM, Lund O. (2018). Rapid and precise alignment of raw reads against redundant databases with KMA. BMC Bioinformatics, 19:307. https://doi.org/10.1186/s12859-018-2336-6
- SAMtools: Danecek P et al. (2021). Twelve years of SAMtools and BCFtools. GigaScience, 10(2):giab008. https://doi.org/10.1093/gigascience/giab008
- BWA-MEM2: Vasimuddin M et al. (2019). Efficient architecture-aware acceleration of BWA-MEM for multicore systems. IEEE IPDPS. https://doi.org/10.1109/IPDPS.2019.00041
- SPAdes: Prjibelski A et al. (2020). Using SPAdes de novo assembler. Current Protocols in Bioinformatics, 70:e102. https://doi.org/10.1002/cpbi.102
- QUAST: Gurevich A et al. (2013). QUAST: quality assessment tool for genome assemblies. Bioinformatics, 29(8):1072–1075. https://doi.org/10.1093/bioinformatics/btt086
- RGI / CARD: Alcock BP et al. (2023). CARD 2023: expanded curation, support for machine learning, and resistome prediction at the Comprehensive Antibiotic Resistance Database. Nucleic Acids Research, 51(D1):D690–D699. https://doi.org/10.1093/nar/gkac920
- MultiQC: Ewels P et al. (2016). MultiQC: summarize analysis results for multiple tools and samples in a single report. Bioinformatics, 32(19):3047–3048. https://doi.org/10.1093/bioinformatics/btw354
- Nextflow: Di Tommaso P et al. (2017). Nextflow enables reproducible computational workflows. Nature Biotechnology, 35:316–319. https://doi.org/10.1038/nbt.3820