Sequence Alignment

BioLang provides pairwise sequence alignment via align(), supporting global (Needleman-Wunsch), local (Smith-Waterman), and semi-global modes. The implementation uses a linear gap score and returns a structured result including the gapped sequences, identity fraction, alignment score, and a CIGAR string.

align(seq1, seq2, mode?, match?, mismatch?, gap?)

Parameter	Type	Default	Description
`seq1`	DNA / RNA / Protein / Str	—	First sequence
`seq2`	DNA / RNA / Protein / Str	—	Second sequence
`mode`	Str	`"global"`	`"global"`, `"local"`, or `"semi"`
`match`	Int	`2`	Score for a matching base or residue
`mismatch`	Int	`-1`	Score for a mismatch
`gap`	Int	`-2`	Linear gap score

Returns a record:

Field	Type	Description
`score`	Int	Alignment score
`aligned_a`	Str	First sequence with gap characters (`-`)
`aligned_b`	Str	Second sequence with gap characters
`identity`	Float	Identical columns ÷ alignment length
`cigar`	Str	CIGAR string: M=match, X=mismatch, I=insertion, D=deletion

Global alignment — Needleman-Wunsch

Aligns two sequences end-to-end, penalizing all gaps including terminal ones.

let a = dna"ACGTACGT"
let b = dna"ACGTTACGT"

let r = align(a, b)
println(r.score)       # alignment score
println(r.aligned_a)   # "ACG-TACGT"
println(r.aligned_b)   # "ACGTTACGT"
println(r.identity)    # 0.888...
println(r.cigar)       # "3M1I5M"

Local alignment — Smith-Waterman

Finds the highest-scoring subsequence alignment, ignoring flanking regions.

let genome = dna"AAACGTACGTGGGGG"
let query  = dna"ACGTACGT"

let hit = align(genome, query, "local")
println(hit.score)          # local best score
println(hit.aligned_a)      # best local region of genome
println(hit.aligned_b)      # corresponding query region
println(hit.identity)       # identity fraction

Protein alignment

let p1 = protein"MKTAYIAKQRQISFVKSHFSRQLEERLGLIEVQAP"
let p2 = protein"MKTAYIAKQRQISFVKSHFSRQLEERLGLIEVQA"

let aln = align(p1, p2, "global", 2, -1, -2)
println("Identity:", aln.identity)
println("CIGAR:   ", aln.cigar)

Batch pairwise alignment

# Align all FASTA records against a single reference
let ref_seq  = dna"ATCGATCGATCGATCG"
let queries  = fasta("queries.fa") |> map(|r| r.seq)

let hits = queries |> map(|q| {
    let r = align(ref_seq, q, "local")
    { score: r.score, identity: r.identity, cigar: r.cigar }
})

hits |> sort(|a, b| b.score - a.score) |> head(10) |> each(println)

consensus(seqs)

Build a plurality-vote consensus string from a list of already-aligned sequences.

let aligned = [
    "ACGT-ACGT",
    "ACGTTACGT",
    "ACGT-ACGG",
]

let result = consensus(aligned)
println(result)                # "ACGTTACGT" (plurality vote)