Tumor Pair Sequencing Pipeline

Version 4.6.1

About

GenPipes Tumor Analysis pipeline is designed to detect somatic variants from a tumor and match normal sample pair more accurately. It offers three protocol options: sv, ensemble (default), and fastpass. See Schema tab for the pipeline workflow.

Fastpass

The fastpass option of the pipeline is meant for quick assessment using exome capture regions and the 1000bp flanking regions. The somatic/germline calls are made using one variant caller VarScan 2 with permissive variant calling thresholds.

Ensemble The pipeline’s ensemble protocol is based on an ensemble approach, which was optimized using both the DREAM3 challenge and the CEPH mixture datasets to select the best combination of callers for both SNV and structural variation detection. For SNVs, multiple callers such as GATK MuTect2, Strelka2, VarScan 2, and VarDict were combined for somatic calls to achieve a sensitivity of 98.1%, precision of 98.4%, and F1 score of 98.3% for variants found in ≥2 callers. For germline calls, Strelka2, VarScan 2 and VarDict calls were combined. Additional annotations are incorporated to the SNV calls using CPSR and/or PCGR, which produce html reports of germline variants implicated in cancer predisposition and somatic variants detected in the tumor, respectively. QC metrics are collected at various stages and visualized using MultiQC.

SV The sv protocol of the pipeline identifies SVs using GridSS to achieve a sensitivity of XX.X%, precision of XX.X%, and F1 score of XX.X% for duplication variants found in the DREAM3 dataset. The pipeline also estimates tumor purity and tumor ploidy of sample pair normal−tumor using specific cancer tool PURPLE, followed by annotation and visualization with Linx.

Usage

Schema

Steps

Fastpass	Ensemble	SV
Picard SAM to FastQ	Picard SAM to FastQ	Picard SAM to FastQ
Skewer Trimming	Skewer Trimming	Skewer Trimming
BWA Mem SAMBAMBA Sort SAM	BWA Mem SAMBAMBA Sort SAM	BWA Mem SAMBAMBA Sort SAM
SAMbamba Sort	SAMbamba Sort	SAMbamba Sort
SamBamba Merge Files	SamBamba Merge Files	SamBamba Merge Files
GATK InDel Realigner	GATK InDel Realigner	GATK InDel Realigner
SamBamba Merge Realigned	SamBamba Merge Realigned	SamBamba Merge Realigned
SamBamba Mark Duplicates	SamBamba Mark Duplicates	SamBamba Mark Duplicates
Recalibration	Recalibration	Recalibration
Manta SV Calls	Conpair Concordance Contamination	Manta SV Calls
Raw Mpileup Panel	Metrics DNA Picard Metrics	Strelka2 Paired Somatic
Paired VarScan 2	Metrics DNA Sample Qualimap	GRIDSS Paired Somatic
Merge VarScan 2 Panel	Metrics DNA FASTQ	Purple Ploidy Estimator
PreProcess VCF Panel	Sequenza Step	Linx Annotations Somatic
SNP Effect Panel	Strelka2 Paired Somatic	Linx Annotations Germline
Gemini Annotations Panel	Strelka2 Paired Germline	Linx Plot
Conpair Concordance Contamination	Strelka2 Paired Germline SnpEff
Metrics DNA Picard Metrics	Purple Ploidy Estimator
Metrics DNA Sample Qualimap	Raw Mpileup
Metrics DNA FASTQ	Paired Var Scan 2
Sequenza Step	Merge Var Scan 2
Run Pair MultiQC	Paired Mutect2
Sym Link Report	Merge Mutect2
Sym Link FASTQ Pair	VarDict Paired
Sym Link Panel	Merge Filter Paired VarDict
	Ensemble Somatic
	GATK Variant Annotator Somatic
	Merge GATK Variant Annotator Somatic
	Ensemble Germline Loh
	GATK Variant Annotator Germline
	Filter Ensemble Somatic
	Report CPSR
	Report PCGR
	Run Pair MultiQC
	Sym Link FASTQ Pair
	Sym Link Final BAM
	Sym Link Report
	Sym Link Ensemble

Step Details

