benchseq.jl

using Bio.Seq


# Utils for Benchmark
# -------------------

function random_seq(n, nts, probs)
    cumprobs = cumsum(probs)
    x = Array(Char, n)
    for i in 1:n
        x[i] = nts[searchsorted(cumprobs, rand()).start]
    end
    return convert(ASCIIString, x)
end

function random_dna(n, probs=[0.24, 0.24, 0.24, 0.24, 0.04])
    return random_seq(n, ['A', 'C', 'G', 'T', 'N'], probs)
end

function gitrevision()
    commit = chomp(readall(`git rev-parse HEAD`))
    branch = chomp(readall(`git rev-parse --abbrev-ref HEAD`))
    return string(branch, '(', commit[1:7], ')')
end

macro bench(ex)
    name = ex.args[1]
    quote
        # JIT compiling
        $ex
        for i in 1:20
            t = @elapsed $ex
            println(revision, '\t', $name, '\t', size, '\t', t)
        end
    end
end


# Benchmark Set
# -------------

function access_all(s)
    for i in 1:endof(s)
        s[i]
    end
end

function fill_A!(s)
    for i in 1:endof(s)
        s[i] = DNA_A
    end
end

function make_copy(s)
    @assert ismutable(s)
    @assert copy(s) !== s
end

function convert_asciistring(s)
    convert(DNASequence, s)
end

function make_reverse(s)
    reverse(s)
end

function make_complement(s)
    complement(s)
end

function each_3mer_step1(s)
    for (i, kmer) in each(DNAKmer{3}, s)
    end
end

function each_3mer_step3(s)
    for (i, kmer) in each(DNAKmer{3}, s, 3)
    end
end

function count_nucs(s)
    c = NucleotideCounts(s)
    @assert c[DNA_A] + c[DNA_C] + c[DNA_G] + c[DNA_T] + c[DNA_N] == length(s)
end

function count_mismatches(s)
    c = mismatches(s, s)
    @assert c == 0
end


# Entry Point
# -----------

revision = gitrevision()
println("revision\tname\tsize\telapsed")
for size in logspace(1, 6, 6)
    s = DNASequence(random_dna(Int(size)))

    @bench access_all(s)
    @bench make_reverse(s)
    @bench make_complement(s)
    @bench each_3mer_step1(s)
    @bench each_3mer_step3(s)
    @bench count_nucs(s)
    @bench count_mismatches(s)

    let
        s = DNASequence(random_dna(Int(size)))
        mutable!(s)
        @bench fill_A!(s)
    end

    let
        s = DNASequence(random_dna(Int(size)))
        mutable!(s)
        @bench make_copy(s)
    end

    let s = random_dna(Int(size))
        @bench convert_asciistring(s)
    end
end