diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/CHANGELOG.md b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/CHANGELOG.md new file mode 100644 index 00000000..cd08132c --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/CHANGELOG.md @@ -0,0 +1,4 @@ +# Workflow change log + +## [1.0.0](https://github.com/nasa/GeneLab_Data_Processing/tree/NF_Amp454IonTor_1.0.0/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor) +- workflow version that conxverted snakemake to nextflow \ No newline at end of file diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/README.md b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/README.md new file mode 100644 index 00000000..30e3bcb1 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/README.md @@ -0,0 +1,194 @@ +# Workflow Information and Usage Instructions + +## General Workflow Info + +### Implementation Tools + +The current GeneLab 454 and IonTorrent amplicon sequencing data processing pipeline (Amp454IonTor), [GL-DPPD-7106.md](../../Pipeline_GL-DPPD-7106_Versions/GL-DPPD-7106.md), is implemented as a [Nextflow](https://nextflow.io/) DSL2 workflow and utilizes [Singularity](https://docs.sylabs.io/guides/3.10/user-guide/introduction.html) containers or [conda](https://docs.conda.io/en/latest/) environments to install/run all tools. This workflow is run using the command line interface (CLI) of any unix-based system. While knowledge of creating workflows in nextflow is not required to run the workflow as is, [the Nextflow documentation](https://nextflow.io/docs/latest/index.html) is a useful resource for users who want to modify and/or extend this workflow. + +## Utilizing the Workflow + +1. [Install Nextflow and Singularity](#1-install-nextflow-and-singularity) + 1a. [Install Nextflow](#1a-install-nextflow) + 1b. [Install Singularity](#1b-install-singularity) + +2. [Download the workflow files](#2-download-the-workflow-files) + +3. [Fetch Singularity Images](#3-fetch-singularity-images) + +4. [Run the workflow](#4-run-the-workflow) + 4a. [Approach 1: Run slurm jobs in singularity containers with OSD accession as input](#4a-approach-1-run-slurm-jobs-in-singularity-containers-with-osd-accession-as-input) + 4b. [Approach 2: Run slurm jobs in singularity containers with a csv file as input](#4b-approach-2-run-slurm-jobs-in-singularity-containers-with-a-csv-file-as-input) + 4c. [Approach 3: Run jobs locally in conda environments and specify the path to one or more existing conda environments](#4c-approach-run-jobs-locally-in-conda-environments-and-specify-the-path-to-one-or-more-existing-conda-environments) + 4d. [Modify parameters and cpu resources in the nextflow config file](#4d-modify-parameters-and-cpu-resources-in-the-nextflow-config-file) + +5. [Workflow outputs](#5-workflow-outputs) + 5a. [Main outputs](#5a-main-outputs) + 5b. [Resource logs](#5b-resource-logs) + +
+ +--- + +### 1. Install Nextflow and Singularity + +#### 1a. Install Nextflow + +Nextflow can be installed either through [Anaconda](https://anaconda.org/bioconda/nextflow) or as documented on the [Nextflow documentation page](https://www.nextflow.io/docs/latest/getstarted.html). + +> Note: If you want to install Anaconda, we recommend installing a Miniconda, Python3 version appropriate for your system, as instructed by [Happy Belly Bioinformatics](https://astrobiomike.github.io/unix/conda-intro#getting-and-installing-conda). +> +> Once conda is installed on your system, you can install the latest version of Nextflow by running the following commands: +> +> ```bash +> conda install -c bioconda nextflow +> nextflow self-update +> ``` + +
+ +#### 1b. Install Singularity + +Singularity is a container platform that allows usage of containerized software. This enables the GeneLab workflow to retrieve and use all software required for processing without the need to install the software directly on the user's system. + +We recommend installing Singularity on a system wide level as per the associated [documentation](https://docs.sylabs.io/guides/3.10/admin-guide/admin_quickstart.html). + +> Note: Singularity is also available through [Anaconda](https://anaconda.org/conda-forge/singularity). + +
+ +--- + +### 2. Download the workflow files + +All files required for utilizing the NF_XXX GeneLab workflow for processing 454 ion torrent data are in the [workflow_code](workflow_code) directory. To get a copy of latest *NF_XXX* version on to your system, the code can be downloaded as a zip file from the release page then unzipped after downloading by running the following commands: + +```bash +wget https://github.com/nasa/GeneLab_Data_Processing/releases/download/NF_Amp454IonTor/NF_Amp454IonTor.zip +unzip NF_Amp454IonTor.zip && cd NF_XXX-X_X.X.X +``` + +OR by using the genelab-utils conda package + +```bash +GL-get-workflow Amplicon-454-IonTorrent +``` + +
+ +--- + +### 3. Fetch Singularity Images + +Although Nextflow can fetch Singularity images from a url, doing so may cause issues as detailed [here](https://github.com/nextflow-io/nextflow/issues/1210). + +To avoid this issue, run the following command to fetch the Singularity images prior to running the NF_Amp454IonTor workflow: + +> Note: This command should be run in the location containing the `NF_Amp454IonTor` directory that was downloaded in [step 2](#2-download-the-workflow-files) above. + +```bash +bash ./bin/prepull_singularity.sh nextflow.config +``` + +Once complete, a `singularity` folder containing the Singularity images will be created. Run the following command to export this folder as a Nextflow configuration environment variable to ensure Nextflow can locate the fetched images: + +```bash +export NXF_SINGULARITY_CACHEDIR=$(pwd)/singularity +``` + +
+ +--- + + +### 4. Run the Workflow + +For options and detailed help on how to run the workflow, run the following command: + +```bash +nextflow run main.nf --help +``` + +> Note: Nextflow commands use both single hyphen arguments (e.g. -help) that denote general nextflow arguments and double hyphen arguments (e.g. --csv_file) that denote workflow specific parameters. Take care to use the proper number of hyphens for each argument. + +
+ +#### 4a. Approach 1: Run slurm jobs in singularity containers with OSD accession as input + +```bash +nextflow run main.nf -resume -profile slurm,singularity --GLDS_accession OSD-72 --target_region 16S --min_bbduk_len 50 --min_bbduk_avg_quality 15 +``` + +
+ +#### 4b. Approach 2: Run slurm jobs in singularity containers with a csv file as input + +```bash +nextflow run main.nf -resume -profile slurm,singularity --csv_file file.csv --target_region 16S --F_primer AGAGTTTGATCCTGGCTCAG --R_primer CTGCCTCCCGTAGGAGT --min_bbduk_len 50 --min_bbduk_avg_quality 15 +``` + +
+ +#### 4c. Approach 3: Run jobs locally in conda environments and specify the path to one or more existing conda environment(s) + +```bash +nextflow run main.nf -resume -profile conda --csv_file file.csv --target_region 16S --F_primer AGAGTTTGATCCTGGCTCAG --R_primer CTGCCTCCCGTAGGAGT --min_bbduk_len 50 --min_bbduk_avg_quality 15 --conda.qc +``` + +
+ +**Required Parameters For All Approaches:** + +* `-run main.nf` - Instructs nextflow to run the NF_XXX workflow +* `-resume` - Resumes workflow execution using previously cached results +* `-profile` – Specifies the configuration profile(s) to load, `singularity` instructs nextflow to setup and use singularity for all software called in the workflow +* `--target_region` – Specifies the amplicon target region to be analyzed, 16S or ITS. +* `--min_bbduk_len` – Specifies the minimum read length to retain after filtering with bbduk. +* `--min_bbduk_avg_quality` – Specifies the minimum average read quality for bbduk read filtering. + + + + *Required only if you would like to pull and process data directly from OSDR* + +* `--GLDS_accession` – A Genelab / OSD accession number e.g. OSD-72. + +*Required only if --GLDS_accession is not passed as an argument* + +* `--csv_file` – A 2-column input file with these headers [sample_id, read]. Please see the sample [file.csv](workflow_code/file.csv) in this repository for an example on how to format this file. + +* `--F_primer` – Forward primer sequence. + +* `--R_primer` – Reverse primer sequence. + +> See `nextflow run -h` and [Nextflow's CLI run command documentation](https://nextflow.io/docs/latest/cli.html#run) for more options and details on how to run nextflow. + +
+ +#### 4d. Modify parameters and cpu resources in the nextflow config file + +Additionally, the parameters and workflow resources can be directly specified in the nextflow.config file. For detailed instructions on how to modify and set parameters in the nextflow.config file, please see the [documentation here](https://www.nextflow.io/docs/latest/config.html). + +Once you've downloaded the workflow template, you can modify the parameters in the `params` scope and cpus/memory requirements in the `process` scope in your downloaded version of the [nextflow.config](workflow_code/nextflow.config) file as needed in order to match your dataset and system setup. For example, you can directly set the the full paths to available conda environments in the `conda` scope within the `params` scope. Additionally, if necessary, you'll need to modify each variable in the [nextflow.config](workflow_code/nextflow.config) file to be consistent with the study you want to process and the machine you're using. + +
+ +--- + +### 5. Workflow outputs + +#### 5a. Main outputs + +The outputs from this pipeline are documented in the [GL-DPPD-7106](../../Pipeline_GL-DPPD-7106_Versions/GL-DPPD-7106.md) processing protocol. + +#### 5b. Resource logs + +Standard nextflow resource usage logs are also produced as follows: + +- Output: + - Resource_Usage/execution_report_{timestamp}.html (an html report that includes metrics about the workflow execution including computational resources and exact workflow process commands) + - Resource_Usage/execution_timeline_{timestamp}.html (an html timeline for all processes executed in the workflow) + - Resource_Usage/execution_trace_{timestamp}.txt (an execution tracing file that contains information about each process executed in the workflow, including: submission time, start time, completion time, cpu and memory used, machine-readable output) + +> Further details about these logs can also found within [this Nextflow documentation page](https://www.nextflow.io/docs/latest/tracing.html#execution-report). + + diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/454-IonTorrent-R-processing.R b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/454-IonTorrent-R-processing.R new file mode 100644 index 00000000..050e848d --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/454-IonTorrent-R-processing.R @@ -0,0 +1,101 @@ +#!/usr/bin/env Rscript +################################################################################## +## R processing script for 454/Ion Torrent amplicon data ## +## Developed by Michael D. Lee (Mike.Lee@nasa.gov) ## +################################################################################## + +# as called from the associated Snakefile, this expects to be run as: Rscript full-R-processing.R + +# setting variables used within R: +args <- commandArgs(trailingOnly = TRUE) + +suppressWarnings(otus_fasta_path <- args[1]) +suppressWarnings(trimmed_dir <- args[2]) +suppressWarnings(filtered_dir <- args[3]) +suppressWarnings(final_outputs_dir <- args[4]) +suppressWarnings(output_prefix <- args[5]) +suppressWarnings(target_region <- args[6]) +suppressWarnings(assay_suffix <- args[7]) + + # loading libraries +library(DECIPHER) +library(biomformat) + +# Set default internet timeout to 1 hour +options(timeout=3600) + + + ### assigning taxonomy ### + # reading OTUs into a DNAStringSet object +dna <- readDNAStringSet(paste0(final_outputs_dir, output_prefix, "OTUs", assay_suffix, ".fasta")) + + + # downloading reference R taxonomy object +cat("\n\n Downloading reference database...\n\n") + +if ( target_region == "16S" ) { + + download.file("https://figshare.com/ndownloader/files/46245217", "SILVA_SSU_r138_2019.RData") + load("SILVA_SSU_r138_2019.RData") + #file.remove("SILVA_SSU_r138_2019.RData") + +} else if ( target_region == "ITS" ) { + download.file("https://figshare.com/ndownloader/files/46245586", "UNITE_v2020_February2020.RData") + load("UNITE_v2020_February2020.RData") + #file.remove("UNITE_v2020_February2020.RData") +} + + +# assigning taxonomy +cat("\n\n Assigning taxonomy...\n\n") + +tax_info <- IdTaxa(dna, trainingSet, strand="both", processors=NULL) + +cat("\n\n Making and writing out tables...\n\n") + + # making and writing out a taxonomy table: + # creating vector of desired ranks +ranks <- c("domain", "phylum", "class", "order", "family", "genus", "species") + + # creating table of taxonomy and setting any that are unclassified as "NA" +tax_tab <- t(sapply(tax_info, function(x) { + m <- match(ranks, x$rank) + taxa <- x$taxon[m] + taxa[startsWith(taxa, "unclassified_")] <- NA + taxa +})) + +colnames(tax_tab) <- ranks +row.names(tax_tab) <- NULL +otu_ids <- names(tax_info) +tax_tab <- data.frame("OTU_ID"=otu_ids, tax_tab, check.names=FALSE) + +write.table(tax_tab, paste0(final_outputs_dir, output_prefix,"taxonomy", assay_suffix, ".tsv"), sep = "\t", quote=F, row.names=FALSE) + + # reading in counts table to generate other outputs +otu_tab <- read.table(paste0(final_outputs_dir, output_prefix, "counts", assay_suffix, ".tsv"), sep="\t", header=TRUE, check.names=FALSE) + + # generating and writing out biom file format +biom_object <- make_biom(data=otu_tab, observation_metadata=tax_tab) +write_biom(biom_object, paste0(final_outputs_dir, output_prefix, "taxonomy-and-counts", assay_suffix, ".biom")) + + # making a tsv of combined tax and counts +tax_and_count_tab <- merge(tax_tab, otu_tab) +write.table(tax_and_count_tab, paste0(final_outputs_dir, output_prefix, "taxonomy-and-counts", assay_suffix, ".tsv"), sep="\t", quote=FALSE, row.names=FALSE) + +# making final count summary table +cutadapt_tab <- read.table(paste0(trimmed_dir, output_prefix, "trimmed-read-counts", assay_suffix, ".tsv"), sep="\t", header=TRUE) +bbduk_tab <- read.table(paste0(filtered_dir, output_prefix, "filtered-read-counts", assay_suffix, ".tsv"), sep="\t", header=TRUE)[,c(1,3)] + # re-reading in counts table to this time set first col as rownames (rather than doing it another way) +otu_tab <- read.table(paste0(final_outputs_dir, output_prefix, "counts", assay_suffix, ".tsv"), sep="\t", header=TRUE, check.names=FALSE, row.names = 1) +mapped_sums <- colSums(otu_tab) +mapped_tab <- data.frame(sample=names(mapped_sums), mapped_to_OTUs=mapped_sums, row.names=NULL) + +t1 <- merge(cutadapt_tab, bbduk_tab) +count_summary_tab <- merge(t1, mapped_tab) +count_summary_tab$final_perc_reads_retained <- round(count_summary_tab$mapped_to_OTUs / count_summary_tab$raw_reads * 100, 2) + +write.table(count_summary_tab, paste0(final_outputs_dir, output_prefix, "read-count-tracking", assay_suffix, ".tsv"), sep="\t", quote=FALSE, row.names=FALSE) + +cat("\n\n Session info:\n\n") +sessionInfo() diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/clean-paths.sh b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/clean-paths.sh new file mode 100644 index 00000000..416758a2 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/clean-paths.sh @@ -0,0 +1,24 @@ +#!/usr/bin/env bash +set -e + +# only built for use on N288 cluster + +# example usage: bash clean-paths.sh + +# making sure by chance we are not overwriting a wanted file called 't' + +if [ -s t ]; then + printf "\n This simple program temporarily writes to a file called 't'\n" + printf " Since that exists already here, we are not going to continue.\n\n" + exit +fi + + +ROOT_DIR=$(echo $2 | awk '{N=split($0,a,"/"); for(i=0; i < N-1; i++) printf "%s/", a[i]}' | sed 's|//|/|') + + +sed -E 's|.*/GLDS_Datasets/(.+)|\1|g' ${1} \ + | sed -E 's|.+/miniconda.+/envs/[^/]*/||g' \ + | sed -E 's|/[^ ]*/GLDS-|GLDS-|g' \ + | sed -E 's|/[a-z]{6}/[^ ]*||g' \ + | sed -E "s|${ROOT_DIR}||g" > t && mv t ${1} \ No newline at end of file diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/create_runsheet.sh b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/create_runsheet.sh new file mode 100644 index 00000000..0ed6b216 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/create_runsheet.sh @@ -0,0 +1,15 @@ +#!/usr/bin/env bash + +# A script to the the input csv file rtequired by this pipeline when a GLDS accession is provided +# Rather than the required input csv file + +ASSAY_TABLE=$1 # target_region_assay_table.txt + +cat ${ASSAY_TABLE} | \ +sed 's/"//g' | \ +awk -v PWD=$PWD 'BEGIN{FS="\t"; OFS=","; print "sample_id,read,f_primer,r_primer,target_region"} \ + NR==1{ for(i=1; i<=NF; i++) header[$i]=i} \ + NR>1{ split ($header["Parameter Value[Primer Info]"], primers, ","); \ + printf "%s,%s/Raw_Sequence_Data/%s,%s,%s,%s\n", $header["Sample Name"], PWD, $header["Raw Data File"], primers[2], primers[4],$header["Parameter Value[Target Molecule]"] }' | \ + sed -E "s/5'-([A-Z]+)-3'/\1/g" | \ + sed -E 's/\s+//g' diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/get_R_package_version.R b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/get_R_package_version.R new file mode 100644 index 00000000..fb5ad6f0 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/get_R_package_version.R @@ -0,0 +1,10 @@ +#!/usr/bin/env Rscript + +# Get versions +VERSIONS <- sprintf("DECIPHER %s\nbiomformat %s\n", + packageVersion("DECIPHER"), + packageVersion("biomformat")) + +# Write versions to file + +write(x= VERSIONS, file="versions.txt", append=TRUE) diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/prepull_singularity.sh b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/prepull_singularity.sh new file mode 100644 index 00000000..125130f1 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/bin/prepull_singularity.sh @@ -0,0 +1,31 @@ +#!/usr/bin/env bash + +# Addresses issue: https://github.com/nextflow-io/nextflow/issues/1210 + +CONFILE=${1:-nextflow.config} +OUTDIR=${2:-./singularity} + +if [ ! -e $CONFILE ]; then + echo "$CONFILE does not exist" + exit +fi + +TMPFILE=`mktemp` + +CURDIR=$(pwd) + +mkdir -p $OUTDIR + +cat ${CONFILE}|grep 'container'|perl -lane 'if ( $_=~/container\s*\=\s*\"(\S+)\"/ ) { $_=~/container\s*\=\s*\"(\S+)\"/; print $1 unless ( $1=~/^\s*$/ or $1=~/\.sif/ or $1=~/\.img/ ) ; }' > $TMPFILE + +cd ${OUTDIR} + +while IFS= read -r line; do + name=$line + name=${name/:/-} + name=${name//\//-} + echo $name + singularity pull ${name}.img docker://$line +done < $TMPFILE + +cd $CURDIR \ No newline at end of file diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/R.yaml b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/R.yaml new file mode 100644 index 00000000..c315c5a2 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/R.yaml @@ -0,0 +1,8 @@ +channels: + - conda-forge + - bioconda + - defaults +dependencies: + - r-base=4.1.1 + - bioconductor-decipher=2.20.0 + - bioconductor-biomformat=1.20.0 diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/bbmap.yaml b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/bbmap.yaml new file mode 100644 index 00000000..034320d5 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/bbmap.yaml @@ -0,0 +1,6 @@ +channels: + - conda-forge + - bioconda + - defaults +dependencies: + - bbmap=38.86 diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/cutadapt.yaml b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/cutadapt.yaml new file mode 100644 index 00000000..f9d073f3 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/cutadapt.yaml @@ -0,0 +1,6 @@ +channels: + - conda-forge + - bioconda + - defaults +dependencies: + - cutadapt=1.16 diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/genelab.yaml b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/genelab.yaml new file mode 100644 index 00000000..9f2ec80b --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/genelab.yaml @@ -0,0 +1,8 @@ +name: genelab-utils +channels: + - conda-forge + - bioconda + - defaults + - astrobiomike +dependencies: + - genelab-utils==1.3.22=py312_1 \ No newline at end of file diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/qc.yaml b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/qc.yaml new file mode 100644 index 00000000..2648f294 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/qc.yaml @@ -0,0 +1,9 @@ +channels: + - conda-forge + - bioconda + - defaults +dependencies: + - fastqc=0.12.1 + - multiqc=1.19 + - zip=3.0 + - python=3.8 diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/vsearch.yaml b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/vsearch.yaml new file mode 100644 index 00000000..cbe534e7 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/envs/vsearch.yaml @@ -0,0 +1,8 @@ +channels: + - conda-forge + - bioconda + - defaults + - astrobiomike +dependencies: + - vsearch=2.15.2 + - bit=1.8.37 diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/file.csv b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/file.csv new file mode 100644 index 00000000..eaa63bc2 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/file.csv @@ -0,0 +1,10 @@ +sample_id,read +SAMN03652399,path/to/Raw_Data/SAMN03652399_raw.fastq.gz +SAMN03652400,path/to/Raw_Data/SAMN03652400_raw.fastq.gz +SAMN03652401,path/to/Raw_Data/SAMN03652401_raw.fastq.gz +SAMN03652402,path/to/Raw_Data/SAMN03652402_raw.fastq.gz +SAMN03652403,path/to/Raw_Data/SAMN03652403_raw.fastq.gz +SAMN03652404,path/to/Raw_Data/SAMN03652404_raw.fastq.gz +SAMN03652405,path/to/Raw_Data/SAMN03652405_raw.fastq.gz +SAMN03652406,path/to/Raw_Data/SAMN03652406_raw.fastq.gz +SAMN03652407,path/to/Raw_Data/SAMN03652407_raw.fastq.gz \ No newline at end of file diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/main.nf b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/main.nf new file mode 100644 index 00000000..cd44823a --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/main.nf @@ -0,0 +1,242 @@ +#!/usr/bin/env nextflow +nextflow.enable.dsl = 2 + +// color defs +c_back_bright_red = "\u001b[41;1m"; +c_bright_green = "\u001b[32;1m"; +c_blue = "\033[0;34m"; +c_reset = "\033[0m"; + +params.help = false +params.debug = false + +/************************************************** +* HELP MENU ************************************** +**************************************************/ + +if (params.help) { + println() + println("Nextflow Amp454IonTor Consensus Pipeline: $workflow.manifest.version") + println("USAGE:") + println("Example 1: Submit and run jobs with slurm in singularity containers with OSD accession as input.") + println(" > nextflow run main.nf -resume -profile slurm,singularity --GLDS_accession OSD-72 --target_region 16S --min_bbduk_len 50 --min_bbduk_avg_quality 15") + println() + println("Example 2: : Submit and run jobs with slurm in conda environments.") + println(" > nextflow run main.nf -resume -profile slurm,singularity --csv_file file.csv --target_region ITS --F_primer TCCGTAGGTGAACCTGCGG --R_primer GCTGCGTTCTTCATCGATGC --min_bbduk_len 50 --min_bbduk_avg_quality 15") + println() + println("Example 3: Run jobs locally in conda environments and specify the path to an existing conda environment") + println(" > nextflow run main.nf -resume -profile conda --csv_file file.csv --target_region 16S --F_primer AGAGTTTGATCCTGGCTCAG --R_primer CTGCCTCCCGTAGGAGT --min_bbduk_len 50 --min_bbduk_avg_quality 15 --conda.qc ") + println() + println() + println("Required arguments:") + println() + println("""-profile [STRING] Specifies the profile to be used to run the workflow. Options are [singularity, docker, conda, slurm]. + singularity, docker and conda will run the pipeline locally using singularity, docker, and conda, respectively. + You can combine profiles by separating them with comma. For example, to submit and run jobs using slurm in singularity containers pass 'slurm,singularity' as argument. """) + println("--csv_file [PATH] Required only if --GLDS_accession is not set. A 2-column input file with these headers [sample_id, read] e.g. file.csv. The sample_id column should contain unique sample ids while the read column should contain the absolute or relative path to the sample's reads.") + println("--target_region [STRING] Specifies the amplicon target region to be analyzed. options are one of [16S, ITS]. Default: 16S") + println() + println("Cutadapt (trimming) parameters:") + println(" --F_primer [STRING] Required only if --GLDS_accession is not set. Forward primer sequence e.g. AGAGTTTGATCCTGGCTCAG") + println(" --R_primer [STRING] Required only if --GLDS_accession is not set. Reverse primer sequence e.g. CTGCCTCCCGTAGGAGT") + println() + println("BBDUK (filtering) parameters:") + println(" --min_bbduk_len [INT] Specifies the minimum read length threshold for bbduk. Default: 50") + println(" --min_bbduk_avg_quality [INT] Specifies the minimum average quality for bbduk. Default: 15") + println() + println() + println("Optional arguments:") + println() + println(" --help Print this help message and exit") + println(" --debug Show a detailed log of the parameters set by the user when the workflow runs.") + println(" --publishDir_mode [STRING] Specifies how nextflow handles file outputs. Options can be found here https://www.nextflow.io/docs/latest/process.html#publishdir Default: link.") + println(" --errorStrategy [STRING] Specifies how nextflow handles errors. Options can be found here https://www.nextflow.io/docs/latest/process.html#errorstrategy. Default: terminate") + println() + println("File Suffixes:") + println(" --primer_trimmed_suffix [STRING] Suffix to use for naming your primer trimmed reads. Default: _trimmed.fastq.gz") + println(" --filtered_suffix [STRING] Suffix to use for naming your quality filtered reads. Default: _filtered.fastq.gz") + println() + println("Output directories:") + println(" --raw_reads_dir [PATH] Specifies where processed raw reads will be published. Default: ../Raw_Sequence_Data/") + println(" --fastqc_out_dir [PATH] Specifies where fastqc and multiqc outputs will be published. Default: ../FastQC_Outputs/") + println(" --trimmed_reads_dir [PATH] Specifies where cutadapt trimmed reads will be published. Default: ../Trimmed_Sequence_Data/") + println(" --filtered_reads_dir [PATH] Specifies where BBDUK filtered reads will be published. Default: ../Filtered_Sequence_Data/") + println() + println("Genelab specific arguements:") + println(" --GLDS_accession [STRING] A Genelab GLDS or OSD accession number if the --csv_file parameter is not set. If this parameter is set, it will ignore the --csv_file parameter.") + println(" --RawFilePattern [STRING] If we do not want to download all files (which we often won't), we can specify a pattern here to subset the total files.") + println(" For example, if we know we want to download just the fastq.gz files, we can say 'fastq.gz'. We can also provide multiple patterns") + println(" as a comma-separated list. For example, If we want to download the fastq.gz files that also have 'Amplicon', and 'raw' in") + println(" their filenames, we can provide '-p fastq.gz,Amplicon,raw'. Default: null.") + println(" --assay_suffix [STRING] Genelabs assay suffix. Default: _GLAmpSeq.") + println(" --output_prefix [STRING] Unique name to tag on to output files. Default: ''") + println() + println("Paths to existing conda environments to use otherwise a new one will be created using the yaml file in envs/.") + println(" --conda.qc [PATH] Path to a conda environment containing fastqc, multiqc, zip and python. Default: null.") + println(" --conda.R [PATH] Path to a conda environment containing R along with the packages decipher and biomformat installed. Default: null.") + println(" --conda.bbmap [PATH] Path to a conda environment containing bbmap. Default: null.") + println(" --conda.cutadapt [PATH] Path to a conda environment containing cutadapt. Default: null.") + println(" --conda.vsearch [PATH] Path to a conda environment containing vsearch and bit. Default: null.") + println() + print("Advanced users can edit the nextflow.config file for more control over default settings such as container choice, number of cpus, memory per task etc.") + exit 0 + } + + +/************************************************ +*********** Show pipeline parameters ************ +*************************************************/ + +if(params.debug){ + +log.info """ + Nextflow Amp454IonTor Consensus Pipeline: $workflow.manifest.version + You have set the following parameters: + Profile: ${workflow.profile} + GLDS_accession : ${params.GLDS_accession} + Input csv file : ${params.csv_file} + Amplicon target region : ${params.target_region} + Foward Primer: ${params.F_primer} + Reverse Primer: ${params.R_primer} + Minimum read length: ${params.min_bbduk_len}bp + Minimum read quality: ${params.min_bbduk_avg_quality} + Directory publishing mode: ${params.publishDir_mode} + Nextflow Error strategy: ${params.errorStrategy} + + File Suffixes: + Primers Trimmed Reads Suffix: ${params.primer_trimmed_suffix} + Filtered Reads Suffix: ${params.filtered_suffix} + + Output directories: + Raw reads: ${params.raw_reads_dir} + FastQC: ${params.fastqc_out_dir} + Trimmed Reads: ${params.trimmed_reads_dir} + Filtered Reads: ${params.filtered_reads_dir} + + Genelab parameters: + Genelab Assay Suffix: ${params.assay_suffix} + Output prefix: ${params.output_prefix} + + Conda environments: + qc: ${params.conda.qc} + R: ${params.conda.R} + bbmap: ${params.conda.bbmap} + cutadapt: ${params.conda.cutadapt} + vsearch: ${params.conda.vsearch} + """.stripIndent() +} + + +// Create GLDS runsheet +include { GET_RUNSHEET } from "./modules/create_runsheet.nf" + +// Read quality check and filtering +include { FASTQC as RAW_FASTQC ; MULTIQC as RAW_MULTIQC } from './modules/quality_assessment.nf' +include { CUTADAPT; COMBINE_CUTADAPT_LOGS_AND_SUMMARIZE } from './modules/quality_assessment.nf' +include { BBDUK; COMBINE_BBDUK_LOGS_AND_SUMMARIZE } from './modules/quality_assessment.nf' +include { FASTQC as FILTERED_FASTQC ; MULTIQC as FILTERED_MULTIQC } from './modules/quality_assessment.nf' +include { pick_otus } from './modules/vsearch.nf' +include { RUN_R} from './modules/assign_taxonomy.nf' +include { ZIP_BIOM } from './modules/zip_biom.nf' + + +workflow { + + // Capture software versions + software_versions_ch = Channel.empty() + + if(params.GLDS_accession){ + + GET_RUNSHEET(params.GLDS_accession, params.target_region) + GET_RUNSHEET.out.input_file + .splitCsv(header:true) + .set{file_ch} + + GET_RUNSHEET.out.params_file + .splitCsv(header:true) + .set{params_ch} + + target_region = params_ch.map{row -> "${row.target_region}"}.first() + primers_ch = params_ch.map{ + row -> ["${row.f_primer}", "${row.r_primer}"] + }.first() + + GET_RUNSHEET.out.version | mix(software_versions_ch) | set{software_versions_ch} + + }else{ + + Channel.fromPath(params.csv_file, checkIfExists: true) + .splitCsv(header:true) + .set{file_ch} + + } + + + file_ch.map{row -> tuple( "${row.sample_id}", [file("${row.read}", checkIfExists: true)] )} + .set{reads_ch} + + // Generating a file with sample ids on a new line + file_ch.map{row -> "${row.sample_id}"} + .collectFile(name: "${baseDir}/unique-sample-IDs.txt", newLine: true) + .set{sample_ids_ch} + + // Read quality check and trimming + RAW_FASTQC(reads_ch) + raw_fastqc_files = RAW_FASTQC.out.html.flatten().collect() + RAW_MULTIQC("raw", raw_fastqc_files) + + // Trim reads + if(!params.GLDS_accession) primers_ch = Channel.value([params.F_primer, params.R_primer]) + CUTADAPT(reads_ch, primers_ch) + trim_counts = CUTADAPT.out.trim_counts.map{ sample_id, count -> file("${count}")}.collect() + trim_logs = CUTADAPT.out.logs.map{ sample_id, log -> file("${log}")}.collect() + COMBINE_CUTADAPT_LOGS_AND_SUMMARIZE(trim_counts, trim_logs) + + // Filter reads + BBDUK(CUTADAPT.out.reads) + filter_counts = BBDUK.out.filter_counts.map{ sample_id, count -> file("${count}")}.collect() + filter_logs = BBDUK.out.logs.map{ sample_id, log -> file("${log}")}.collect() + COMBINE_BBDUK_LOGS_AND_SUMMARIZE(filter_counts, filter_logs) + + FILTERED_FASTQC(BBDUK.out.reads) + filtered_fastqc_files = FILTERED_FASTQC.out.html.flatten().collect() + FILTERED_MULTIQC("filtered", filtered_fastqc_files) + + // Pick OTUs with vsearch + pick_otus(BBDUK.out.reads) + + // Assign taxonomy + RUN_R(pick_otus.out.otus, pick_otus.out.counts, + COMBINE_CUTADAPT_LOGS_AND_SUMMARIZE.out.counts, + COMBINE_BBDUK_LOGS_AND_SUMMARIZE.out.counts) + + // Zip biom file + ZIP_BIOM(RUN_R.out.biom) + + + + // Software Version Capturing - combining all captured sofware versions + RAW_FASTQC.out.version | mix(software_versions_ch) | set{software_versions_ch} + RAW_MULTIQC.out.version | mix(software_versions_ch) | set{software_versions_ch} + CUTADAPT.out.version | mix(software_versions_ch) | set{software_versions_ch} + BBDUK.out.version | mix(software_versions_ch) | set{software_versions_ch} + FILTERED_FASTQC.out.version | mix(software_versions_ch) | set{software_versions_ch} + FILTERED_MULTIQC.out.version | mix(software_versions_ch) | set{software_versions_ch} + pick_otus.out.versions | mix(software_versions_ch) | set{software_versions_ch} + RUN_R.out.version | mix(software_versions_ch) | set{software_versions_ch} + + + nf_version = "Nextflow Version: ".concat("${nextflow.version}\n<><><>\n") + nextflow_version_ch = Channel.value(nf_version) + + // Write software versions to file + software_versions_ch | map { it.text + "\n<><><>\n"} + | unique + | mix(nextflow_version_ch) + | collectFile(name: "${params.metadata_dir}/software_versions.txt", newLine: true, cache: false) + | set{final_software_versions_ch} +} + +workflow.onComplete { + log.info ( workflow.success ? "\nDone! Workflow completed without any error\n" : "Oops .. something went wrong" ) +} diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/assign_taxonomy.nf b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/assign_taxonomy.nf new file mode 100644 index 00000000..22856521 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/assign_taxonomy.nf @@ -0,0 +1,44 @@ +#!/usr/bin/env nextflow +nextflow.enable.dsl = 2 + +process RUN_R { + + tag "Assigning taxonomy to OTUs using decipher..." + + input: + path(otus) + path(counts) + path(trimmed_read_counts) + path(filtered_read_counts) + output: + path("Final_Outputs/${params.output_prefix}taxonomy${params.assay_suffix}.tsv"), emit: taxonomy + path("Final_Outputs/${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.biom"), emit: biom + path("Final_Outputs/${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.tsv"), emit: tsv + path("Final_Outputs/${params.output_prefix}read-count-tracking${params.assay_suffix}.tsv"), emit: read_count + path("versions.txt"), emit: version + script: + """ + mkdir Trimmed_Sequence_Data/ && mv ${trimmed_read_counts} Trimmed_Sequence_Data/ + mkdir Filtered_Sequence_Data/ && mv ${filtered_read_counts} Filtered_Sequence_Data/ + mkdir Final_Outputs/ && \\ + cp ${otus} Final_Outputs/ && \\ + mv ${counts} Final_Outputs/ + + 454-IonTorrent-R-processing.R \ + "${otus}" \ + "Trimmed_Sequence_Data/" \ + "Filtered_Sequence_Data/" \ + "Final_Outputs/" \ + "${params.output_prefix}" \ + "${params.target_region}" \ + "${params.assay_suffix}" + + # Sort the taxonomy count by ASV id + (head -n 1 "Final_Outputs/${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.tsv"; \\ + awk 'NR>1{print}' "Final_Outputs/${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.tsv" | sort -V -k1) \\ + > temp_tax_cont.tsv && mv temp_tax_cont.tsv "Final_Outputs/${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.tsv" + + R --vanilla --version |grep "R version" > versions.txt + get_R_package_version.R + """ +} diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/create_runsheet.nf b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/create_runsheet.nf new file mode 100644 index 00000000..03babf25 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/create_runsheet.nf @@ -0,0 +1,70 @@ +#!/usr/bin/env nextflow +nextflow.enable.dsl = 2 + +//params.GLDS_accession = "OSD-72" +//params.RawFilePattern = null // Pattern of files on OSDR for the OSD accession you want to process + +process GET_RUNSHEET { + + beforeScript "chmod +x ${baseDir}/bin/*" + + input: + val(GLDS_accession) + val(target_region) + output: + tuple path("a_*amplicon*.txt"), path("target_assay_table.txt"), path("runsheet.csv"), emit: tables + path("*.zip"), emit: zip + path("GLfile.csv"), emit: input_file + path("GLparams_file.csv"), emit: params_file + path("versions.txt"), emit: version + script: + """ + # Download ISA zip file for the GLDS_accession then unzip it + GL-download-GLDS-data -g ${GLDS_accession} -p ISA -f && unzip *-ISA.zip + + (head -n1 a_*amplicon*.txt ; \\ + grep "${target_region}" a_*amplicon*.txt) > target_assay_table.txt + + if [ ${params.RawFilePattern} == null ];then + + # Attempt to download the sequences using the assay table, if that fails then + # attempt retrieving all fastq.gz files + GL-download-GLDS-data -f -g ${GLDS_accession} -a target_assay_table.txt -o Raw_Sequence_Data || \\ + GL-download-GLDS-data -f -g ${GLDS_accession} -p ".fastq.gz" -o Raw_Sequence_Data + + else + + + GL-download-GLDS-data -f -g ${GLDS_accession} -p ${params.RawFilePattern} -o Raw_Sequence_Data + + fi + + # Handle case where URLs contain the "+" sign and replaces it with %2B + if grep -q '+' *wanted-file-download-commands.sh;then + grep '+' *wanted-file-download-commands.sh | \\ + sort -u | \\ + awk '{gsub(/\\+/,"%2B", \$NF);print}' \\ + > plus_containing_${GLDS_accession}-wanted-file-download-commands.sh + cat plus_containing_${GLDS_accession}-wanted-file-download-commands.sh | parallel -j $task.cpus + fi + + # Create runsheet, input and parameter files from the target assay table + create_runsheet.sh target_assay_table.txt > runsheet.csv + cut -d "," -f1-2 runsheet.csv > GLfile.csv + cut -d "," -f3- runsheet.csv | uniq > GLparams_file.csv + + GL-version 2>&1 | grep "GeneLab utils"| sed -E 's/^\\s+//' > versions.txt + """ +} + + +workflow { + + GET_RUNSHEET(params.GLDS_accession, params.target_region) + file_ch = GET_RUNSHEET.out.input_file + .splitCsv(header:true) + + params_ch = GET_RUNSHEET.out.params_file + .splitCsv(header:true) + +} diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/quality_assessment.nf b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/quality_assessment.nf new file mode 100644 index 00000000..54085137 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/quality_assessment.nf @@ -0,0 +1,207 @@ +#!/usr/bin/env nextflow +nextflow.enable.dsl = 2 + +/**************************************************************************************** +********************* Sequence quality assessment and control processes **************** +****************************************************************************************/ + +// A 2-column (single-end) or 3-column (paired-end) file +//params.csv_file = "${baseDir}/file.csv" +//params.prefix = "raw" + +// FastQC performed on reads +process FASTQC { + + tag "Running fastqc on ${sample_id}..." + beforeScript "chmod +x ${baseDir}/bin/*" + label "fastqc" + + input: + tuple val(sample_id), path(reads) + output: + tuple path("*.html"), path("*.zip"), emit: html + path("versions.txt"), emit: version + script: + """ + fastqc -o . \\ + -t ${task.cpus} -q \\ + ${reads} + + fastqc --version > versions.txt + """ +} + + +process MULTIQC { + + tag "Running multiqc on the ${prefix} files..." + beforeScript "chmod +x ${baseDir}/bin/*" + + input: + val(prefix) + path(files) + output: + path("${params.output_prefix}${prefix}_multiqc${params.assay_suffix}_data.zip"), emit: data + path("${params.output_prefix}${prefix}_multiqc${params.assay_suffix}_report.html"), emit: html + path("versions.txt"), emit: version + script: + """ + multiqc -z -q -o . \\ + -n "${params.output_prefix}${prefix}_multiqc${params.assay_suffix}" . \\ + > /dev/null 2>&1 + + + # Renaming html file + mv ${params.output_prefix}${prefix}_multiqc${params.assay_suffix}.html \\ + ${params.output_prefix}${prefix}_multiqc${params.assay_suffix}_report.html + + multiqc --version > versions.txt + """ + } + + + + +// This process runs cutadapt +process CUTADAPT { + + tag "Trimming off primers for ${sample_id} using cutadapt..." + beforeScript "chmod +x ${baseDir}/bin/*" + + input: + tuple val(sample_id), path(reads) + tuple val(F_primer), val(R_primer) + output: + tuple val(sample_id), path("${sample_id}${params.primer_trimmed_suffix}"), emit: reads + tuple val(sample_id), path("${sample_id}-cutadapt.log"), emit: logs + tuple val(sample_id), path("${sample_id}-trimmed-counts.tsv"), emit: trim_counts + path("versions.txt"), emit: version + script: + """ + cutadapt -g ${F_primer} \\ + -a ${R_primer} \\ + -o ${sample_id}${params.primer_trimmed_suffix} \\ + ${reads[0]} > ${sample_id}-cutadapt.log 2>&1 + + paste <( printf "${sample_id}" ) \\ + <( grep "Total reads processed:" ${sample_id}-cutadapt.log | tr -s " " "\\t" | cut -f 4 | tr -d "," ) \\ + <( grep "Reads written (passing filters):" ${sample_id}-cutadapt.log | tr -s " " "\\t" | cut -f 5 | tr -d "," ) \\ + > ${sample_id}-trimmed-counts.tsv + + VERSION=`cutadapt --version` + echo "cutadapt \${VERSION}" > versions.txt + """ +} + + + + +// This process combines the cutadapt logs and summarizes them. +process COMBINE_CUTADAPT_LOGS_AND_SUMMARIZE { + + tag "Combining the logs generated by cutadapt..." + + input: + path(counts) + path(logs) + output: + path("${params.output_prefix}cutadapt${params.assay_suffix}.log"), emit: logs + path("${params.output_prefix}trimmed-read-counts${params.assay_suffix}.tsv"), emit: counts + script: + """ + cat ${logs} > ${params.output_prefix}cutadapt${params.assay_suffix}.log + + cat <( printf "sample\\traw_reads\\tcutadapt_trimmed\\n" ) \\ + <( cat ${counts} ) > ${params.output_prefix}trimmed-read-counts${params.assay_suffix}.tsv + """ +} + + + +// This process runs quality filtering/trimming on input fastq files. +process BBDUK { + + + tag "Quality filtering ${sample_id}s reads.." + beforeScript "chmod +x ${baseDir}/bin/*" + + input: + tuple val(sample_id), path(reads) + output: + tuple val(sample_id), path("${sample_id}${params.filtered_suffix}"), emit: reads + tuple val(sample_id), path("${sample_id}-bbduk.log"), emit: logs + tuple val(sample_id), path("${sample_id}-filtered-counts.tsv"), emit: filter_counts + path("versions.txt"), emit: version + script: + """ + bbduk.sh in=${reads[0]} out1=${sample_id}${params.filtered_suffix} \\ + qtrim=r trimq=10 mlf=0.5 \\ + minavgquality=${params.min_bbduk_avg_quality} \\ + minlength=${params.min_bbduk_len} \\ + > ${sample_id}-bbduk.log 2>&1 + + paste <( printf "${sample_id}" ) <( grep "Input:" ${sample_id}-bbduk.log | \\ + tr -s " " "\\t" | cut -f 2 ) <( grep "Result:" ${sample_id}-bbduk.log | \\ + tr -s " " "\\t" | cut -f 2 ) > ${sample_id}-filtered-counts.tsv + + VERSION=`bbversion.sh` + echo "bbtools \${VERSION}" > versions.txt + """ +} + + + +// This process combines the bbduk logs and summarizes them. +process COMBINE_BBDUK_LOGS_AND_SUMMARIZE { + + tag "Combining the logs generated by bbduk..." + + input: + path(counts) + path(logs) + output: + path("${params.output_prefix}bbduk${params.assay_suffix}.log"), emit: logs + path("${params.output_prefix}filtered-read-counts${params.assay_suffix}.tsv"), emit: counts + script: + """ + cat ${logs} > ${params.output_prefix}bbduk${params.assay_suffix}.log + + cat <( printf "sample\\tinput_reads\\tfiltered_reads\\n" ) \\ + <( cat ${counts} ) > ${params.output_prefix}filtered-read-counts${params.assay_suffix}.tsv + """ +} + + + +workflow quality_check { + + take: + prefix_ch + multiqc_config + reads_ch + + + main: + FASTQC(reads_ch) + fastqc_ch = FASTQC.out.html.flatten().collect() + MULTIQC(prefix_ch, multiqc_config, fastqc_ch) + + software_versions_ch = Channel.empty() + FASTQC.out.version | mix(software_versions_ch) | set{software_versions_ch} + MULTIQC.out.version | mix(software_versions_ch) | set{software_versions_ch} + + emit: + versions = software_versions_ch +} + +workflow { + + Channel.fromPath(params.csv_file) + .splitCsv(header:true) + .map{row -> tuple( "${row.sample_id}", [file("${row.read}", checkIfExists: true)] )} + .set{reads_ch} + + + quality_check(Channel.of(params.prefix), params.multiqc_config, reads_ch) + CUTADAPT(reads_ch) +} diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/vsearch.nf b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/vsearch.nf new file mode 100644 index 00000000..e3abcac8 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/vsearch.nf @@ -0,0 +1,139 @@ +#!/usr/bin/env nextflow +nextflow.enable.dsl = 2 + +process VSEARCH_DEREP_SAMPLE { + + tag "dereplicating ${sample_id}s sequences..." + label "vsearch" + + input: + tuple val(sample_id), path(reads) + output: + path("${sample_id}-derep.fa.tmp"), emit: derep_temp + path("versions.txt"), emit: version + script: + """ + vsearch --derep_fulllength ${reads} \\ + --strand both \\ + --output "${sample_id}-derep.fa.tmp" \\ + --sizeout \\ + --relabel "sample=${sample_id};seq_" > /dev/null + + vsearch --version 2>&1 |grep "vsearch" |head -n1 | sed -E 's/(vsearch v.+?)_linux.+/\\1/' > versions.txt + """ +} + + +process VSEARCH_COMBINE_DEREPD_SAMPLES { + + tag "Combining all dereplicated samples..." + label "vsearch" + + input: + path(derepd_reads) + output: + path("all-samples.fa.tmp") + script: + """ + cat ${derepd_reads} > all-samples.fa.tmp + """ +} + + process VSEARCH_PROCESS_ALL { + + tag "Clustering your sequences to OTUs using vsearch..." + label "vsearch" + + input: + path(all_samples_fasta) + output: + path("${params.output_prefix}OTUs.fasta"), emit: fasta + path("${params.output_prefix}counts${params.assay_suffix}.tsv"), emit: counts + path("versions.txt"), emit: version + script: + """ + # Dereplicate all + vsearch --derep_fulllength ${all_samples_fasta} \\ + --strand both \\ + --output all-samples_derep.fa.tmp \\ + --sizein --sizeout + + # Clustering to get rep seqs + vsearch --cluster_size all-samples_derep.fa.tmp \\ + --id 0.97 \\ + --strand both \\ + --sizein \\ + --sizeout \\ + --relabel "OTU_" \\ + --centroids rep-seqs.fa.tmp + + # Removing singletons + vsearch --sortbysize rep-seqs.fa.tmp \\ + --minsize 2 \\ + --output rep-seqs-no-singletons.fa.tmp + + # Chimera check and removal + vsearch --uchime_denovo rep-seqs-no-singletons.fa.tmp \\ + --sizein \\ + --nonchimeras ${params.output_prefix}OTUs.fasta \\ + --relabel "OTU_" + + # Mapping seqs to OTUs to get OTU abundances per sample + vsearch --usearch_global ${all_samples_fasta} \\ + -db ${params.output_prefix}OTUs.fasta \\ + --sizein --id 0.97 \\ + --otutabout counts.tmp + + sed 's/^#OTU ID/OTU_ID/' counts.tmp \\ + > ${params.output_prefix}counts${params.assay_suffix}.tsv + + vsearch --version 2>&1 |grep "vsearch" |head -n1 | sed -E 's/(vsearch v.+?)_linux.+/\\1/' > versions.txt + """ + } + + +process REMOVE_LINE_WRAPS { + + tag "Removing line wraps from OTU fasta file..." + + input: + path(temp_fasta) + output: + path("${params.output_prefix}OTUs${params.assay_suffix}.fasta"), emit: fasta + path("versions.txt"), emit: version + script: + """ + # Removing line wraps from fasta file + bit-remove-wraps ${temp_fasta} \\ + > ${params.output_prefix}OTUs${params.assay_suffix}.fasta.tmp && \\ + mv ${params.output_prefix}OTUs${params.assay_suffix}.fasta.tmp \\ + ${params.output_prefix}OTUs${params.assay_suffix}.fasta + + bit-version 2>&1 |grep "Bioinformatics Tools"|sed -E 's/^\\s+//' > versions.txt + """ +} + + +workflow pick_otus { + + take: + reads_ch + + main: + VSEARCH_DEREP_SAMPLE(reads_ch) + derep_temp_ch = VSEARCH_DEREP_SAMPLE.out.derep_temp.collect() + VSEARCH_COMBINE_DEREPD_SAMPLES(derep_temp_ch) | VSEARCH_PROCESS_ALL + REMOVE_LINE_WRAPS(VSEARCH_PROCESS_ALL.out.fasta) + + // capture software versions + software_versions_ch = Channel.empty() + VSEARCH_DEREP_SAMPLE.out.version | mix(software_versions_ch) | set{software_versions_ch} + VSEARCH_PROCESS_ALL.out.version | mix(software_versions_ch) | set{software_versions_ch} + REMOVE_LINE_WRAPS.out.version | mix(software_versions_ch) | set{software_versions_ch} + + emit: + otus = REMOVE_LINE_WRAPS.out.fasta + counts = VSEARCH_PROCESS_ALL.out.counts + versions = software_versions_ch + +} diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/zip_biom.nf b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/zip_biom.nf new file mode 100644 index 00000000..e1da0fa0 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/modules/zip_biom.nf @@ -0,0 +1,25 @@ +#!/usr/bin/env nextflow +nextflow.enable.dsl = 2 + + +/**************************************************************************************** +********************* Zip Biom ******************************************************** +****************************************************************************************/ + +process ZIP_BIOM { + + tag "Zipping the taxonomy counts...." + + input: + path(taxonomy_and_counts_biom) + output: + path("${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.biom.zip"), emit: biom + path("versions.txt"), emit: version + script: + """ + zip -q ${params.output_prefix}taxonomy-and-counts${params.assay_suffix}.biom.zip \\ + ${taxonomy_and_counts_biom} + + zip -h | grep "Zip" | sed -E 's/(Zip.+\\)).+/\\1/' > versions.txt + """ +} diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/nextflow.config b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/nextflow.config new file mode 100644 index 00000000..040d66a3 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/nextflow.config @@ -0,0 +1,231 @@ +//******** Global parameters *****************// +params { + + //---------- Required parameters -----------------------------// + csv_file = "file.csv" // A 2-column input file ["sample_id", "read"] + publishDir_mode = "link" // "link" , "copy" + // Cutadapt parameters + F_primer = "" + R_primer = "" + + // BBDUK parameters + min_bbduk_len = 50 // minimum length threshold for bbduk + min_bbduk_avg_quality = 15 // bbduk minimum average quality + + target_region = "16S" // options are ["16S", "ITS"] + + //----------------------- Optional parameters ------------------------------------// + // Suffixes + primer_trimmed_suffix = "_trimmed.fastq.gz" + filtered_suffix = "_filtered.fastq.gz" + // Directories + raw_reads_dir = "../Raw_Sequence_Data/" + fastqc_out_dir = "../FastQC_Outputs/" + trimmed_reads_dir = "../Trimmed_Sequence_Data/" + filtered_reads_dir = "../Filtered_Sequence_Data/" + final_outputs_dir = "../Final_Outputs/" + metadata_dir = "../Metadata/" + genelab_dir = "../GeneLab/" + + + + // Specify paths to existing conda environments (/path/to/envs/qc) else leave as is so that + // new conda environements will be created if `-profile conda` is used to run the pipeline + conda{ + qc = null // /path/to/envs/qc + R = null // /path/to/envs/R + bbmap = null // /path/to/envs/bbmap + cutadapt = null // /path/to/envs/cutadapt + vsearch = null // /path/to/envs/vsearch + genelab = null // /path/to/envs/genelab + } + + // Genelab specific parameters + assay_suffix = "_GLAmpSeq" + output_prefix = "" + GLDS_accession = false // OSD acession number for the data to be processed + // Pattern of files on OSDR for the GLDS_accession you want to process. + RawFilePattern = null // e.g. "_Amplicon_" + + + errorStrategy = "terminate" // how should errors be handled by nextflow. + debug = false // should a detailed log of set parameters be shown before workflow execution? +} + + +// Setting the default container engine as singularity +params.containerEngine = "singularity" +// Conda shouldn't be used by default except when using conda-based profiles +params.use_conda = false + +/******************************************************************************************************* +*************************************** Workflow Profiles ********************************************** +********************************************************************************************************/ + +profiles { + + slurm { + process.executor = 'slurm' + } + + conda { + conda.enabled = true + params.use_conda = true + } + + singularity { + singularity.enabled = true + singularity.autoMounts = true + singularity.cacheDir = "singularity/" + params.containerEngine = "singularity" + } + + docker { + docker.enabled = true + docker.runOptions = '-u $(id -u):$(id -g)' + docker.userEmulation = true + params.containerEngine = "docker" + } + +} + +// Number of jobs to run in parallel +executor.queueSize = 20 + +/****************************************************************************************************************** +***************** Tune process specific resources (cpu, container, memory etc.) *********************************** +*******************************************************************************************************************/ + +process { + + //******************* Default process settings ************************// + + // "ignore" will ignore errors while "retry" will retry the failed task as many times as specified by maxRetries below + errorStrategy = { params.errorStrategy ? params.errorStrategy : "ignore"} + //maxRetries = 2 // uncomment if you'd like to retry a failed task. + cpus = 2 + memory = '5 GB' + cache = 'lenient' + //debug = true // uncomment to see what is being emitted to the standard output + + /********************************************************************************************* + ******************************** Process Specific Settings ********************************** + *********************************************************************************************/ + +//************************* GLDS_accession runsheet and input file retrieval **************************************// + + withName: GET_RUNSHEET { + conda = {params.conda.genelab != null ? params.conda.genelab : "envs/genelab.yaml"} + container = "olabiyi/genelab-utils:1.3.22" + publishDir = [path: params.genelab_dir, mode: params.publishDir_mode] + } + +//********************************** Read quality control and assesment ********************************************// + withLabel: fastqc { + conda = {params.conda.qc != null ? params.conda.qc : "envs/qc.yaml"} + container = "staphb/fastqc:0.12.1" + } + + withName: RAW_FASTQC { + publishDir = [path: params.raw_reads_dir, mode: params.publishDir_mode] + } + + withName: "RAW_MULTIQC|FILTERED_MULTIQC" { + conda = {params.conda.qc != null ? params.conda.qc : "envs/qc.yaml"} + container = "staphb/multiqc:1.19" + publishDir = [path: params.fastqc_out_dir, mode: params.publishDir_mode] + } + + withName: "CUTADAPT|COMBINE_CUTADAPT_LOGS_AND_SUMMARIZE" { + conda = {params.conda.cutadapt != null ? params.conda.cutadapt : "envs/cutadapt.yaml"} + container = "zavolab/cutadapt:1.16" + publishDir = [path: params.trimmed_reads_dir, mode: params.publishDir_mode] + } + + withName: "BBDUK|COMBINE_BBDUK_LOGS_AND_SUMMARIZE" { + conda = {params.conda.bbmap != null ? params.conda.bbmap : "envs/bbmap.yaml"} + container = "staphb/bbtools:38.86" + memory = "20 GB" + cpus = 5 + publishDir = [path: params.filtered_reads_dir, mode: params.publishDir_mode] + } + + withName: FILTERED_FASTQC { + publishDir = [path: params.filtered_reads_dir, mode: params.publishDir_mode ] + } + +//******************************** OTU Picking with VSEARCH **************************************// + withLabel: vsearch { + conda = {params.conda.vsearch != null ? params.conda.vsearch : "envs/vsearch.yaml"} + container = "quay.io/biocontainers/vsearch:2.15.2--h2d02072_0" + memory = '100 GB' + cpus = 10 + } + + withName: VSEARCH_PROCESS_ALL { + publishDir = [path: params.final_outputs_dir, + pattern: "${params.output_prefix}counts${params.assay_suffix}.tsv", + mode: params.publishDir_mode] + } + + withName: REMOVE_LINE_WRAPS { + conda = {params.conda.vsearch != null ? params.conda.vsearch : "envs/vsearch.yaml"} + container = "olabiyi/bit-astrobiomike:1.0" + memory = "5 GB" + cpus = 2 + publishDir = [path: params.final_outputs_dir, mode: params.publishDir_mode] + } + + withName: RUN_R { + conda = {params.conda.R != null ? params.conda.R : "envs/R.yaml"} + container = "olabiyi/r-dada-decipher-biomformat:1.0" + memory = "100 GB" + cpus = 10 + publishDir = [path: params.final_outputs_dir , pattern: "Final_Outputs/*.{tsv,biom}", + mode: params.publishDir_mode, saveAs: { fn -> fn.substring(fn.lastIndexOf('/')+1)} ] + } + + withName: ZIP_BIOM { + conda = {params.conda.qc != null ? params.conda.qc : "envs/qc.yaml"} + container = "staphb/multiqc:1.19" + publishDir = [path: params.final_outputs_dir, mode: params.publishDir_mode] + } + +} + + +/***************************************************************************** +********************** Workflow Resource Usage Capturing ********************* +******************************************************************************/ + +// Adapted from : https://github.com/nf-core/rnaseq/blob/master/nextflow.config +def trace_timestamp = new java.util.Date().format( 'yyyy-MM-dd_HH-mm-ss') +timeline { + enabled = true + file = "../Resource_Usage/execution_timeline_${trace_timestamp}.html" +} +report { + enabled = true + file = "../Resource_Usage/execution_report_${trace_timestamp}.html" +} +trace { + enabled = true + file = "../Resource_Usage/execution_trace_${trace_timestamp}.txt" +} + + +/****************************************************************************** +**************************** Workflow Metadata ******************************** +*******************************************************************************/ + +manifest { + author = 'Olabiyi Aderemi Obayomi, Mike Douglas Lee' + homePage = 'https://github.com/nasa/GeneLab_Data_Processing/blob/master/Amplicon/' + description = '454 Ion Torrent workflow for pipeline document GL-DPPD-7106' + mainScript = 'main.nf' + defaultBranch = 'main' + nextflowVersion = '>=22.10.1' + version = '1.0.0' +} + + diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/slurm_submit.slurm b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/slurm_submit.slurm new file mode 100644 index 00000000..fcddf897 --- /dev/null +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/NF_Amp454IonTor/workflow_code/slurm_submit.slurm @@ -0,0 +1,64 @@ +#!/bin/bash + +#SBATCH --job-name="nf_master" ## Replace job_name with the name of the job you are running ## +#SBATCH --output=nf_master.o.%j ## Replace job_name with the name of the job you are running ## +#SBATCH --error=nf_master.e.%j ## Replace job_name with the name of the job you are running ## +#SBATCH --partition=normal ## Specifies the job queue to use, for urgent jobs change normal to priority ## +#SBATCH --mem=2G ## Memory required to run the job in MB, this example is showing 10,000 MB or 10GB, change this number based on how much RAM you need ## +#SBATCH --cpus-per-task=1 ## Number of CPUs to run the job, this example is showing 5 CPUs, change this number based on how many CPUs you need ## +#SBATCH --mail-user=email@domain.com ## Specifies the e-mail address to e-mail when the job is complete, replace this e-mail address with your NASA e-mail address ## +#SBATCH --mail-type=END ## Tells slurm to e-mail the address above when the job has completed ## + +. ~/.profile + + +echo "nf_master" ## Replace job_name with the name of the job you are running ## +echo "" + + +## Add a time-stamp at the start of the job ## +start=$(date +%s) +echo "start time: $start" + +## Print the name of the compute node executing the job ## +echo $HOSTNAME + + +## Activate the conda environemnt containing the tools you need to run your job ## +## You can see a list of all available environments by running the command: conda env list ## +## If you need a conda envrionment installed request it using JIRA ## + +source activate /path/to/envs/nextflow ## Replace conda_env_name with the name of the environment with nextflow installed ## + + +## Print the version of the tool you are using to ensure the tool version is recorded ## +echo "" +echo "Nextflow version: " ## Replace Tool with the name of the tool you are using ## +nextflow -v ## Replace this command with the command the tool uses to print its version ## +echo "" + + +## The command(s) that you want to run in this slurm job ## +export NXF_SINGULARITY_CACHEDIR=singularity/ +## Replace command with the command(s) you want to run ## +#nextflow run main.nf -profile slurm,singularity -resume --csv_file file.csv --target_region 16S --F_primer AGAGTTTGATCCTGGCTCAG --R_primer CTGCCTCCCGTAGGAGT --min_bbduk_len 50 --min_bbduk_avg_quality 15 +nextflow run main.nf -profile slurm,singularity -resume --GLDS_accession OSD-72 --target_region 16S --min_bbduk_len 50 --min_bbduk_avg_quality 15 + +## Add a time-stamp at the end of the job then calculate how long the job took to run in seconds, minutes, and hours ## +echo "" +end=$(date +%s) +echo "end time: $end" +runtime_s=$(echo $(( end - start ))) +echo "total run time(s): $runtime_s" +sec_per_min=60 +sec_per_hr=3600 +runtime_m=$(echo "scale=2; $runtime_s / $sec_per_min;" | bc) +echo "total run time(m): $runtime_m" +runtime_h=$(echo "scale=2; $runtime_s / $sec_per_hr;" | bc) +echo "total run time(h): $runtime_h" +echo "" + + +## Print the slurm job ID so you have it recorded and can view slurm job statistics if needed ## +echo "slurm job ID: ${SLURM_JOB_ID}" + diff --git a/Amplicon/454-and-IonTorrent/Workflow_Documentation/README.md b/Amplicon/454-and-IonTorrent/Workflow_Documentation/README.md index 3c46d116..a5253a52 100644 --- a/Amplicon/454-and-IonTorrent/Workflow_Documentation/README.md +++ b/Amplicon/454-and-IonTorrent/Workflow_Documentation/README.md @@ -6,8 +6,8 @@ |Pipeline Version|Current Workflow Version (for respective pipeline version)| |:---------------|:---------------------------------------------------------| -|*[GL-DPPD-7106.md](../Pipeline_GL-DPPD-7106_Versions/GL-DPPD-7106.md)|[1.0.0](SW_Amp454IonTor)| +|*[GL-DPPD-7106.md](../Pipeline_GL-DPPD-7106_Versions/GL-DPPD-7106.md)|[1.0.0](NF_Amp454IonTor)| *Current GeneLab Pipeline/Workflow Implementation -> See the [workflow change log](SW_Amp454IonTor/CHANGELOG.md) to access previous workflow versions and view all changes associated with each version update. +> See the [workflow change log](NF_Amp454IonTor/CHANGELOG.md) to access previous workflow versions and view all changes associated with each version update.