MGnify amplicon analysis pipeline

This repository contains the beta version of the v6.0 MGnify amplicon analysis pipeline. It is, first and foremost, a refactor of the existing v5.0 amplicon analysis pipeline, replacing CWL with Nextflow as its workflow management system. This pipeline re-implements all existing closed-reference v5.0 features, and makes multiple significant changes and additions.

Pipeline description

Features

The amplicon analysis pipeline v6.0 re-implements all of the existing features from v5.0:

• Reads quality control
• rRNA sequence extraction using Infernal/cmsearch
• Closed-reference-based taxonomic classification and visualiation of rRNA using MAPseq and Krona

The amplicon analysis pipeline v6.0 also contains multiple significant changes:

• Refactoring from CWL to Nextflow for pipeline definition
• Simplification the reads quality control using fastp
• Automatic amplified region inference for 16S and 18S rRNA
• Automatic primer identification, trimming, and validation
• Addition of Amplicon Sequence Variant (ASV) calling using DADA2
• Taxonomic classification and visualisation of ASVs using MAPseq and Krona to complement the existing closed-reference analysis
• Addition of PR2 as a reference database
• Updating of existing reference databases (SILVA, UNITE, ITSoneDB, Rfam)

Valid amplicons

At this stage, the only sequence amplicons that this pipeline is built for are:

Amplicon	Closed-reference analysis	ASV analysis
16S	✓	✓
18S	✓	✓
LSU	✓	✗
ITS	✓	✗

Tools

Tool	Version	Purpose
fastp	0.23.4	Read quality control
seqtk	1.3-r106	FASTQ file manipulation
easel	0.49	FASTA file manipulation
bedtools	2.30.0	FASTA sequence masking
Infernal/cmsearch	1.1.5	rRNA sequence searching
cmsearch_tblout_deoverlap	0.09	Deoverlapping of cmsearch results
MAPseq	2.1.1b	Reference-based taxonomic classification of rRNA
Krona	2.8.1	Krona chart visualisation
cutadapt	4.6	Primer trimming
R	4.3.3	R programming language (runs DADA2)
DADA2	1.30.0	ASV calling
mgnify-pipelines-toolkit	0.1.5	Toolkit containing various in-house processing scripts

Reference databases

This pipeline uses five different reference databases. The files the pipeline uses are processed from the raw files available on each database's website, for use by MAPseq and cmsearch. We provide ready-made versions of these processed files on our FTP, which you can find here:

Reference database	Version	Purpose	Processed file paths
SILVA	138.1	16S+18S+LSU rRNA database	https://ftp.ebi.ac.uk/pub/databases/metagenomics/pipelines/tool-dbs/silva-ssu/ https://ftp.ebi.ac.uk/pub/databases/metagenomics/pipelines/tool-dbs/silva-lsu/
PR2	5.0	Protist-focused 18S+16S rRNA database	https://ftp.ebi.ac.uk/pub/databases/metagenomics/pipelines/tool-dbs/pr2/
UNITE	9.0	ITS database	https://ftp.ebi.ac.uk/pub/databases/metagenomics/pipelines/tool-dbs/unite/
ITSoneDB	1.141	ITS database	https://ftp.ebi.ac.uk/pub/databases/metagenomics/pipelines/tool-dbs/itsonedb/
Rfam	14.10	RNA family profile database	https://ftp.ebi.ac.uk/pub/databases/metagenomics/pipelines/tool-dbs/rfam/

How to run

Requirements

At the moment the only prerequisites for running it are Nextflow and Docker/Singularity, since all of the Nextflow processes use pre-built containers.

Input shape

The input data for the pipeline is amplicon sequencing reads (either paired-end or single-end) in the form of FASTQ files. These files should be specified using a .csv samplesheet file with this format:

sample,fastq_1,fastq_2,single_end
SRR9674618,/path/to/reads/SRR9674618.fastq.gz,,true
SRR17062740,/path/to/reads/SRR17062740_1.fastq.gz,/path/to/reads/SRR17062740_2.fastq.gz,false

Execution

You can run the current version of the pipeline on SLURM like this:

nextflow run -profile codon_slurm main.nf --input /path/to/samplesheet.csv --outdir /path/to/outdir main.nf

Output directory structure

Example output directory structure for one run:

├── amplified-region-inference
│   ├── ERR2356289_16S-V4-V5_extracted.fastq.gz
│   ├── ERR2356289.16S.V4-V5.txt
│   └── ERR2356289.tsv
├── asv
│   ├── 16S-V4-V5
│   │   ├── ERR2356289_16S.V4-V5_asv_read_counts.tsv
│   │   ├── ERR2356289_16S.V4-V5_pr2_asv_krona_counts.txt
│   │   └── ERR2356289_16S.V4-V5_silva_asv_krona_counts.txt
│   ├── ERR2356289_asv_seqs.fasta
│   ├── ERR2356289_DADA2-PR2_asv_taxa.tsv
│   ├── ERR2356289_DADA2-SILVA_asv_taxa.tsv
│   ├── ERR2356289_proportion_chimeric.txt
│   └── ERR2356289_proportion_matched.txt
├── primer-identification
│   ├── ERR2356289_16S-V4-V5_auto_primers.fasta
│   ├── ERR2356289_16S-V4-V5_std_primers.fasta
│   ├── ERR2356289_1.trim.fastq.gz
│   ├── ERR2356289_2.trim.fastq.gz
│   ├── ERR2356289_final_concat_primers.fasta
│   ├── ERR2356289_rev_comp_se_primers.fasta
│   └── ERR2356289_trimming_conductor_out_16S-V4-V5.txt
├── qc
│   ├── ERR2356289.merged.fastq.gz
│   └── ERR2356289.seqtk-seq.fasta.gz
├── sequence-categorisation
│   ├── ERR2356289_SSU.fasta
│   └── ERR2356289.tblout.deoverlapped
└── taxonomy-summary
    ├── DADA2-PR2
    │   └── ERR2356289_16S-V4-V5.html
    ├── DADA2-SILVA
    │   └── ERR2356289_16S-V4-V5.html
    ├── ITSoneDB
    │   ├── ERR2356289.html
    │   └── ERR2356289.txt
    ├── PR2
    │   ├── ERR2356289.html
    │   └── ERR2356289.txt
    ├── SSU
    │   ├── ERR2356289.html
    │   └── ERR2356289.txt
    └── UNITE
        ├── ERR2356289.html
        └── ERR2356289.txt