From 42c7eb67dffdaa62137f76ff4064e29a345b48a4 Mon Sep 17 00:00:00 2001 From: bioaddict Date: Thu, 2 Nov 2023 14:14:03 +0100 Subject: [PATCH] Rename the pipeline and add the README file --- README.md | 57 +++++++++++++++++- fast5pipeline | 164 ++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 220 insertions(+), 1 deletion(-) create mode 100755 fast5pipeline diff --git a/README.md b/README.md index 6b10510..1f9ad7d 100644 --- a/README.md +++ b/README.md @@ -1 +1,56 @@ -# FAST5toGenome \ No newline at end of file + +# Fast5pipeline : Bacterial Isolate Genome Sequencing Data Processing pipeline + +fast5pipeline is designed to process bacterial isolate whole genome sequencing data obtained using Nanopore sequencing technology (ONT). The pipeline takes fast5 files as input and includes basecalling, demultiplexing, assembling, and polishing of the genome sequences. This pipeline has been design and tested on ubuntu 22.04. You may face some issues while trying to run it in other linux distributions. + +``` +OPTIONS: + +-c: Path to the Guppy configuration file. +-i: Path to the directory containing the raw nanopore sequencing data (fast5 files) +-o: Path to the directory where processed data will be stored. This directory will be created if it doesn't exist. +-k: Barcode kits to be used for barcoding. +-m: medaka model for consensus polishing. +-b: GPU memory control for medaka (optional) +-v: Enable verbosity. +``` + +## Instructions to running fast5pipeline +### Setup environment +Before using this pipeline, you need to setup the required environment by installing the following software, if they are not already installed: + +- [Guppy Basecalling Software](https://community.nanoporetech.com/docs/prepare/library_prep_protocols/Guppy-protocol/v/gpb_2003_v1_revax_14dec2018/linux-guppy), (C) Oxford Nanopore Technologies plc. Version 6.5.7+ca6d6af, minimap2 version 2.24-r1122 or above +- [GNU gzip](https://www.gnu.org/software/gzip/) +- [flye 2.9.2-b1786](https://github.com/fenderglass/Flye) or above +- [conda](https://conda.io/projects/conda/en/latest/user-guide/install/index.html) or [mamba](https://mamba.readthedocs.io/en/latest/mamba-installation.html#mamba-install) +- [medaka 1.8.0 or above](https://anaconda.org/bioconda/medaka) (To be installed with conda or mamba in an environment named medaka) + +### Run fast5pipeline + +To run this pipeline, kindly follow this instructions: + +```bash +git clone https://github.com/AlbertRockG/fast5toGenome +cd fast5Genome/ +export PATH="$PWD/fast5Genome:$PATH" + +# Example of running +./fast5pipeline -c dna_r9.4.1_450bps_hac.cfg \\ + -i data/IN_DIR \\ + -o data/OUT_DIR \\ + -k SQK-RBK110-96 \\ + -m r941_min_hac_g507 \\ + -v +``` + +## Acknowledgments + +I would like to express my gratitude to the following individuals and organizations for their valuable contributions and support: + +- Marco, a senior bioinformatician, for providing guidance and assistance throughout the process. +- FAO (Food and Agriculture Organization), SeqAfrica, and FlemingFund for their support in funding and resources. +- The Noguchi Memorial Institute for Medical Research for providing the infrastructure and facilities for this project. + +## License + +This project is licensed under the MIT License - see the [LICENSE](LICENSE) file for details. diff --git a/fast5pipeline b/fast5pipeline new file mode 100755 index 0000000..aee8038 --- /dev/null +++ b/fast5pipeline @@ -0,0 +1,164 @@ +#!/usr/bin/env bash +# -*- coding: utf-8 -*- +# Created on Fri June 23 14:31:04 2023 +# @author: AlbertRockG + +export LC_ALL="C" +set -euop pipefail +### FUNCTIONS -------------------------------------------------------------------------------------------------------------------------- +# Function to display usage information +usage() { + printf ' +Usage: %s [OPTIONS] + + Basecalls, trims, joins, assembles and polish bacteria isolates genome sequencing + data. Take fast5 files as inputs. + + OPTIONS + -c Guppy config file + -i Input directory + -o Output directory: different from the input directory + -k Barcode kits + -m Medaka model: only high accuracy models supported + -b GPU memory control: controls memory use (default: 100) + -v Enable verbosity + ' "$(basename "$0")" && + exit 1 +} + +err_msg() { + echo "${0##*/}: $*" >&2; +} + +emit() { + (( !VERBOSE )) || err_msg "$*"; +} + +err_exit() { + err_msg "$*"; + exit 1; +} + +# Function to concatenate fastq.gz files in a subfolder +concatenate_fastq() { + local subfolder="$1" + + emit "Processing subfolder: $subfolder" + + # Create the "joined" folder in the parent directory + parent_dir="${subfolder%/*/*}" + joined_folder="$parent_dir/joined" + mkdir -p "$joined_folder" + + # Create output file name based on subfolder name + output_file="$joined_folder/$(basename "$subfolder").fastq.gz" + + # Concatenate fastq.gz files into the output file + gunzip -c "$subfolder"/*.fastq.gz | gzip -c > "$output_file" + + emit "Concatenated fastq.gz files into: $output_file" + emit "" +} + +# Function to run flye on a joined_fastq file +run_flye() { + local input_file="$1" + + emit "Processing file: $input_file" + + # Create the "assembled" folder in the parent directory + parent_dir="${input_file%/*/*}" + assembly_folder="$parent_dir/assembled" + mkdir -p "$assembly_folder" + + # Create output file name based on subfolder name + output="$assembly_folder/$(basename "$input_file")" + + # Run flye on the input file + flye -t "$(nproc)" --out-dir "$output" --nano-hq "$input_file" + + emit "Flye analysis completed for: $input_file" + emit "" +} + +# Function to run flye on a joined_fastq file +run_medaka() { + local input_file="$1" + + emit "Processing file: $input_file" + + # Create the "polished" folder in the parent directory + main_dir="${input_file%/*/*}" + assembly_dir="$(basename "$input_file")" + polished_folder="$main_dir/polished" + mkdir -p "$polished_folder" + + # Create output dir name based on the input file name + output="$polished_folder/$(basename "$input_file")" + + # Run medaka on the input file + medaka_consensus -i "$input_file" \ + -d "$main_dir"/assembled/"$assembly_dir"/assembly.fasta \ + -o "$output" -m "$medaka_model" -t "$(nproc)" -b "$BATCHSIZE" + + emit "Medaka analysis completed for: $input_file" + emit "" +} + +### VARIABLE DECLARATION --------------------------------------------------------------------------------------------------------------- + +VERBOSE=0 +BATCHSIZE=100 + +### ARGS PARSING ----------------------------------------------------------------------------------------------------------------------- +# Parse command-line options + +while getopts "c:i:o:k:m:bv" opt; do + case "${opt}" in + c) guppy_config=${OPTARG};; + i) input_dir=${OPTARG};; + o) output_dir=${OPTARG};; + k) barcode_kits=${OPTARG};; + m) medaka_model=${OPTARG};; + b) BATCHSIZE=${OPTARG};; + v) VERBOSE=1;; + *) usage;; + esac +done + +# Check if all required options are provided + +if [[ -z "${guppy_config:-}" || -z "${input_dir:-}" || -z "${output_dir:-}" || -z "${barcode_kits:-}" || -z "${medaka_model:-}" ]]; then + usage +fi + +### DATA PROCESSING -------------------------------------------------------------------------------------------------------------------- + +# Step 1: Run guppy_basecaller +guppy_basecaller --disable_pings -x auto \ + -c "$guppy_config" \ + -i "$input_dir" --recursive \ + -s "$output_dir/basecalled" + +# Step 2: Run guppy_barcoder +guppy_barcoder -t "$(nproc)" --disable_pings -x auto \ + --barcode_kits "$barcode_kits" --enable_trim_barcodes \ + -i "$output_dir/basecalled" --recursive \ + -s "$output_dir/trimmed" --compress_fastq + +# Step 3: Join fastq.gz files + +# Loop through subfolders in the input folder +for subfolder in "$output_dir"/trimmed/*; do + [[ ! -d "$subfolder" ]] || concatenate_fastq "$subfolder" || emit "no such directory: $subfolder" +done + +# Step 4: Run flye +for input_file in "$output_dir"/joined/*.fastq.gz; do + [[ ! -f $input_file ]] || run_flye "$input_file" || emit "no such file: $input_file" +done + +# Step 6: Run medaka_consensus +for input_file in "$output_dir"/joined/*.fastq.gz; do + [[ ! -f $input_file ]] || run_medaka "$input_file" || emit "no such file: $input_file" +done