trust1/ubuntu

•Updated over 5 years ago

Target Capture Cloud Platform (TCCP) Oz Mammals Genomics Phylogenomics aka OMG pipeline

Image

266

trust1/ubuntu repository overview

⁠Target Capture Cloud Platform (TCCP aka OMG Pipeline) Brief Overview

The Target Capture Cloud Platform is a cloud computing solution for Phylogenetics with targeted capture DNA sequencing protocols (e.g. exon capture). It comprises several related pipelines and additional software tools, and is being built by the ANU Bioinformatics Consultancy in consultation with researchers across Australia. The requirements were specified by the Oz Mammals Genomics (OMG) Project, and funded by Professor Craig Moritz (ANU) and Bio Platforms Australia (BPA).

⁠Credits (in chronological order):

Dr Jason Bragg
Dr Cameron Jack
Dr Zhi-Ping Feng
Bob Buckley
Dr Rust Turakulov :: roustam_t[at]yahoo.com
Dr Terry Bertozzi Senior Research Scientist - Bioinformatics // South Australian Museum //

==========================================================================

This docker image is attempt to package all OMG pipeline tools into one container. Please be aware not all of the scripts from original Ansible distribution by ABC group had been tested. Work is in progress. Detailed protocol ⁠ for the quick test run on inbuilt small data-set run.

⁠Instructions for the new user.

First three steps may required In case you do not have Docker installed locally. If your system already have docker proceed straight to command 4): to generate new container on the system.

Create docker account: www.docker.com⁠ (optional)
Install docker on your computer. To check if it is already there, type: .

   docker -v

Run command from terminal or command prompt [windows]: (optional)

   docker login

Start OMG image. In case of the first time start, the image will not be available locally, however it will be downloaded automatically from docker repository. If you have old image stored you may want to check "docker pull" command to update it. Please replace directory path 'Z:/dockerdata' to your local directory path in following command in terminal (command prompt):

   docker run -v Z:/dockerdata:/mnt/OMG  --name YourLbl4Container   -it trust1/ubuntu:OMGv001

Input path: Z:/dockerdata and container name: Your Lbl4Container are the only variable in this command to adjust to your local system. Now you will be inside container which is indicated by changed location of cursor path like this: root@f283a9697c1e:/#
Container is based on ubuntu v18. You are welcome to use available bash tools inside. Modification are not saved for the next use of the image unless you commit changes. See the docker documentation how to save your changes within container.

Start core OMG pipeline Target Capture De-novo Phylogenetic Pipeline with the local FASTQ files.
1. Create a new working folder, the name must match your sample name. Inside that folder create sub folder named 'raw' for example
 
 Outside docker on host system: Z:/dockerdata/Sample_123/raw
 
 Inside the docker container path /mnt/OMG/Sample_123/raw
2. The fastq files placed in 'raw' folder must have names with the pattern: /mnt/OMG//raw/_L00[1-8]_R[12]_0001.fastq.gz
for example paired read sequencing data for sample ' Sample_123 ':
```
 /mnt/OMG/Sample_123/raw/Sample_123_L001_R1_0001.fastq.gz

 /mnt/OMG/Sample_123/raw/Sample_123_L001_R2_0001.fastq.gz
```
Data file names must end in _L001_R1_0001.fastq.gz (they've been compressed with gzip)
1. To submit sample files for processing, use the following command (inside the docker container):

      nohup /home/ubuntu/omg/code/run-omgpipeline.sh Sample_123 qtest >nohup_OMG_pipeline.txt&

where the 'qtest' options could be replaced to 'marsupial', 'bat', or 'skink'. The default enrichment kit is 'marsupial'. The 'qtest' kit has only nine exons from the bat targets. This feature is for the quick pipeline run and diagnostics purposes. The 'qtest' targets are better be used with small training sample set. Training datasets with 3 samples saved within container at directory: ~/omg/data/quicktest move those fastq files to /mnt/OMG if you want to run diagnostic or training on them.

You can monitor progress inside docker container with the following commands:

        tail -f /mnt/OMG/Sample_123/batch.*

6) Building the phylogenetic tree from multiple samples Once several samples have been processed, they can be used to build a phylogenetic tree. The script for doing this is run-mktree.sh: which will run muscle alignment of assembled fasta files for all samples it can find in /mnt/OMG/ folder and then generate iqtree. This tree file can be taken elsewhere for plotting and adding graphical decoration.

    ~/omg/code/run-mktree.sh

The exon data for each processed sample is collected and the results stored in the mktree directory. For docker system if you run data in multiple containers most likely you would have to combine all data together and make it accessible from one of the containers. Avoid to use soft links for raw fastq of pipeline generated files it may not work inside different container.

==========================================================================

⁠ Run batch of samples

The primary source for raw sequencing data for this pipeline is Bio Platforms Australia (BPA) Data Portal. Here is instructions⁠ how to find and download bulk of raw reads to your local drive from BPA repo. You can always use your own data but make sure your fastq name system matched BPA (see above). At the host system copy batch of conventionally named sequence files ( *_L01_R1_0001.fastq.gz ) to the folder you wish to attach to docker container ( Z:/dockerdata ). You can have samples split across multiple files for the forward or reverse read. For example sample can be sequenced across multiple lines and produce line split fastq files. You do not have to combine such files to single one for each read.
Fire docker image (from the host system):

        docker run -v Z:/dockerdata:/mnt/OMG -it trust1/ubuntu:OMGv001

