A Snakemake workflow to process single samples or cohorts of Illumina paired-end sequencing data (WGS or WES) using trim galore/bwa/GATK4/parabricks.
public
1yr ago
Version:
v2.1.0
0 bookmarks
A Snakemake workflow to process single samples (unrelated individuals) or cohort samples (related individuals) of paired-end sequencing data (WGS or WES) using bwa and GATK4 . Quality control checks are also undertaken. The fastq files can be optionally trimmed with Trim Galore and the pipeline can run on NVIDIA GPU's where nvidia clara parabricks software is available for significant speedups in analysis times. This workflow is designed to follow the GATK best practice workflow for germline short variant discovery (SNPs + Indels) . This pipeline is designed to be followed by vcf_annotation_pipeline and the data ingested into scout for clinical interpretation. However, this pipeline also stands on it's own, taking the data from fastq to vcf (raw sequencing data to called variants). This pipeline has been developed with human genetic data in mind, however we designed it to be species agnostic. Genetic data from other species can be analysed by setting a species-specific reference genome and variant databases in the configuration file (but not all situations have been tested).
Pipeline summary - single samples
-
Adapter trimming ( Trim Galore ) ( optional )
-
Alignment against reference genome ( Burrows-Wheeler Aligner )
-
Mark duplicates ( GATK MarkDuplicates )
-
Base recalibration ( GATK BaseRecalibrator and GATK ApplyBQSR )
-
Haplotype calling ( GATK HaplotypeCalller )
Pipeline summary - single samples - GPU accelerated
-
Adapter trimming ( Trim Galore ) ( optional )
-
Alignment against reference genome, mark duplicates, base recalibration and haplotype calling ( parabricks germline pipeline )
- Equivilant to Burrows-Wheeler Aligner , GATK MarkDuplicates , GATK BaseRecalibrator , GATK ApplyBQSR and GATK HaplotypeCalller
Pipeline summary - cohort samples
-
Adapter trimming ( Trim Galore ) ( optional )
-
Alignment against reference genome ( Burrows-Wheeler Aligner )
-
Mark duplicates ( GATK MarkDuplicates )
-
Base recalibration ( GATK BaseRecalibrator and GATK ApplyBQSR )
-
Haplotype calling ( GATK HaplotypeCalller )
-
Combine GVCF into multi-sample GVCF ( GATK CombineGVCFs )
-
Genotyping ( GATK GenotypeGVCFs )
Pipeline summary - cohort samples - GPU accelerated
-
Adapter trimming ( Trim Galore ) ( optional )
-
Alignment against reference genome, mark duplicates, base recalibration and haplotype calling ( parabricks germline pipeline )
- Equivilant to Burrows-Wheeler Aligner , GATK MarkDuplicates , GATK BaseRecalibrator , GATK ApplyBQSR and GATK HaplotypeCalller
-
Combine GVCF into multi-sample GVCF ( parabricks trio combine gvcf )
- Equivalent to GATK CombineGVCFs
-
Genotyping ( GATK GenotypeGVCFs )
Main output files
Single samples:
-
results/qc/multiqc_report.html -
results/mapped/sample1_recalibrated.bam -
results/called/sample1_raw_snps_indels.vcf
Cohort samples:
-
results/qc/multiqc_report.html -
results/mapped/sample1_recalibrated.bam -
results/mapped/sample2_recalibrated.bam -
results/mapped/sample3_recalibrated.bam -
results/called/proband1_raw_snps_indels.vcf
Prerequisites
-
Prerequisite hardware: NVIDIA GPUs (for GPU accelerated runs)
-
Prerequisite software: NVIDIA CLARA parabricks and dependencies (for GPU accelerated runs), Git (tested with version 2.7.4), Mamba (tested with version 0.4.4) with Conda (tested with version 4.8.2), gsutil (tested with version 4.52), gunzip (tested with version 1.6)
Test human_genomics_pipeline
The provided test dataset can be used to test running this pipeline on a new machine, or test pipeline developments/releases.
Run human_genomics_pipeline
See the docs for a walkthrough guide for running human_genomics_pipeline on:
Contribute back!
-
Raise issues in the issues page
-
Create feature requests in the issues page
-
Start a discussion in the discussion page
-
Contribute your code! Create your own branch from the development branch and create a pull request to the development branch once the code is on point!
Contributions and feedback are always welcome! :blush:
Code Snippets
21 22 | shell: "bwa mem -R {params.readgroup} -t {threads} -K 10000000 {input.refgenome} {input.fastq} | gatk SortSam --java-options {params.maxmemory} {params.sortsam} -O={output} --TMP_DIR={params.tdir} &> {log}" |
16 17 | shell: "fastqc {input} -o ../results/qc/fastqc/ -t {threads} &> {log}" |
20 21 | shell: "gatk ApplyBQSR --java-options {params.maxmemory} -I {input.bam} -bqsr {input.recal} -R {input.refgenome} -O {output} {params.padding} {params.intervals}" |
21 22 | shell: "gatk BaseRecalibrator --java-options {params.maxmemory} -I {input.bams} -R {input.refgenome} -O {output} --tmp-dir {params.tdir} {params.padding} {params.intervals} {params.recalibration_resources} &> {log}" |
22 23 | shell: "gatk CombineGVCFs -R {input.refgenome} {params.command} -O {output.vcf} --tmp-dir {params.tdir} {params.other} &> {log}" |
21 22 | shell: "gatk GenotypeGVCFs --java-options {params.maxmemory} -R {input.refgenome} -V {input.gvcf} -O {output} {params.padding} {params.intervals} {params.other} &> {log}" |
23 24 | shell: "gatk HaplotypeCaller --java-options {params.maxmemory} -I {input.bams} -R {input.refgenome} -D {input.dbsnp} -O {output.vcf} --tmp-dir {params.tdir} {params.padding} {params.intervals} {params.other} &> {log}" |
21 22 | shell: "gatk HaplotypeCaller --java-options {params.maxmemory} -I {input.bams} -R {input.refgenome} -D {input.dbsnp} -O {output} --tmp-dir {params.tdir} {params.padding} {params.intervals} &> {log}" |
18 19 | shell: "gatk MarkDuplicates --java-options {params.maxmemory} -I {input} -O {output.bam} -M {output.metrics} --TMP_DIR {params.tdir} &> {log}" |
12 13 | shell: "multiqc {input} -o ../results/qc/ &> {log}" |
14 15 | shell: "multiqc {input.fastqc} {input.trimming1} {input.trimming2} -o ../results/qc/ &> {log}" |
26 27 | shell: "pbrun germline --ref {input.refgenome} --in-fq {input.fastq} {params.recalibration_resources} --out-bam {output.bam} --out-variants {output.vcf} --out-recal-file {output.recal} --num-gpus {resources.gpu} {params.readgroup} --tmp-dir {params.tdir} {params.padding} {params.intervals} {params.other_params} --num-cpu-threads {threads} &> {log}" |
17 18 | shell: "pbrun triocombinegvcf --ref {input.refgenome} {params.command} --out-variants {output} &> {log}" |
22 23 | shell: "trim_galore {input} -o ../results/trimmed/ {params.adapters} {params.other} -j {threads} &> {log}" |
Support
- Share:
-
-
-
- Future updates
Related Workflows
