samblaster

repository·master·Indexed 42 days ago

https://github.com/gregoryfaust/samblaster

A high-performance tool for marking duplicate reads in read-id grouped SAM files and extracting structural variant-related reads—including discordant, split, and unmapped/clipped reads—into separate files. It is designed for use in pipelines such as BWA MEM and follows duplicate identification logic similar to Picard MarkDuplicates.

Tokens
2.7K
Snippets
4
Records
6
Agent score
4%

What's inside samblaster

  1. How samblaster identifies duplicates

    master

    A duplicate read pair is defined by a signature consisting of the combination of:

    1. Sequence name (QNAME)
    2. Strand
    3. Reference offset where the 5' end of the read would fall if fully aligned (calculated via POS, strand, and CIGAR).

    This signature matches the logic used by Picard MarkDuplicates.

    Key Rules:

    • For mapped pairs, both signatures must match.
    • For orphans (one side mapped), the signature of the mapped read must match a previously seen orphan. Since version 0.1.25, forward and reverse strand orphans/singletons can be duplicates of each other.
    • No doubly unmapped pair is marked as a duplicate.
    • Any secondary or supplementary alignment associated with a duplicate primary alignment is also marked as a duplicate.

    Requirement: Input SAM files must be read-id grouped (all alignments for a single QNAME must be in adjacent lines).

  2. How split reads and discordants are identified

    master

    Discordant Read Pairs

    A discordant pair meets these criteria:

    1. Both sides are mapped (neither FLAG 0x4 nor 0x8 is set).
    2. The properly paired FLAG (0x2) is not set.
    3. Secondary or supplementary alignments are not output as discordant.
    4. Duplicate pairs are output as discordant unless --excludeDups (-e) is used.

    Split Reads

    Split reads are identified when a single read has between two and --maxSplitCount primary and supplementary alignments. They are output if they are adjacent on the read and meet these criteria:

    • Each covers at least --minNonOverlap base pairs not covered by the other.
    • They map to different reference sequences/strands, OR they map to the same sequence/strand but represent a structural variant (SV) at least --minIndelSize in length with at most --maxUnmappedBases between them.

    Unmapped/Clipped Reads

    • Unmapped: A primary alignment with the unmapped read FLAG (0x4) set.
    • Clipped: A mapped primary alignment with a CIGAR string starting or ending with at least --minClipSize unaligned bases (S or H), provided it does not have supplementary alignments.

    Note: To output the full sequence for clipped reads, the input SAM must have soft-clipped primary alignments.

  3. Use samblaster to mark duplicates in a BWA MEM pipeline

    master

    The most common usage is to pipe alignments from bwa mem directly into samblaster, and then into samtools to create a compressed BAM file.

    Note: If you use the -M option in bwa mem, you must also use the -M option in samblaster for compatibility.

    Standard pipeline:

    bwa mem <idxbase> samp.r1.fq samp.r2.fq | samblaster | samtools view -Sb - > samp.out.bam

    With BWA -M compatibility:

    bwa mem -M <idxbase> samp.r1.fq samp.r2.fq | samblaster -M | samtools view -Sb - > samp.out.bam
    bwa mem <idxbase> samp.r1.fq samp.r2.fq | samblaster | samtools view -Sb - > samp.out.bam
  4. Install samblaster

    master

    samblaster is self-contained and requires only g++ and make. You can install it by cloning the repository, building with make, and moving the binary to your PATH.

    git clone git://github.com/GregoryFaust/samblaster.git
    cd samblaster
    make
    cp samblaster /usr/local/bin/.
    git clone git://github.com/GregoryFaust/samblaster.git
    cd samblaster
    make
    cp samblaster /usr/local/bin/.
  5. Extract discordant, split, and unmapped reads

    master

    samblaster can output specific types of alignments to separate files while processing the main stream.

    • -d or --discordantFile: Output discordant read pairs.
    • -s or --splitterFile: Output split reads (useful for SV detection).
    • -u or --unmappedFile: Output unmapped/clipped reads as FASTQ (if QUAL is available) or FASTA.
    • -e or --excludeDups: Exclude reads marked as duplicates from these extra files.

    Example: Extracting discordants and split reads from BWA MEM:

    bwa mem <idxbase> samp.r1.fq samp.r2.fq | samblaster -e -d samp.disc.sam -s samp.split.sam | samtools view -Sb - > samp.out.bam
    bwa mem <idxbase> samp.r1.fq samp.r2.fq | samblaster -e -d samp.disc.sam -s samp.split.sam | samtools view -Sb - > samp.out.bam
  6. Reference: samblaster CLI options

    master

    Input/Output Options

    OptionLong FlagDescription
    -i--inputInput sam file [stdin].
    -o--outputOutput sam file for all input alignments [stdout].
    -d--discordantFileOutput discordant read pairs to this file. [no discordant file output]
    -s--splitterFileOutput split reads to this file. [no splitter file output]
    -u--unmappedFileOutput unmapped/clipped reads as FASTQ/FASTA to this file. [no unmapped file output]

    Other Options

    OptionLong FlagDescription
    -a--acceptDupMarksAccept duplicate marks already in input file instead of looking for duplicates.
    -e--excludeDupsExclude reads marked as duplicates from discordant, splitter, and/or unmapped file.
    -r--removeDupsRemove duplicate reads from all output files (Implies --excludeDups).
    --addMateTagsAdd MC and MQ tags to all output paired-end SAM lines.
    --ignoreUnmatedSuppress abort on unmated alignments. Use only for singleton/filtered data.
    -MCompatibility mode for BWA MEM -M. Treats 0x100 and 0x800 as supplementary.
    --maxReadLengthINT Maximum allowed length of SEQ/QUAL string. [500]
    --maxSplitCountINT Max split alignments for a read to be in splitter file. [2]
    --maxUnmappedBasesINT Max un-aligned bases between alignments in splitter file. [50]
    --minIndelSizeINT Min SV feature size for split alignments in splitter file. [50]
    --minNonOverlapINT Min non-overlapping base pairs for split reads. [20]
    --minClipSizeINT Min bases a mapped read must be clipped to be in unmapped file. [20]
    -h--helpPrint help to stderr.
    -q--quietOutput fewer statistics.
    --versionPrint version to stderr.
    -i --input           FILE Input sam file [stdin].
    -o --output          FILE Output sam file for all input alignments [stdout].
    -d --discordantFile  FILE Output discordant read pairs to this file. [no discordant file output]
    -s --splitterFile    FILE Output split reads to this file abiding by paramaters below. [no splitter file output]
    -u --unmappedFile    FILE Output unmapped/clipped reads as FASTQ to this file abiding by parameters below. [no unmapped file output].
                              Requires soft clipping in input file.  Will output FASTQ if QUAL information available, otherwise FASTA.
    
    Other Options:
    -a --acceptDupMarks       Accept duplicate marks already in input file instead of looking for duplicates in the input.
    -e --excludeDups          Exclude reads marked as duplicates from discordant, splitter, and/or unmapped file.
    -r --removeDups           Remove duplicates reads from all output files. (Implies --excludeDups).
       --addMateTags          Add MC and MQ tags to all output paired-end SAM lines.
       --ignoreUnmated        Suppress abort on unmated alignments. Use only when sure input is read-id grouped, 
                              and either paired-end alignments have been filtered or the input file contains singleton reads.
                              <b>--ignoreUnmated is not recommended for general use on paired-end data. It disables checks that detect incorrectly sorted input.</b>
    -M                        Compatibility mode (details below); both FLAG 0x100 and 0x800 denote supplementary (chimeric). Similar to <i>bwa mem</i> <b>-M</b> option.
       --maxReadLength    INT Maximum allowed length of the SEQ/QUAL string in the input file. [500]
       --maxSplitCount    INT Maximum number of split alignments for a read to be included in splitter file. [2]
       --maxUnmappedBases INT Maximum number of un-aligned bases between two alignments to be included in splitter file. [50]
       --minIndelSize     INT Minimum structural variant feature size for split alignments to be included in splitter file. [50]
       --minNonOverlap    INT Minimum non-overlaping base pairs between two alignments for a read to be included in splitter file. [20]
       --minClipSize      INT Minumum number of bases a mapped read must be clipped to be included in unmapped file. [20]
    
    -h --help                 Print samblaster help to stderr.
    -q --quiet                Output fewer statistics.
       --version              Print samblaster version number to stderr.