Python Bio Generators

Generators for Bioinformatics

When to Use

• Parsing a FASTA/FASTQ file too large to fit in memory as a list of records.
• Chaining several sequence transforms (quality filter -> N filter -> translate) where you don't want an intermediate list at every stage.
• Scanning a genome/sequence with a sliding window (GC content, k-mers, motif positions) without materializing every window up front.
• Hitting a MemoryError, slow startup, or high RSS on a script that reads a whole FASTA/FASTQ into a list before processing.
• Finding ORFs, k-mers, or homopolymer runs where results should stream out as found rather than accumulate.

Version Compatibility

Pure Python standard library only (itertools, generator/yield syntax). Works unchanged on Python >=3.7; f-strings and walrus-free code shown here run on Python >=3.8. No third-party packages required.

Prerequisites

• Comfortable with Python functions, for loops, and basic file I/O.
• Familiarity with FASTA/FASTQ format (see python-bio-sequences) helps but isn't required — the readers below are self-contained.

Goal: Stream FASTA/FASTQ records and derived sequence data instead of loading them into lists.

Approach: Every stage is a generator function (yield) or generator expression; nothing runs until the final stage is consumed (e.g., by list(), a for loop, or sum()). Each function takes an iterable in and yields items out, so stages compose into a pipeline with O(1) memory per stage.

def read_fasta(filename):
    """Yield (header, sequence) tuples one record at a time. Memory: O(single record)."""
    header, parts = None, []
    with open(filename) as f:
        for line in f:
            line = line.strip()
            if not line:
                continue
            if line.startswith('>'):
                if header is not None:
                    yield header, ''.join(parts)
                header, parts = line[1:], []
            else:
                parts.append(line)
    if header is not None:
        yield header, ''.join(parts)


def read_fastq(filename):
    """Yield one FASTQ record dict at a time. Memory: O(single record)."""
    with open(filename) as f:
        while True:
            header = f.readline().strip()
            if not header:
                break
            seq = f.readline().strip()
            f.readline()  # '+' separator line
            qual = f.readline().strip()
            if header.startswith('@'):
                yield {
                    'id': header[1:].split()[0],
                    'sequence': seq,
                    'quality': qual,
                    'avg_qual': sum(ord(c) - 33 for c in qual) / len(qual) if qual else 0,
                }


def filter_quality(records, min_avg_qual=30):
    """Pass through only records whose average Phred quality meets the threshold."""
    for rec in records:
        if rec['avg_qual'] >= min_avg_qual:
            yield rec


def filter_no_n(records):
    """Pass through only records with no ambiguous 'N' bases."""
    for rec in records:
        if 'N' not in rec['sequence'].upper():
            yield rec


# Pipeline: nothing is read from disk or computed until list() consumes it.
passing = list(filter_no_n(filter_quality(read_fastq('reads.fastq'))))

Goal: Scan a sequence in overlapping windows and translate DNA to protein lazily, without building intermediate lists at each pipeline stage.

Approach: A sliding-window generator yields (position, subsequence) pairs; a codon/translate pipeline chains codons -> translate -> filter generators so a stop codon can short-circuit the whole pipeline via return inside a generator (raises StopIteration).

CODON_TABLE = {
    'TTT': 'F', 'TTC': 'F', 'TAA': '*', 'TAG': '*', 'TGA': '*',
    'ATG': 'M', 'GCT': 'A', 'CGT': 'R', 'AAT': 'N', 'GAT': 'D',
    # ... full 64-codon table in practice; use Bio.Seq.translate for production code
}


def sliding_window(sequence, window_size, step=1):
    """Yield (start_index, window) for every overlapping window."""
    for i in range(0, len(sequence) - window_size + 1, step):
        yield i, sequence[i:i + window_size]


def gc_content(seq):
    s = seq.upper()
    return (s.count('G') + s.count('C')) / len(s) if s else 0.0


def codons(sequence):
    """Yield successive non-overlapping 3-mers (reading frame 0)."""
    for i in range(0, len(sequence) - 2, 3):
        yield sequence[i:i + 3]


def translate(codon_gen, table=CODON_TABLE):
    """Yield amino acids; stop (raise StopIteration) at the first stop codon."""
    for codon in codon_gen:
        aa = table.get(codon.upper(), 'X')
        if aa == '*':
            return
        yield aa


dna = "GCGCGCATATATATGCGCGCATATATATGCGCGCGC"

# High-GC windows found without ever storing the full window list.
high_gc = sum(1 for _, w in sliding_window(dna, window_size=8, step=4) if gc_content(w) > 0.6)

# Lazy translate pipeline: nothing runs until list() below.
protein = list(translate(codons("ATGAAGCGCGATGAAATCGATGAAGTGGTTGAAATCGAA")))

Goal: Generate k-mers across multiple k values, or process results with itertools, without allocating them all up front.

Approach: yield from delegates to a sub-generator per k; itertools.takewhile/groupby/product/chain provide lazy building blocks for common streaming patterns.

import itertools


def all_kmers_multiK(sequence, k_values):
    """Yield every k-mer for each k in k_values, in order, via delegation."""
    for k in k_values:
        yield from (sequence[i:i + k] for i in range(len(sequence) - k + 1))


# takewhile: stop consuming quality scores as soon as one drops below threshold
quality_scores = [40, 39, 38, 35, 30, 20, 10]
high_quality = list(itertools.takewhile(lambda q: q >= 30, quality_scores))

# groupby: detect homopolymer runs (groups only consecutive equal elements)
dna_run = "AAATTTGGGGCCAATTGCCCCCC"
runs = [(base, sum(1 for _ in group)) for base, group in itertools.groupby(dna_run)]

# product: all possible k-mers of length k over the DNA alphabet (4**k of them)
k = 3
all_kmers = (''.join(c) for c in itertools.product('ATGC', repeat=k))

# chain: splice exons into one mRNA sequence without concatenating strings first
mrna = ''.join(itertools.chain("ATGAAAGCC", "TTTGGGTGA"))

Pitfalls

• Generators are single-pass. Once exhausted (drained by a for loop or list()), re-iterating yields nothing — no error, just empty output. Call list(gen) once and reuse the list, or re-call the generator function to restart.
• return inside a generator stops iteration, it doesn't return a value to the caller — it raises StopIteration internally. Don't rely on a generator's return <value> being retrievable except via .send()/StopIteration.value internals.
• itertools.groupby only groups consecutive equal items — sort first if you need global grouping, or you'll silently get one group per run instead of one group per distinct value.
• Debugging is harder: a traceback inside a generator points to the yield line, not where the value is eventually consumed. Materialize with list(gen) temporarily when debugging a pipeline stage.
• sys.getsizeof() on a generator object is tiny regardless of what it will produce — don't use it to estimate the size of the eventual output.

See Also

• python-bio-iterators — the iterator protocol (__iter__/__next__) that generators implement under the hood.
• python-bio-sequences — Biopython Seq/SeqIO FASTA/FASTQ parsing and translation (production-grade alternative to the hand-rolled readers here).
• python-bio-file-operations — file handling and context-manager basics used inside the streaming readers.
• bio-core-sequence-motifs — motif/k-mer finding with regex, complementary to the generator-based scans here.