Motif Finding & Scanning

Motifs are short, recurring sequence patterns that often correspond to functional elements such as transcription factor binding sites, splice signals, or regulatory elements. BioLang provides tools for motif discovery, position weight matrix (PWM) construction, and genome-wide motif scanning.

motif_find

Search for exact or degenerate motif patterns in a sequence. Uses IUPAC ambiguity codes for flexible matching:

# Exact motif search
let seq = dna"ATCGATCGAAATCGATCG"
let hits = motif_find(seq, "ATCG")
# Returns a list of match records with start, end positions

print(len(hits), "matches found")

# IUPAC degenerate motif
# R = A|G, Y = C|T, N = any, W = A|T, S = C|G
hits = motif_find(seq, "RTCG")   # Matches ATCG or GTCG

# Count occurrences of a motif pattern
let n = motif_count(seq, "ATCG")
print(n, "occurrences")

motif_count

Count occurrences of an IUPAC motif pattern in a sequence:

let seq = dna"ATCGATCGATCGATCG"

# Count exact matches
let n = motif_count(seq, "ATCG")
print(n, "occurrences")

# IUPAC degenerate patterns
n = motif_count(seq, "RTCG")   # R = A or G
print(n, "degenerate matches")

# Scan a FASTA file for a regulatory element
read_fasta("data/sequences.fasta")
  |> map(|r| {
    gene: r.id,
    tata_boxes: motif_count(r.seq, "TATAWAW")
  })
  |> to_table()
  |> filter(|r| r.tata_boxes > 0)
  |> print()

consensus

Derive a consensus sequence from a collection of aligned sequences. Each position uses the most frequent base:

# Build consensus from aligned sequences
let aligned = [
  dna"ATCGATCG",
  dna"ATCAATCG",
  dna"ATCGATCG",
  dna"ATCGATCA",
  dna"ATCGATCG",
]

let cons = consensus(aligned)
print(cons)   # => consensus sequence string

pwm (Position Weight Matrix)

Build a position weight matrix from a set of aligned sequences. The PWM captures the position-specific base preferences of a motif:

# Build a PWM from binding site sequences
let sites = [
  dna"ATCGATCG",
  dna"ATCAATCG",
  dna"ATCGATCG",
  dna"ATCGATCA",
  dna"ATCGATCG",
]

let matrix = pwm(sites)

# The PWM is a list of per-position frequency records {A, C, G, T}
matrix |> map(|pos| print(pos))
# Each entry shows base frequencies at that position

pwm_scan

Scan a sequence with a PWM to find high-scoring matches. Returns positions where the PWM score exceeds a threshold:

# Build a PWM from known binding sites
let binding_sites = read_fasta("data/sequences.fasta")
  |> map(|r| r.seq)
  |> collect()

let matrix = pwm(binding_sites)

# Scan a promoter sequence — pwm_scan(seq, pwm, threshold?)
let promoter = dna"ATCGATCGATCGATCGATCGATCGATCGATCGATCGATCG"

# With default threshold
let hits = pwm_scan(promoter, matrix)

# With custom threshold (third argument)
hits = pwm_scan(promoter, matrix, 0.8)

hits |> map(|h| print(h.pos, h.score))

# Genome-wide scan
read_fasta("data/sequences.fasta")
  |> map(|chr| {
    chrom: chr.id,
    sites: pwm_scan(chr.seq, matrix, 0.85)
  })
  |> filter(|r| len(r.sites) > 0)
  |> map(|r| r.sites |> map(|s| {
    chrom: r.chrom, pos: s.pos, score: s.score
  }))
  |> flatten()
  |> to_table()
  |> print()

Combining Motif Functions

The motif functions compose naturally to build analysis workflows:

# Find motifs, build a PWM, then scan other sequences
let sites = [
  dna"TATAAATA",
  dna"TATAATTA",
  dna"TATAAATA",
  dna"TATAAGTA",
]