Following are the various steps that are part of GenPipes Tumor Pair Pipeline:

Picard SAM to FastQ

This step converts SAM/BAM files from the input readset file into FASTQ format, if FASTQ files are not already specified in the readset file. Otherwise, it does nothing.

BWA Mem SAMBAMBA Sort SAM

In this step, the filtered reads are aligned to a reference genome. The alignment is done per sequencing readset. The BWA alignment software is with algorithm: bwa mem. BWA output BAM files are then sorted by coordinate using SamBamba.

This step uses trimmed FASTQ files, or FASTQ files from the readset file as input. If neither of these are available, then FASTQ output files from previous picard_sam_to_fastq conversion of BAM files is used as input.

SamBamba Merge Files

This step takes as input files:

Aligned and sorted BAM output files from previous BWA SAMBAMBA Sort SAM step, if available
Else, BAM files from the readset file

BAM readset files are merged into one file per sample. Merge is done using Picard.

GATK InDel Realigner

Insertion and deletion realignment is performed on regions where multiple base mismatches are preferred over indels by the aligner since it can appear to be less costly by the algorithm. Such regions will introduce false positive variant calls which may be filtered out by realigning those regions properly. Realignment is done using GATK. The reference genome is divided by a number regions given by the nb_jobs parameter.

Note: modified to use both normal and tumor BAMs to reduce FPs around indels

SamBamba Merge Realigned

BAM files of regions of realigned reads are merged per sample using SamBamba.

Sambamba Mark Duplicates

Mark duplicates. Aligned reads per sample are duplicates if they have the same 5’ alignment positions (for both mates in the case of paired-end reads). All but the best pair (based on alignment score) will be marked as a duplicate in the BAM file. Marking duplicates is done using Picard.

Recalibration

This step re-calibrates base quality scores of sequencing-by-synthesis reads in an aligned BAM file. After recalibration, the quality scores in the QUAL field in each read in the output BAM are more accurate in that the reported quality score is closer to its actual probability of mismatching the reference genome. Moreover, the recalibration tool attempts to correct for variation in quality with machine cycle and sequence context, and by doing so, provides not only more accurate quality scores but also more widely dispersed ones.

Conpair Concorance Contamination

Conpair is a concordance and contamination estimator for tumor–normal pairs. This step is a quality control process to ensure the normal sample and cancer sample come from the same patient. It estimates this by determining the amount of normal sample in the tumor and the amount of tumor in normal sample.

Metrics DNA Picard Metrics

This step runs specific QC metrics on DNA data Functions:

collect_multiple_metrics,
collectOxoGMetrics,
collect_sequencing_artifacts_metrics

using Picard Metric Definition Tools.

Metrics DNA Sample Qualimap

In this step, QC alignment metrics generated by Qualimap.

Metrics DNA FASTQ

In this step, quality controlling metrics are generated on the read level using FastQC.

Sequenza Step

Sequenza is a novel set of tools providing a fast python script to genotype cancer samples, and an R package to estimate cancer cellularity, ploidy, genome wide copy number profile and infer for mutated alleles.

Strelka2 Paired Somatic

Strelka2 is a fast and accurate small variant caller optimized for analysis of germline variation in small cohorts and somatic variation in tumor/normal sample pairs This implementation is optimized for somatic calling.

Strelka2 Paired Germline

Strelka2 is a fast and accurate small variant caller optimized for analysis of germline variation in small cohorts and somatic variation in tumor/normal sample pairs This implementation is optimized for germline calling in cancer pairs.

Strelka2 Paired Germline SnpEff

Strelka2 is a fast and accurate small variant caller optimized for analysis of germline variation in small cohorts and somatic variation in tumor/normal sample pairs This implementation is optimized for germline calling in cancer pairs.

Purple Ploidy Estimator

Purple is a purity ploidy estimator for whole genome sequenced (WGS) data. It combines B-allele frequency (BAF) from AMBER, read depth ratios from COBALT, somatic variants and structural variants to estimate the purity and copy number profile of a tumor sample.

Raw Mpileup

