Edit me

Step 2.1a Alignment

Reads need to be aligned to the reference genome in order to identify the similar and polymorphic regions in the sample. As of 2016, the GATK team recommends their b37 bundle as the standard reference for Whole Exome and Whole Genome Sequencing analyses pending the completion of the GRcH38/Hg38 bundle 1. However, the 2018 functional equivalence specifications recommends the GRCh38DH from the 1000 Genomes project 2. Either way, a number of aligners can perform the alignment task.

Among these, BWA MEM 3 and bowtie2 4 have become trusted tools for short reads Illumina data, because they are accurate, fast, well supported, and open-source. Combined with variant callers, different aligners can offer different performance advantages with respect to SNPs, InDels and other structural variants, benchmarked in works like 3,5,6. Functional equivalence specifications recommends BWA-MEM v0.7.15 in particular (with at least the following parameters -K 100000000 -Y, and without -M so that split reads are marked as supplementary reads in congruence with BAM specification ).

The output file is usually in a binary BAM format 7, still taking tens or hundreds of Gigabytes of hard disk space. The alignment step tends to be I/O intensive, so it is useful to place the reference onto an SDD, as opposed to HDD, to speed up the process. The alignment can be easily parallelized by chunking the data into subsets of reads and aligning each subset independently, then combining the results.

Step 2.1b: De-duplication

The presence of duplicate reads in a sequencing project is a notorious problem. The causes are discussed in a blog post by Eric Vallabh Minikel (2012) 8. Duplicately sequenced molecules should not be counted as additional evidence for or against a putative variant – they must be removed prior to the analysis. A number of tools can be used including: samblaster 9, sambamba 10, the commercial novosort from the novocraft suit <11, Picard, and FASTX-Toolkit has fastx_collapser for this purpose. Additionally, MarkDuplicates is shipped as part of GATK4, but is called from Picard tools 12 in older GATK releases. For functional equivalence 2, it is recommended to use Picard tools v2.4.1.

De-duplication can also be performed by a simple in-house written Perl script.

Bibliography

  1. GATK Doc Article #1213. at https://software.broadinstitute.org/gatk/documentation/article.php?id=1213 

  2. Regier, A. A. et al. Functional equivalence of genome sequencing analysis pipelines enables harmonized variant calling across human genetics projects. BioRxiv (2018). doi:10.1101/269316  2

  3. Li, H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. (2013).  2

  4. Langmead, B. & Salzberg, S. L. Fast gapped-read alignment with Bowtie 2. Nat. Methods 9, 357–359 (2012). 

  5. Li, H. & Homer, N. A survey of sequence alignment algorithms for next-generation sequencing. Brief. Bioinformatics 11, 473–483 (2010). 

  6. Thankaswamy-Kosalai, S., Sen, P. & Nookaew, I. Evaluation and assessment of read-mapping by multiple next-generation sequencing aligners based on genome-wide characteristics. Genomics 109, 186–191 (2017). 

  7. Li, H. et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078–2079 (2009). 

  8. How PCR duplicates arise in next-generation sequencing. at http://www.cureffi.org/2012/12/11/how-pcr-duplicates-arise-in-next-generation-sequencing/ 

  9. Faust, G. G. & Hall, I. M. SAMBLASTER: fast duplicate marking and structural variant read extraction. Bioinformatics 30, 2503–2505 (2014). 

  10. Tarasov, A., Vilella, A. J., Cuppen, E., Nijman, I. J. & Prins, P. Sambamba: fast processing of NGS alignment formats. Bioinformatics 31, 2032–2034 (2015). 

  11. NOVOCRAFT TECHNOLOGIES SDN BHD. Novocraft. NOVOCRAFT TECHNOLOGIES SDN BHD at http://www.novocraft.com/ 

  12. Picard Tools - By Broad Institute. at https://broadinstitute.github.io/picard/