Edge Cases & Gotchas

Bioinformatics data is messy. This page covers common pitfalls, edge cases, and how BioLang handles them. Understanding these will save you hours of debugging when working with real-world data.

Empty Sequences

Zero-length sequences in FASTA

# Some FASTA files contain empty sequences. BioLang reads them as empty strings.
let records = read_fasta("data/sequences.fasta")

# This will NOT crash, but may give unexpected results
let lengths = records |> map(|r| len(r.seq))
println(f"Min length: {lengths |> min}")  # Could be 0

# Always filter empty sequences before analysis
let valid = records |> filter(|r| len(r.seq) > 0)
println(f"Sequences with data: {valid |> len}")

# GC content of an empty sequence is defined as 0.0
let empty = dna""
println(empty |> gc_content)  # 0.0

Empty FASTQ records

# FASTQ files with zero-length reads (can happen after trimming)
let reads = read_fastq("data/reads.fastq")

# mean_phred on an empty quality value returns 0.0
println(mean_phred(qual""))  # 0.0

# Filter empty reads when they should not contribute to the analysis
let safe_reads = reads |> filter(|r| len(r.seq) > 0)
let avg_qual = safe_reads
  |> map(|r| mean_phred(r.quality))
  |> mean
println(f"Average quality: {avg_qual |> round(1)}")

Nil Handling

Nil must be handled explicitly

# nil is a value; functions do not silently propagate it
println(to_string(nil))       # "nil"
println(len(to_string(nil)))  # 3

# Numeric aggregations reject nil values, so filter them first
let values = [1.5, nil, 2.5]
let present = values |> filter(|x| x != nil)
println(present |> sum)   # 4.0
println(present |> mean)  # 2.0

# Invalid conversions and missing record keys raise errors
# float("")              # TypeError
# {score: 1.0}["pvalue"] # KeyError

Nil in conditionals

# nil is falsy, but not the same as false
let x = nil

if x { println("truthy") } else { println("falsy") }  # prints "falsy"
if x == false { println("equal") } else { println("not equal") }  # prints "not equal"
if x == nil { println("is nil") }  # prints "is nil"

# Use ?? (nil coalescing) for defaults
let name = nil
println(f"Hello, {name ?? "unknown"}")  # "Hello, unknown"

# Use ?. (nil-safe access) for chained lookups
let record = { info: nil }
println(record.info?.gene)  # nil, no crash
# record.info.gene        # ERROR: cannot access field on nil

Large Files

Streaming vs collecting

# Eager: read_fastq returns a reusable in-memory Table
let all_reads = read_fastq("data/reads.fastq")
# This can exhaust memory for large files!

# Lazy: fastq returns a one-pass Stream
let passing_count = fastq("data/reads.fastq")
  |> filter(|r| mean_phred(r.quality) >= 30.0)
  |> count()
println(f"Passing reads: {passing_count}")

# If you need to collect, sample first
let sample = fastq("data/reads.fastq") |> take(10000) |> collect
let avg_len = sample |> map(|r| len(r.seq)) |> mean
println(f"Estimated average read length: {avg_len}")

Memory-safe aggregation

Requires CLI: this operation reads directly from the local filesystem.

# read_stats scans the path without building an in-memory read table
let stats = read_stats("data/reads.fastq")

# read_stats computes quality and length metrics in constant memory
println(f"Total reads: {stats.total_reads}")
println(f"Mean length: {stats.mean_length |> round(1)}")
println(f"Mean quality: {stats.mean_quality |> round(1)}")

Encoding Issues

Quality score encoding

# BioLang reads modern Sanger/Illumina 1.8+ FASTQ as Phred+33
let reads = read_fastq("data/reads.fastq")
let first_read = reads |> first
println(f"Mean Phred: {mean_phred(first_read.quality) |> round(1)}")

# Legacy Illumina 1.3-1.7 files may use Phred+64. Convert those files to
# Phred+33 before reading them; the reader does not auto-detect +64.

Sequence characters

# DNA sequences may contain ambiguity codes (IUPAC)
let seq = dna"ATCGNNRYSWKM"

# Standard functions handle ambiguity codes
println(seq |> len)              # 12
println(seq |> gc_content)       # Counts only definite G/C bases
println(contains(to_string(seq), "N"))    # true

# Reverse complement preserves ambiguity codes
println(seq |> reverse_complement)  # MKWSYRNNCGAT

# But k-mer counting may give unexpected results
let kmer_counts = seq |> kmers(3) |> frequencies
# K-mers containing N are valid k-mers in BioLang
println(kmer_counts)  # includes "ATC", "TCG", "CGN", "GNN", etc.

# Filter out ambiguous k-mers if needed
let clean_kmers = seq |> kmers(3) |> filter(|k| !contains(str(k), "N"))

Numeric Precision

Floating point comparison

# Classic floating point trap
let a = 0.1 + 0.2
println(a == 0.3)          # false!
println(a)                 # 0.30000000000000004

# Use abs() for floating point comparison
println(abs(a - 0.3) < 1e-10)  # true

# This matters for p-value filtering
let pval = 0.05
let computed = 0.01 + 0.04
# WRONG: filter(|v| { v.pval == 0.05 })
# RIGHT: filter(|v| { v.pval <= 0.05 })

Integer overflow

# Genome coordinates can be large
let pos = 2147483647  # Max 32-bit int

# BioLang uses 64-bit integers by default, so this is fine:
let big_pos = pos + 1000
println(big_pos)  # 2147484647

# But be careful with multiplication
let genome_size = 3_000_000_000
let coverage = 30
let total_bases = genome_size * coverage  # 90 billion -- fits in i64
println(f"Total bases: {total_bases}")

File Format Gotchas

VCF multi-allelic sites

# Multi-allelic VCF records have comma-separated ALT alleles
let vcf = read_vcf("variants.vcf")

for v in vcf |> take(5) {
  if v.alt |> contains(",") {
    println(f"Multi-allelic: {v.chrom}:{v.pos} {v.ref} => {v.alt}")
    # Split ALT on comma to get individual alleles
    for allele in v.alt |> split(",") {
      println(f"  Allele: {allele}")
    }
  }
}

# Common mistake: treating ALT as a single string
# WRONG: filter(|v| { v.alt == "A" })
# RIGHT: filter(|v| { v.alt |> split(",") |> contains("A") })

BED coordinate system

# BED is 0-based, half-open: [start, end)
# VCF/GFF are 1-based, closed: [start, end]
# This off-by-one difference causes many bugs

let bed_start = 100   # First base is position 100
let bed_end = 200     # Last base is position 199
let bed_length = bed_end - bed_start  # 100 bases

# Converting BED to 1-based:
let one_based_start = bed_start + 1  # 101
let one_based_end = bed_end          # 200

# interval() stores numeric half-open coordinates; convert before constructing it
let bed_region = interval("chr1", 100, 200)
let vcf_as_bed = interval("chr1", 100, 200)  # VCF position 101 converted to BED
let tree = interval_tree([bed_region])
println(len(query_overlaps(tree, "chr1", 100, 200)) > 0)  # true

Windows line endings

# Files from Windows have \r\n line endings
# BioLang handles this automatically for standard formats (FASTQ, VCF, etc.)

# But for plain text, trailing \r can cause problems
let lines = read_lines("windows_file.txt")
let first = lines |> first
println(first |> ends_with("gene"))      # might be false!
println(first |> trim |> ends_with("gene"))  # true

# read_lines automatically strips \r, but read_text does not
let raw = read_text("windows_file.txt")
let clean = raw |> replace("\r\n", "\n")