Full pileup (optional). A raw Mpileup file is created using samtools Mpileup and compressed in gzipped format. One packaged Mpileup file is created per sample/chromosome.

Paired Var Scan 2

Variant calling and somatic mutation/CNV detection for next-generation sequencing data. Uses VarScan 2 for Somatic mutation and copy number alteration discovery in cancer by exome sequencing VarScan 2 thresholds based on DREAM3 results generated by SC INFO field remove to prevent collision with Samtools output during ensemble.

Merge Var Scan 2

Merge Mpileup files per sample/chromosome into one compressed gzip file per sample.

Paired Mutect2

GATK Mutect2 Overview caller for SNVs and Indels.

Merge Mutect2

Merge SNVs and indels for GATK MuTect2 Overview Replace TUMOR and NORMAL sample names in VCF to the exact tumor/normal sample names Generate a somatic VCF containing only PASS variants

VarDict Paired

In this step, VarDict caller is used for SNVs and Indels. Note: variants are filtered to remove instance where REF == ALT and REF modified to ‘N’ when REF is AUPAC nomenclature

Merge Filter Paired VarDict

The fully merged VCF is filtered using following steps:

Retain only variants designated as somatic by VarDict: either StrongSomatic or LikelySomatic
Somatics identified in step 1 must have PASS filter.

Ensemble Somatic

Apply Bcbio.variations ensemble approach for GATK MuTect2 Overview, VarDict, Samtools and VarScan 2 calls Filter ensemble calls to retain only calls overlapping 2 or more callers

GATK Variant Annotator Somatic

Add VCF annotations to ensemble VCF: Standard and Somatic annotations.

Merge GATK Variant Annotator Somatic

Merge annotated somatic VCFs.

Compute Cancer Efects Somatic

Variant effect annotation. The .vcf files are annotated for variant effects using the SnpEff software. SnpEff annotates and predicts the effects of variants on genes (such as amino acid changes). Modified arguments to consider paired cancer data.

Ensemble Somatic dbNSFP Annotation

This step performs additional Single Nucleotide Variants (SVN) annotations. It provides extra information about SVN by using numerous published databases. This is applicable to human samples. Databases that are available include Biomart (adds GO annotations based on gene information) and dbNSFP (an integrated database of functional annotations from multiple sources for the comprehensive collection of human non-synonymous SNPs. It compiles prediction scores from four prediction algorithms (SIFT, Polyphen2, LRT and MutationTaster), three conservation scores (PhyloP, GERP++ and SiPhy) and other function annotations).

Sample Gemini Annotations Somatic

This step loads functionally annotated .vcf file into a MySQL lite annotation database - Gemini.

Ensemble Germline Loh

This step applies Bcbio.variations ensemble approach for VarDict, Samtools and VarScan 2 calls. It also filters ensemble calls to retain only calls overlapping 2 or more callers.

GATK Variant Annotator Germline

Add VCF annotations to ensemble VCF: most importantly the AD field.

Merge GATK Variant Annotator Germline

Merge annotated germline and LOH VCFs.

Compute Cancer Effects Germline

Variant effect annotation. The .vcf files are annotated for variant effects using the SnpEff software. SnpEff annotates and predicts the effects of variants on genes (such as amino acid changes). Modified arguments to consider paired cancer data.

Ensemble Germline dbNSFP Annotation

Additional SVN annotations. Provides extra information about SVN by using numerous published databases. Applicable to human samples. Databases available include Biomart (adds GO annotations based on gene information) and dbNSFP (an integrated database of functional annotations from multiple sources for the comprehensive collection of human non-synonymous SNPs. It compiles prediction scores from four prediction algorithms (SIFT, Polyphen2, LRT and MutationTaster), three conservation scores (PhyloP, GERP++ and SiPhy) and other function annotations).

Sample Gemini Annotations Germline

Load functionally annotated vcf file into a MySQL lite annotation database : Gemini

Run Pair MultiQC

Aggregate results from bioinformatics analyses across many samples into a single report. MultiQC searches a given directory for analysis logs and compiles a HTML report. It’s a general use tool, perfect for summarizing the output from numerous bioinformatics tools. For details, refer to MultiQC Info.