Run perl preparation script (inside container).

      perl ~/prep_omg.pl

This script will sort fastq files to sample directories and creates job cards for all samples inside "/mnt/OMG/DONE" folder. It also makes "/mnt/OMG/SUBMIT" and "/mnt/OMG/BATCH" folder for cards movement. Now you have move job sample cards from "/mnt/OMG/DONE" folder to "/mnt/OMG/SUBMIT" folder manually.

Run perl batch queuing script .

      perl ~/qstaffer.pl

This script takes all sample cards from the "/mnt/OMG/SUBMIT" and relocate them to the "/mnt/OMG/BATCH" directory and submit them to the job queue. This directory is indication of what samples are currently in task-spooler queue. Cards of of finished samples should be returned to the "/mnt/OMG/DONE" (under construction).

You can trace jobs in task-spooler queue by issuing:

      tsp -l

Task spooler can flawlessly handle paralleled jobs you had send to the queue. However note that one sample is required up to 4 logical threads from processor and approximately up to 4-8 GB of memory in most pipeline steps however there is HaplotypeCaller step in GATK which uses up to 13GB of virtual memory per sample. In case you may have sufficient computer power to run multiple job in parallel you can change default available job slot from N=1 to the number you like with command:

      tsp -S  [N]

Safe [N] number for the machine with 36 CPUs and 64GB memory machine is 7 (seven samples in parallel). Such setup makes ETA for batch of 50 samples in about 3 days. Refer to the task-spooler documentation (tsp --help) for the additional queue management commands.

⁠Some valuable updates from 2021 by Terry Bertozzi

2/08/2019 : qc_collector.pl - Rust's script for summarizing Heterozygosity across all samples for omg pipe (docker) 19/3/2021 : Modified (and commented) by Terry Bertozzi to log errors rather than halting execution on the first error 1/4/2021 : Modified by Terry Bertozzi to add option of specifying input list for smaller runs and timestamping output files in ISO 8601 format to avoid over-writing Usage: perl qc_collector_tb.pl ... where

is the sample name directory created in /mnt/omg for each processed sample

For example

nohup  ~/qc_collector_tb.pl 36* > nohup.txt&

⁠tccpCoverage.pl by Terry Bertozzi

This script [ ~/tccpCoverage.pl ]takes the aligned exons output from the TCCP pipeline and will output; (1) a table of coverage stats (how much of target) for each sample (including alignment length for each locus and total concatenated alignment) (2) a summary table of numbers of loci above 10, 20, 30..90% thresholds for each sample

The script should be run from the directory containing the samples.

~/tccpCoverage.pl

==========================================================================

⁠Target Capture De-novo Phylogenetic Pipeline Major Stages:

Trimming and deduplication with trimmomatic.
Process sequences trough FASTQC to generate statistics for raw and trimmed data.
Map the nucleotide read sequences to protein targets container location: ~/omg/data. This is done to ensure that synonymous nucleotide changes do not alter the mapability of reads to what should be fairly conserved protein sequences. Reads mappped to the protein targets being collected to separate bins.
Paired-end reads from each protein target bin merged where possible (typically 20-45% in OMG samples) using bbmerge (BBtools).
De novo assembly with merged reads performed separately for each target using SPAdes. As an optimisation, each target gets a single processor thread and multiple assemblies are launched in parallel (4 cores by default). SPAdes creates assemblies with fixed kmer sizes of 25, 35, 45, 55, 65, and 75 bases, and combines the graphs produced for each. We leave error correction switched off and cut the tails of sequences where they have less than 5 reads coverage. SPAdes typically takes between 45 seconds and 5 minutes per target. However this is still the slowest stage of the pipeline may taking up to 10 hrs per sample.
Denovo assembled targets used as a genome reference for mapping reads back with 'bwa mem'.
Produced bam file is used for variant calling performed with GATK v3.8 and HaplotypeCaller.
VCF generated by HaplotypeCaller processed with another GATK tool ReadBackedPhasing and then annotated with VariantAnnotator and filtered to the targets.
Annotated and filtered VCF processed with GATK DepthOfCoverage and finally zipped.
Haplotypes and diplotypes reconstructed from VCF with "bcftools consensus" command.
All assembled nucleotide contigs are matched against their protein targets with Exonerate and the best matching contig with at least 300 nucleotides of perfect matches is mapped against the complete set of targets with Diamond. If the top target isn't the one we've been trying to find, then we assume that is a homologous and this target discarded from further consideration.

The next part of the pipeline replicates some of previous steps but dedicated to the mitochondrial targets:

Mitochondrial reference alignment done on reads which were unmapped to the original protein targets. Those unmapped reads extracted from bam files with sam2fq2 and shell scripts and mapped with 'bwa mem' to the reference fasta at location: ~ubuntu/omg/data/MitoDNA.fa.gz
Mitochondrial mapped reads merged with BBMerge version 37.50 followed with 'diamond blastx' for getting mitochondrial proteins targets.
Mitochondrial enriched reads de novo assembled with Spades on k 25,35,45,55,65,75
Denovo assembled mitochondrial reference reused for mapping paired end reads with 'bwa mem'.
Produced mitochondrial BAM file is used to get coverage statistics.

==========================================================================

Workflow Process Chart ⁠by Johan Gustafsson.

Workflow Process Chart with more focus on tools versionning and commands parameters ⁠by Cameron Jack.