Author: Khalid Mahmood ([email protected])

hiplexpipe is based on the Ruffus library for writing bioinformatics pipelines. Its features include:

• Job submission on a cluster using DRMAA (currently only tested with SLURM).
• Job dependency calculation and checkpointing.
• Pipeline can be displayed as a flowchart.
• Re-running a pipeline will start from the most up-to-date stage. It will not redo previously completed tasks.

See LICENSE.txt in source repository.

hiplexpipe depends on the following programs and libraries:

• python (version 2.7.5)
• java (version 1.8)
• DRMAA for submitting jobs to the cluster (it uses the Python wrapper to do this). You need to install your own libdrama.so for your local job submission system. There are versions available for common schedulers such as Torque/PBS, SLURM and so on.

• SAMtools (version 1.3.1)
• bwa for aligning reads to the reference genome (version 0.7.15)
• GATK for calling variants and genotyping (version 3.6)
• BEDTools for calculating sequencing coverage statistics (version 2.26.0)

hiplexpipe assumes the tools above are installed by the users themselves.

hiplexpipe depends on the following python libraries, tools and wrappers.

• Python 2.7.5
• PyVCF
• Biopython
• pybedtools
• cyvcf2

We recommend using a python virtual environment. Following is an examples of how to setup a hiplexpipe virtual environment ready for analysis:

The following instruction are based on the Melbourne Bioinformatics computing infrastructure.

module load Python/2.7.12-GCC-4.9.3
export DRMAA_LIBRARY_PATH=/usr/local/slurm_drmaa/1.0.7-GCC/lib/libdrmaa.so
virtualenv --system-site-packages hiplexpipe
source hiplexpipe/bin/activate
pip install jupyter
pip install plotly
pip install pybedtools
pip install -U https://github.com/khalidm/undr_rover/archive/master.zip
pip install -U https://github.com/khalidm/hiplexpipe/archive/simple.zip
pip install -U https://github.com/khalidm/offtarget/archive/master.zip
mkdir references
mkdir coverage

hiplexpipe --config pipeline.config --use_threads --log_file pipeline.log --jobs 10 --verbose 3 --just_print

You should have two target interval files for every Hi-Plex experiment.

• rover.txt - this contains the amplicon regions and primer sequences.
• idt.txt - this file contains the primer sequences and their names matching the names in the above rover.txt file.

Make sure heel sequences are removed from rover.txt file

Follow instructions below to prepare the intervals files for the pipeline. (We are working on a tool to automate this task).

Each interval in this bed file is the midpoint of each amplicon. This file is used to calculate alignment and coverage statistics.

cut -f1,2,3,4,5 rover.txt > rover.bed

or

awk ' BEGIN{FS="\t";OFS="\t"}; { print $1,int($2+($3-$2)/2),int($3-($3-$2)/2),$4,$5} ' rover.txt > rover.bed

This file is used to clip primer sequences from the alignments.

awk ' BEGIN{FS="\t";OFS="\t"}; { print $1,$7,$8,$1,$12,$11} ' rover.txt > primer.bedpe

This creates a bed file of intervals for GATK variant calling. Note this is different from rover.bed as this merges overlapping targets and mainly functions to provide a target for variant calling.

cut -f1,2,3 rover.txt | bedtools slop -i - -b 10 -g hg19.genome | bedtools merge -i - > rover.gatk.bed

module load Python/2.7.12-GCC-4.9.3
export DRMAA_LIBRARY_PATH=/usr/local/slurm_drmaa/1.0.7-GCC/lib/libdrmaa.so
module load BEDTools/2.26.0-vlsci_intel-2015.08.25
module load SAMtools
module load VCFtools

I have created a template config file (pipeline_template.config) for all these analysis.

1. Create a new directory for the analysis
2. Make a copy of pipeline_template.config in the new analysis directory.
3. Make relevant changes to the new config file.

• change pipeline_id
• add fastq file paths
• under the comment "hiplex files" - amend paths to files relevant to the design

Log into snowy (HPC)

Run following commands:

module load Python/2.7.12-GCC-4.9.3
screen -S new_analysis
module purge
module load vlsci
module load Python/2.7.12-GCC-4.9.3
module load SAMtools
module load VCFtools
source hiplexpipe/bin/activate

hiplexpipe --config pipeline.config --use_threads --log_file pipeline.log --jobs 50 --verbose 2

From within the virtualenv, run the following command:

python alignment_stats.py > stats.txt

This will generate a table containing various alignment statistics for each sample.

jupyter nbconvert --ExecutePreprocessor.timeout=6000 --to html --execute coverage_analysis_main.ipynb

This will output coverage_analysis_main.html file.

Generates statistics on which amplicons are mapping to incorrect regions of the genome, or not mapping at all.

Run for a few samples picked at random.

offtarget --primers <rover.txt> --fastq1 <fastq_read_1> --fastq2 <fastq_read_2> --bam <sorted bam file> --log offtarget.log > output.txt