Sym Link FASTQ Pair

Generates symbolic links for FASTQ Pair output files.

Sym Link Final BAM

Generates symbolic links for Final BAM output files.

Sym Link Report

Generates symbolic links for Tumor Pair pipeline report.

Sym Link Ensemble

Generates symbolic links for Ensemble processing output.

Sym Link Panel

This step create sym links of panel variants for generating reports and deliverable to the clients.

Manta SV Calls

In this step`Manta`_ is deployed for calling structural variants (SVs) and indels from mapped paired-end sequencing reads. It is optimized for the following analysis:

Germline variation in small sets of individuals.
Somatic variation in tumor/normal sample pairs.

Manta discovers, assembles and scores large-scale SVs, medium-sized indels and large insertions within a single efficient workflow. The Manta step accepts input read mappings from BAM or CRAM files and reports all SV and indel inferences in VCF 4.1 format.

Raw Mpileup Panel

Full pileup (optional). A raw Mpileup file is created using samtools Mpileup and compressed in gzipped format. One packaged Mpileup file is created per sample/chromosome.

Paired VarScan 2

This step involves variant calling and somatic mutation/CNV detection for next-generation sequencing data. For details see - VarScan 2 Somatic mutation and copy number alteration discovery in cancer by exome sequencing.

Merge VarScan 2 Panel

Merge Mpileup files per sample/chromosome into one compressed gzip file per sample.

PreProcess VCF Panel

Preprocess VCF for loading into a annotation database called Gemini. Processes include normalization and decomposition of MNPs by Vt and VCF modification for correct loading into Gemini.

SNP Effect Panel

Variant effect annotation. The .vcf files are annotated for variant effects using the SnpEff software. SnpEff annotates and predicts the effects of variants on genes (such as amino acid changes).

Gemini Annotations Panel

Load functionally annotated VCF file into a MySQL lite annotation database Gemini.

Delly Call Filter

Delly2 is an integrated structural variant (SV) prediction method. It can discover, genotype and visualize deletions, tandem duplications, inversions and translocations at single-nucleotide resolution, in case of short-read, massively parallel sequencing data.

Delly2 uses paired-ends and split-reads to sensitively and accurately delineate genomic rearrangements throughout the genome. SVs can be visualized using Delly-maze and Delly-suave.

This step takes as input the normal and tumor final BAM files. It returns output in the form of .bcf file.

Delly SV Annotation

The Delly SV Annotation or SV Genotype is typically a time consuming step. This exhaustive and comprehensive step interrogates every read overlapping putative SV breakpoints so that Delly can yield accurate read counts for reference-supporting and alternative allele-supporting reads and paired-ends.

Manta SV Annotation

This step uses Manta results and performs SV filtering and annotations to enable data visualizations. It uses Manta SV call step output VCF file and performs SV annotations to compares types and breakpoints for candidate SVs from different callsets and enables fast comparison of SVs to genomic features such as genes and repetitive regions, as well as to previously established SV datasets.

Lumpy Paired SV

This step uses Lumpy for structural variant discovery in the given input file. The output is available in BAM format.

Comprehensive discovery of structural variation (SV) from whole genome sequencing data requires multiple detection signals including read-pair, split-read, read-depth and prior knowledge. Owing to technical challenges, extant SV discovery algorithms either use one signal in isolation, or at best use two sequentially. Lumpy is a novel SV discovery framework that naturally integrates multiple SV signals jointly across multiple samples. It yields improved sensitivity, especially when SV signal is reduced owing to either low coverage data or low intra-sample variant allele frequency.

Lumpy SV Annotation

This step performs LumPy SV Annotation for mapping and characterization of SVs.

Wham Call SV

WHAM (Whole-genome Alignment Metrics) provides a single, integrated framework for both structural variant calling and association testing, thereby bypassing many of the difficulties that currently frustrate attempts to employ SVs in association testing. This step returns a VCF file.

Wham SV Annotation

This step uses the VCF file generated in the previous step Wham Call SV and performs SV annotations.

CNVKit Batch

