White Papers

ManualsBrandsDell ManualsConverged InfrastructureHigh Performance Computing Solution Resources

Variant Calling Benchmark

Not Only Human

Variant call process refers to the identification of a nucleotide difference from reference sequences at a

given position in an individual genome or transcriptome. It includes single nucleotide polymorphism

(SNPs), insertion/deletions (indels) and structural variants. One of the most popular variant calling

applications is GenomeAnalysisTK (GATK) from Broad Institute. Often this GATK is used with BWA to

compose a variant calling workflow focusing on SNPs and indels. After we published Dell HPC System for

Genomics White Paper last year, there were significant changes in GATK. The key process, variant call

step UnifiedGenotyper is no longer recommended in their best practice. Hence, here we recreate BWA-

GATK pipeline according to the recommended practice to test whole genome sequencing data from

mammals and plants in addition to human’s whole genome sequencing data. This is a part of Dell’s

effort to help customers estimating their infrastructure needs for their various genomics data loads by

providing a comprehensive benchmark.

Variant Analysis for Whole Genome Sequencing data

System

The detailed configuration is in Dell HPC System for Genomics White Paper, and the summary of system

configuration and software is in Table 2.

Table 1 Server configuration and software

Component

Detail

Server

40x PowerEdge FC430 in FX2 chassis

Processor

Total of 1120 cores: Intel® Xeon® Dual E5-2695 v3 - 14 cores

Memory

128GB - 8x 16GB RDIMM, 2133 MT/s, Dual Rank, x4 Data Width

Storage

480TB IEEL (Lustre)

Interconnect

InfiniBand FDR

Red Hat Enterprise 6.6

Cluster Management tool

Bright Cluster Manager 7.1

Short Sequence Aligner

BWA 0.7.2-r1039

Variant Analysis

GATK 3.5

Utilities

sambamba 0.6.0, samtools 1.2.1

BWA-GATK pipeline

The current version of GATK is 3.5, and the actual workflow tested obtained from the workshop, ‘GATK

Best Practices and Beyond’. In this workshop, they introduce a new workflow with three phases.

• Best Practices Phase 1: Pre-processing

• Best Practices Phase 2A: Calling germline variants

• Best Practices Phase 2B: Calling somatic variants

• Best Practices Phase 3: Preliminary analyses

Here we tested out phase 1, phase 2A and phase3 for germline variant call pipeline. The details of

commands used in benchmark are listed below.

Summary of content (6 pages)

PAGE 1
Variant Calling Benchmark Not Only Human Variant call process refers to the identification of a nucleotide difference from reference sequences at a given position in an individual genome or transcriptome. It includes single nucleotide polymorphism (SNPs), insertion/deletions (indels) and structural variants. One of the most popular variant calling applications is GenomeAnalysisTK (GATK) from Broad Institute. Often this GATK is used with BWA to compose a variant calling workflow focusing on SNPs and indels.
PAGE 2
Phase 1. Pre-processing Step 1. Aligning and Sorting bwa mem -c 250 -M -t [number of threads] -R ‘@RG\tID:noID\tPL:illumine\tLB:noLB\tSM:bar’ [reference chromosome] [read fastq 1] [read fastq 2] | samtools view -bu - | sambamba sort -t [number of threads] -m 30G --tmpdir [path/to/temp] -o [sorted bam output] /dev/stdin Step 2. Mark and Remove Duplicates sambamba markdup -t [number of threads] --remove-duplicates --tmpdir=[path/to/temp] [input: sorted bam output] [output: bam without duplicates] Step 3.
PAGE 3
not required, but these steps are useful for the reporting purpose. If it is necessary, it can be added as a side workflow to the main procedure. For each job, 9 cores were assigned when 120 concurrent jobs were processed concurrently and 13 cores were used for the test of 80 concurrent jobs.
PAGE 4
Table 3 Mapping qualities of sequence reads data; obtained by using ‘samtools flagstat’. ‘Total QC-passed reads’ is the number of reads passed the criteria of sequencing quality. Among all QC-passed reads, the number of reads actually mapped on a reference genome and its percentage is on ‘Mapped reads (%)’ column. ‘Paired in sequencing’ column is the number of paired reads properly paired by a sequencer.
PAGE 5
Time Measurement Total run time is the elapsed wall time from the earliest start of Phase 1, Step 1 to the latest completion of Phase 3, Step 2. Time measurement for each step is from the latest completion time of the previous step to the latest completion time of the current step as described in Figure 1. Figure 1 measuring the running time The running time for each data set is summarized in Table 4.
PAGE 6
other mammals, and it became worse in plants. It is known that plants genomes contain large number of repeat sequences which make mapping process difficult. It is important to note that the shorter running time for HaplotypeCaller in rice does not reflect a real time since the size of the reference variant file was reduced significantly due to the chromosome length/position mismatches in the data.