let matrix = pwm(sites)
let cons = consensus(sites)
print("Consensus:", cons)

# Scan a new sequence
let target = dna"AAAAATATAAATAGGGGTATAAATACCCC"
let hits = pwm_scan(target, matrix, 0.7)
hits |> map(|h| print("Hit at", h.pos, "score:", h.score))

Practical Example: Transcription Factor Binding Analysis

# Analyze motif enrichment in sequences
let seqs = read_fasta("data/sequences.fasta") |> collect()

# Count occurrences of a known motif across sequences
seqs
  |> map(|r| {
    id: r.id,
    hits: motif_count(r.seq, "TATAWAW")
  })
  |> to_table()
  |> filter(|r| r.hits > 0)
  |> print()

# Build a PWM from the sequences and scan for matches
let matrix = seqs |> map(|r| r.seq) |> pwm()
print("PWM built from", len(seqs), "sequences")

# Scan a target sequence
let target_seq = read_fasta("data/sequences.fasta") |> first()
let hits = pwm_scan(target_seq.seq, matrix, 0.8)
print(len(hits), "PWM hits found")

find_pattern

find_pattern is the enhanced IUPAC-aware pattern search that scans both strands and supports mismatch tolerance. It supersedes find_motif for most use cases, returning rich records rather than bare integer positions.

Return type: List<{pos, strand, matched, mismatches}>

pos — 0-based start position on the forward strand
strand — "+" or "-"
matched — the actual substring that matched (reverse-complemented for "-" hits)
mismatches — number of mismatches relative to the pattern (0 for exact matches)

Simple exact search

let seq = dna"ATCGATCGTTCGATCG"

# Exact match on both strands, max_mismatches defaults to 0
let hits = find_pattern(seq, "ATCG")

hits |> map(|h| print(h.pos, h.strand, h.matched))
# 0  +  ATCG
# 4  -  ATCG   (reverse-complement hit)
# 8  +  ATCG
# ...

Mismatch-tolerant search

let seq = dna"ATCGATCGTTCGATCG"

# Allow up to 1 mismatch
let hits = find_pattern(seq, "ATCG", 1)

hits |> map(|h| {
  pos:        h.pos,
  strand:     h.strand,
  matched:    h.matched,
  mismatches: h.mismatches
}) |> to_table() |> print()

Biological example: transcription factor binding sites with a degenerate consensus

# NF-κB canonical site: GGGRNNTCC
# R = A|G, N = any base
let promoter = dna"TTAGGGAAATCCAAAGGGCGGTCCTTAG"

# Exact IUPAC match
let nfkb_hits = find_pattern(promoter, "GGGRNNTCC")
print(len(nfkb_hits), "canonical NF-κB sites found")

# Relaxed search — tolerate a single mismatch
let nfkb_fuzzy = find_pattern(promoter, "GGGRNNTCC", 1)
print(len(nfkb_fuzzy), "sites found with ≤1 mismatch")

nfkb_fuzzy |> map(|h| {
  position:   h.pos,
  strand:     h.strand,
  sequence:   h.matched,
  mismatches: h.mismatches
}) |> to_table() |> print()

Difference from find_motif

find_motif returns a flat list of integer start positions and searches the forward strand only. find_pattern returns structured records that include strand, the actual matched substring, and a mismatch count, and it searches both strands automatically. Use find_motif when you only need forward-strand positions; use find_pattern when you need dual-strand coverage, IUPAC-tolerant fuzzy matching, or richer result metadata.

let seq = dna"ATCGATCGTTCGATCG"

# find_motif: forward strand only, returns List<Int>
let positions = find_motif(seq, "ATCG")
# [0, 4, 12]  — integer positions

# find_pattern: both strands, returns List<{pos, strand, matched, mismatches}>
let records = find_pattern(seq, "ATCG")
# [{pos:0, strand:"+", matched:"ATCG", mismatches:0},
#  {pos:4, strand:"-", matched:"ATCG", mismatches:0},
#  ...]