CNVkit is a Python library and command-line software toolkit to infer and visualize copy number from high-throughput DNA sequencing data. It is designed for use with hybrid capture, including both whole-exome and custom target panels, and short-read sequencing platforms such as Illumina and Ion Torrent.

In this step, the CNVKit library used to run ‘batch’ command. Following steps are automated with the CNVKit batch command. #. Create a pooled reference of per-bin copy number estimates from several normal samples; then #. Use this reference in processing all tumor samples that were sequenced with the same platform and library prep.

CNVkit uses the bait BED file, reference genome sequence, and (optionally) sequencing-accessible regions along with your BAM files as input for the batch command.

CNVKit SV Annotation

In this step, structural variant annotation is carried out after CNVKit batch processing.

Scones

This step aims to estimate somatic copy number variation using BVAtools and SCONE.

BVAtools generate the binned Depth ratio values from the tumor and normal BAM files. SCONE is tool to deconvolution the logR Signal of the tumor-normal coverage into a mixture of Bayesian sub-signal for each copy number state. The result is a set of several deconvolution using 0-7 sub-signal. As each tumor sample is unique the choice of the best final model (number of sub-signal) needs to be manually evaluated using the log ratio graphical representation.

SvABA SV Assemble

This step performs structural variation and indel analysis by assembly. It uses SvABA tool to perform this analysis. SvABA is an efficient and accurate method for detecting SVs from short-read sequencing data using genome-wide local assembly with low memory and computing requirements.

SvABA SV Annotation

In this step, the output generated by SvABA SV Assemble step as VCF file is used as input and SV annotations are done.

Ensemble MetaSV Somatic

This step uses MetaSV Tool to perform an accurate and integrative structural-variant calling analysis. MetaSV, an integrated SV caller which leverages multiple orthogonal SV signals for high accuracy and resolution. MetaSV proceeds by merging SVs from multiple tools for all types of SVs. It also analyzes soft-clipped reads from alignment to detect insertions accurately since existing tools underestimate insertion SVs. Local assembly in combination with dynamic programming is used to improve breakpoint resolution. Paired-end and coverage information is used to predict SV genotypes.

Ensemble MetaSV Germline

Performs MetaSV Germline analysis.

MetaSV SV Annotation

In this step the output of MetaSV analysis earlier is annotated for structural variant detection.

Sym Link Sequenza

In this step, sym links for Sequenza output are generated.

Sym Link MetaSV

In this step, sym links for MetaSV output are generated.

Sym Link Manta

In this step, sym links for Manta output are generated.

Sym Link Delly

In this step, sym links for Delly output are generated.

Sym Link Lumpy

In this step, sym links for Lumpy output are generated.

Sym Link Wham

In this step, sym links for Wham output are generated.

Sym Link CNVKit

In this step, sym links for CNVKit output are generated.

SAMbamba Sort

BWA output BAM files are sorted by coordinate using Sambamba This step takes as input files: BAM file from bwa mem alignment compressed by SAMbamba.

Filter Ensemble Germline

Applies custom script to inject FORMAT information - tumor/normal DP and VAP into the INFO field the filter on those generated fields.

Filter Ensemble Germline

Applies custom script to inject FORMAT information - tumor/normal DP and VAP into the INFO field the filter on those generated fields.

Report CPSR

Creates a CPSR Germline Report where input is the filtered ensemble germline vcf and output is the html report and additional flat files.

Report PCGR

Creates a PCGR Somatic Germline Report where input is the filtered ensemble germline vcf and the output is the html report and additional flat files.

GRIDSS Paired Somatic

GRIDSS is a module software suite containing tools useful for the detection of genomic rearrangements. GRIDSS includes a genome-wide break-end assembler, as well as a structural variation caller for Illumina sequencing data. GRIDSS calls variants based on alignment-guided positional de Bruijn graph genome-wide break-end assembly, split read, and read pair evidence.

Linx Annotations Somatic

The LINX algorithm classifies somatic structural variation in tumors. It is an annotation, interpretation and visualization tool for structural variants. The primary function is grouping together individual SV calls into distinct events and properly classify and annotating the event to understand both its mechanism and genomic impact.

LINX helps in the interpretation of structural variant and copy number data derived from short-read, whole-genome sequencing. It classifies raw structural variant calls into distinct events and predicts their effect on the local structure of the derivative chromosome and the functional impact on affected genes. Visualizations facilitate further investigation of complex rearrangements. LINX allows insights into a diverse range of structural variation events and can reliably detect pathogenic rearrangements, including gene fusions, immunoglobulin enhancer rearrangements, intragenic deletions, and duplications.

Linx Annotations Germline

Though LINX is designed primarily for somatic SV, it can also be run in a more limited germline mode to annotate and interpret germline rearrangements

Linx Plot

LINX Plot offers visualization methods that provide insight into complex genomic rearrangements. It leverages the integrated structural variant and copy number calling to cluster individual structural variants into higher order events and chains them together to predict local derivative chromosome structure.

Skewer Trimming

Skewer is used for adapter trimming in GenPipes pipelines. Adapter trimming is a prerequisite step for analyzing next-generation sequencing (NGS) data when the reads are longer than the target DNA/RNA fragments. Although typically used in small RNA sequencing, adapter trimming is also used widely in other applications, such as genome DNA sequencing and transcriptome RNA/cDNA sequencing, where fragments shorter than a read are sometimes obtained because of the limitations of NGS protocols. For the newly emerged Nextera long mate-pair (LMP) protocol, junction adapters are located in the middle of all properly constructed fragments; hence, adapter trimming is essential to gain the correct paired reads.

BWA SAMbamba Sort SAM

The input for this step is the trimmed FASTQ files if available. Otherwise, it uses the FASTQ files from the readset. If those are not available then it uses FASTQ output files from the previous ‘Picard SAM to FASTQ`_ step where BAM files are converted to FASTQ format.

In this step, filtered reads are aligned to a reference genome. The alignment is done per sequencing readset. The alignment software used is BWA with algorithm: bwa mem. BWA output BAM files are then sorted by coordinate using SAMBAMBA.

SAMBAM Merge SAM Files

This step takes as input files:

Aligned and sorted BAM output files from previous bwa_mem_picard_sort_sam step if available
Else, BAM files from the readset file

In this step, BAM readset files are merged into one file per sample. Merge is done using Picard.

SAMBAM Merge Realigned

In this step, BAM files of regions of realigned reads are merged per sample using SAMBAMBA.

Details

Tumor Pair pipeline helps in inferring the cancer cell copy number to normal cell copy number.

Human genome comprises of a set of chromosome pairs. One chromosome in each pair, called homolog, is derived from each parent. It is typically referred to as diploid whereas the set of chromosomes from a single parent is called haploid genome. For a given gene on a given chromosome, there is a comparable, if not identical, gene on the other chromosome in the pair, known as an allele. Large structural alterations in chromosomes can change the number of copies of affected genes on those chromosomes. This is one of the key reasons for causing tumors or cancer. In cancer cells, instead of having a homologous pair of alleles for a given gene, there may be deletions or duplications of those genes.

Such alterations leads to unequal contribution of one allele over the other, altering the copy number of a given allele. These variations in copy number indicated by the ratio of cancer cell copy number to normal cell copy number can provide information regarding the structure and history of cancer. However, when DNA is extracted, there is a mix of cancer and normal cells and this information regarding absolute copy number per cancer cell is lost in DNA extraction process. Hence it must be inferred.

Inferring absolute copy number is difficult for three reasons:

Cancer cells are nearly always intermixed with an unknown fraction of normal cells; the measure for this is tumor purity.
The actual quantity of DNA in the cancer cells after gross numerical and structural chromosomal changes is unknown; the measure for this is tumor ploidy.
The cancer cell population may be heterogeneous, possibly because of ongoing mutations and changes.

Tumor purity and ploidy have a substantial impact on next-gen sequence analyses of tumor samples and may alter the biological and clinical interpretation of results.

More information

For the latest implementation and usage details see pipeline implementation README.md.

MuTect2 Tool for calling somatic SNVs and indels via local assembly of haplotypes - See here.
A three-caller pipeline for variant analysis of cancer whole-exome sequencing data.