Skip to content

Instantly share code, notes, and snippets.

@znorgaard
Created June 9, 2026 17:02
Show Gist options
  • Select an option

  • Save znorgaard/7cd1494d492f57aea14e9c22953afbef to your computer and use it in GitHub Desktop.

Select an option

Save znorgaard/7cd1494d492f57aea14e9c22953afbef to your computer and use it in GitHub Desktop.
bioinformatics-python — Claude Code skill for Python bioinformatics (fgpyo, pybedlite, fgmetric). See install.sh to install into ~/.claude/skills/.
#!/usr/bin/env sh
# Install the `bioinformatics-python` Claude Code skill from this gist.
#
# Gists are flat, so the skill's reference docs are stored with a
# `references__` prefix. This script clones the gist and rebuilds the
# real directory layout:
#
# <dest>/
# SKILL.md
# references/
# fastx.md fgmetric.md io.md pybedlite.md sam.md util.md vcf.md
#
# Usage:
# sh install.sh # installs to ~/.claude/skills/bioinformatics-python
# sh install.sh <dest> # installs to a custom directory
#
# One-liner (requires git):
# curl -fsSL https://gist.githubusercontent.com/znorgaard/7cd1494d492f57aea14e9c22953afbef/raw/install.sh | sh
set -eu
GIST_ID="7cd1494d492f57aea14e9c22953afbef"
DEST="${1:-$HOME/.claude/skills/bioinformatics-python}"
command -v git >/dev/null 2>&1 || { echo "error: git is required" >&2; exit 1; }
TMP="$(mktemp -d)"
trap 'rm -rf "$TMP"' EXIT INT TERM
echo "Cloning gist $GIST_ID ..."
git clone --depth 1 "https://gist.github.com/${GIST_ID}.git" "$TMP/gist" >/dev/null 2>&1
mkdir -p "$DEST/references"
cp "$TMP/gist/SKILL.md" "$DEST/SKILL.md"
for f in "$TMP/gist"/references__*.md; do
[ -e "$f" ] || continue
base="$(basename "$f")"
cp "$f" "$DEST/references/${base#references__}"
done
echo "Installed bioinformatics-python skill to:"
echo " $DEST"
echo "Restart Claude Code (or start a new session) to pick it up."

fgpyo.fasta and fgpyo.fastx reference

FastaBuilder — fluent FASTA test fixture

fgpyo.fasta.FastaBuilder. Use this in tests instead of writing FASTA strings to disk by hand. Produces .fasta, .fai, and .dict files in one call.

from pathlib import Path
from fgpyo.fasta import FastaBuilder

builder = FastaBuilder()

# `.add(name)` returns a ContigBuilder; `.add(bases, repeats)` appends to it.
builder.add("chr1").add("ACGT", 100)              # 400 bp of ACGT-repeats
builder.add("chr2").add("N", 1000).add("CGCG", 50)
builder.to_file(path=Path("ref.fasta"))            # writes .fasta + .fasta.fai + .dict

ContigBuilder API:

contig = builder.add("chrom_name")
contig.add(bases: str, repeats: int = 1) -> ContigBuilder    # chainable
contig.bases                                                  # accumulated sequence (str)

builder.to_file(path) produces all three sidecar files; subsequent pysam.FastaFile(path) opens cleanly.

SequenceDictionary — parse / construct .dict files

fgpyo.fasta.sequence_dictionary.SequenceDictionary. Use this when working with sequence dictionaries (@SQ lines) instead of poking at SAM headers directly.

from fgpyo.fasta.sequence_dictionary import SequenceDictionary, Keys

sd = SequenceDictionary.from_sam(sam_header)             # from a pysam AlignmentHeader
sd = SequenceDictionary.from_fasta(Path("ref.fasta"))    # via the .dict / .fai sidecar

# Lookup by index, name, or alias
entry = sd[0]
entry = sd["chr1"]
entry[Keys.MD5]                                         # any standard SAM @SQ field
entry[Keys.ASSEMBLY] = "GRCh38"
del entry[Keys.SPECIES]

# Round-trip back to a SAM header
header = sd.to_sam_header()

Keys enum contains the standard @SQ fields (SN, LN, M5, AS, SP, UR, etc.) plus their long-name aliases.

There's also a sequence_dictionary(fasta_path) helper that loads (or generates) the .dict and returns a list of dicts — fine for ad-hoc use, but SequenceDictionary is the typed entry point.

FastxZipped — zipped iteration over paired FASTX streams

fgpyo.fastx.FastxZipped. Replaces ad-hoc zip() over two pysam.FastxFiles, with mismatched-name detection and proper cleanup. Works with FASTA and FASTQ, and any number of input streams.

from fgpyo.fastx import FastxZipped

with FastxZipped("R1.fastq.gz", "R2.fastq.gz") as zipped:
    for r1, r2 in zipped:                       # both pysam.FastxRecord
        assert r1.name == r2.name               # FastxZipped enforces this
        seq1, seq2 = r1.sequence, r2.sequence
        qual1, qual2 = r1.quality, r2.quality   # str of phred+33 chars (FASTQ only)

Constructor:

FastxZipped(*paths: Path | str, persist: bool = False)
  • persist=False (default) — memory-efficient. Each yielded FastxRecord is invalidated on the next iteration; copy out anything you need to keep.
  • persist=True — yielded records remain valid; uses more memory.

Three-or-more-way zips work the same way (e.g. R1/R2/I1 of an Illumina lane):

with FastxZipped("R1.fq.gz", "R2.fq.gz", "I1.fq.gz") as z:
    for r1, r2, i1 in z:
        ...

Common mistakes

  • Calling FastaBuilder.to_file(path).fasta etc. to_file returns the path written; the indexes are siblings (<path>.fai, <path>.replace(.fasta, .dict)).
  • Iterating a FastxZipped without a with block. It opens lazily in __enter__; using outside a context manager skips proper cleanup.
  • Storing yielded FastxRecord instances when persist=False. They get invalidated on the next iteration. Copy .name / .sequence / .quality out, or pass persist=True.
  • Reaching for raw pysam.FastaFile to mock a small reference. FastaBuilder produces fully-indexed FASTA in 3 lines and is the convention in tests across the codebase.

fgmetric reference

fgmetric (PyPI: fgmetric, on GitHub at fg-labs/fgmetric) is a Pydantic-backed Metric class for typed delimited files. Prefer fgmetric.Metric over fgpyo.util.Metric — fgmetric's validation, error messages, list/Counter handling, and Pydantic ecosystem (custom validators / serializers, aliases) are stronger.

Requires Python ≥ 3.12. Pydantic ≥ 2.11.4 is pulled in transitively.

Subclass and go

from pathlib import Path
from fgmetric import Metric, MetricWriter


class AlignmentMetric(Metric):
    read_name: str
    mapping_quality: int
    is_duplicate: bool = False


# Write — context manager only; header is emitted on enter.
with MetricWriter(AlignmentMetric, Path("out.tsv")) as w:
    w.write(AlignmentMetric(read_name="r1", mapping_quality=60))
    w.write(AlignmentMetric(read_name="r2", mapping_quality=30, is_duplicate=True))

# Read — streaming iterator. Strings are coerced to int/float/bool by Pydantic.
for m in AlignmentMetric.read(Path("out.tsv")):
    ...

No @attr.s, no @dataclassMetric inherits from Pydantic BaseModel, so fields are declared the Pydantic way. Field(alias=...), custom validators, and @field_serializer all work.

API signatures

class Metric(BaseModel, ABC):
    collection_delimiter: ClassVar[str] = ","         # used for list[T] fields

    @classmethod
    def read(
        cls,
        path: Path,
        delimiter: str = "\t",
        fieldnames: Sequence[str] | None = None,      # supply for headerless files
    ) -> Iterator[Self]: ...
class MetricWriter[T: Metric](AbstractContextManager):
    def __init__(
        self,
        metric_class: type[T],
        filename: Path | str,
        delimiter: str = "\t",
        lineterminator: str = "\n",
    ) -> None: ...
    def write(self, metric: T) -> None: ...          # one row
    def writeall(self, metrics: Iterable[T]) -> None: # iterable of rows
    def close(self) -> None: ...
  • File handle is opened in __init__; the header is written on enter. Always use MetricWriter as a context manager.
  • Reads use csv.DictReader; the file is opened with utf-8-sig so a BOM in the header is tolerated.
  • Writes use csv.DictWriter with model_dump(mode="json"), so enums/datetimes serialize as JSON-compatible values.

Headerless files

Pass fieldnames= to .read() for a file with no header row:

for m in AlignmentMetric.read(
    Path("out.tsv"),
    fieldnames=["read_name", "mapping_quality", "is_duplicate"],
):
    ...

(Order must match the file's column order. With a header present, omit fieldnames= and let validation surface mismatches.)

CSV (or any single-char delimiter)

for m in MyMetric.read(Path("data.csv"), delimiter=","):
    ...

with MetricWriter(MyMetric, Path("out.csv"), delimiter=",") as w:
    w.writeall(rows)

list[T] fields auto-split / join

Fields annotated as list[T] are read from a delimited substring and serialized back the same way. The within-field delimiter defaults to , and is controlled by the collection_delimiter class variable (must be a single character).

class TaggedRead(Metric):
    read_id: str
    tags: list[str]                  # "A,B,C" ↔ ["A", "B", "C"]
    scores: list[int]                # "1,2,3" ↔ [1, 2, 3]
    optional_tags: list[str] | None  # ""      ↔ None
    nullable_elements: list[int | None]  # "1,,3" ↔ [1, None, 3]


class SemicolonMetric(Metric):
    collection_delimiter = ";"
    values: list[int]                # "1;2;3" ↔ [1, 2, 3]

Round-trips are lossy if list elements themselves contain the delimiter — pick a delimiter your data never contains.

Counter[StrEnum] fields pivot wide ↔ long

A single field typed Counter[E] (where E: StrEnum) is pivoted out to one column per enum member on write, and rolled back up into a Counter on read. Missing members default to count 0. No CounterPivotTable mixin is neededMetric already inherits it.

from collections import Counter
from enum import StrEnum
from fgmetric import Metric


class Base(StrEnum):
    A = "A"
    C = "C"
    G = "G"
    T = "T"


class BaseCountMetric(Metric):
    position: int
    counts: Counter[Base]


# Input TSV:
#   position  A   C   G   T
#   1         10  5   3   2
#
# Reads as: BaseCountMetric(position=1,
#                           counts=Counter({Base.A: 10, Base.C: 5, Base.G: 3, Base.T: 2}))
# Writes back to the same wide layout.

Constraints (enforced at class-init):

  • At most one Counter field per model.
  • The type parameter must be a StrEnum subclass (not int, not a Literal, not a non-Str enum).
  • Counter[E] | None is not allowed — required field.

Optional fields, empty strings, defaults

  • Empty string in an input column whose field annotation is T | None is converted to None automatically.
  • Empty string in a required field surfaces as a pydantic.ValidationError (with field name + reason).
  • Defaults work as in Pydantic — count: int = 0 etc. — and missing optional columns fall back to those defaults on read.

Common mistakes

  • MyMetric.write(path, *rows). fgmetric has no class-level write method. Use MetricWriter as a context manager. (This is the most common drift from fgpyo's Metric.)
  • MetricWriter(MyMetric, Path("out.tsv.gz")) expecting auto-gzip. It will write uncompressed text into a file with a .gz suffix. fgmetric does not detect or apply compression. If you need gzip output: write to a .tsv temp path and gzip after, or open gzip.open(...)/fgpyo.io.to_writer(...) yourself and write rows via the underlying csv.DictWriter (not via MetricWriter).
  • @attr.s / @dataclass on the subclass. Metric is a Pydantic BaseModel — declare fields directly with type annotations and (optionally) Field(...). Adding attrs decorators breaks Pydantic's field discovery.
  • MyMetric(name="x", count="3")-style coercion test in code. Pydantic does coerce, but lean on .read() (which is fed strings by csv.DictReader) rather than constructing models from string-typed Python values yourself.
  • Forgetting with MetricWriter(...) as w:. Header is written on __enter__; using the writer outside a with block (or forgetting close()) can leave the header unwritten or the file open.
  • from fgpyo.util.metric import Metric. That's the legacy import. Use from fgmetric import Metric, MetricWriter.
  • Stacking CounterPivotTable ahead of Metric in the bases list. Unnecessary — Metric already inherits the pivot mixin.

Differences vs fgpyo.util.Metric (quick map)

You're used to (fgpyo) fgmetric equivalent
@attr.s(frozen=True, kw_only=True, auto_attribs=True, slots=True) on the subclass Just class M(Metric): name: str; n: int (Pydantic)
MyMetric.write(path, *rows) No class-level write — with MetricWriter(M, path) as w: w.writeall(rows)
MyMetric.read(path, ignore_extra_fields=True, strip_whitespace=False, threads=N) MyMetric.read(path, delimiter="\t", fieldnames=None) — no threads, no whitespace flag, no ignore-extra
Gzip auto-detect on .gz / .bgz None. Plain text only.
cls._parsers() override for custom field parsing Pydantic validators (@field_validator, @model_validator(mode="before"))
cls.format_value(value) override for serialization Pydantic @field_serializer / @model_serializer
MyMetric.header() No public equivalent. The internal cls._header_fieldnames() is private; if you need the header before writing, instantiate a MetricWriter and the header is written automatically.
MetricWriter(append=True, include_fields=..., exclude_fields=...) Not supported — write a fresh file, or filter the iterable yourself.

Legacy fallback: fgpyo.util.Metric

If a project pins fgmetric out or forbids it, fall back to fgpyo.util.Metric — same conceptual shape, attrs-based, with gzip support and Metric.write(path, *rows) / MetricWriter. The fgpyo Metric pattern lives in the project's existing code; reach for it only when fgmetric is unavailable.

fgpyo.io reference

Auto-gzip-aware text I/O plus path assertions. Default to fgpyo.io instead of open / gzip.open switched on suffix.

For tabular typed records (Metric subclasses), see fgmetric.mdfgmetric.Metric is preferred over fgpyo.util.Metric.

Readers and writers (auto-detect .gz / .bgz)

import fgpyo.io as fio
from pathlib import Path

# Auto-detects .gz / .bgz; returns a TextIOWrapper.
with fio.to_reader(Path("data.tsv.gz")) as r:
    for line in r:
        ...

with fio.to_writer(Path("out.txt.gz"), append=False) as w:
    w.write("hello\n")

# Streaming line iteration (handles .gz transparently)
for line in fio.read_lines(Path("data.txt.gz"), strip=True):
    ...

fio.write_lines(Path("out.txt"), ["row1", "row2"])

Signatures:

to_reader(path: Path, threads: int | None = None) -> TextIOWrapper
to_writer(path: Path, append: bool = False, threads: int | None = None) -> TextIOWrapper
read_lines(path: Path, strip: bool = False, threads: int | None = None) -> Iterator[str]
write_lines(path: Path, lines_to_write: Iterable[Any], append: bool = False,
            threads: int | None = None) -> None

threads=N enables parallel zlib decompression/compression for .gz paths (uses zlib-ng).

Path assertions — replace custom existence checks

fio.assert_path_is_readable(p)                    # exists, is file, is readable
fio.assert_path_is_writable(p, parent_must_exist=True)
fio.assert_directory_exists(p)
fio.assert_fasta_indexed(fasta, dictionary=False, bwa=False)  # checks .fai/.dict/.bwt

(assert_path_is_writeable is a deprecated spelling kept for backwards compatibility — prefer assert_path_is_writable.)

Misc

with fio.suppress_stderr():
    noisy_call()                                  # stderr → /dev/null inside block

Pairing fgpyo.io with fgmetric.MetricWriter for gzip output

fgmetric does not auto-detect .gz paths. If you need gzipped typed output, the simplest path is to write uncompressed and gzip after:

import shutil, gzip
from pathlib import Path
from fgmetric import MetricWriter

tmp = Path("out.tsv")
with MetricWriter(MyMetric, tmp) as w:
    w.writeall(rows)
with tmp.open("rb") as fin, gzip.open("out.tsv.gz", "wb") as fout:
    shutil.copyfileobj(fin, fout)
tmp.unlink()

For plain (non-Metric) gzip text, fgpyo.io.to_writer is the right tool.

pybedlite reference

Lightweight typed BED + interval-overlap library from Fulcrum Genomics. Install: pip install pybedlite. All coordinates 0-based, half-open (BED convention).

Top-level surface

import pybedlite as pybed
from pybedlite import BedSource, BedWriter, BedRecord, BedStrand
from pybedlite.overlap_detector import OverlapDetector, Interval, Span, StrandedSpan

OverlapDetector and Interval are not re-exported at top level — import them from pybedlite.overlap_detector.

BedRecord — typed BED record

pybedlite.BedRecord. Frozen attrs class, kw-only, supports BED3–BED12 (the present field count is auto-detected via bed_field_num).

from pybedlite import BedRecord, BedStrand

rec = BedRecord(
    chrom="chr1", start=1000, end=2000,
    name="gene1", score=500, strand=BedStrand.Positive,
)
rec.bed_field_num                   # 6
rec.as_bed_line()                   # "chr1\t1000\t2000\tgene1\t500\t+"
rec.as_bed_line(number_of_output_fields=3)   # "chr1\t1000\t2000"
rec.refname                          # "chr1"  (alias for chrom — satisfies Span)
rec.negative                         # False  (False if strand is Positive or None)

Full constructor:

BedRecord(*,
    chrom: str, start: int, end: int,
    name: Optional[str] = None,
    score: Optional[int] = None,                     # plain int (not clamped to 0..1000)
    strand: Optional[BedStrand] = None,              # None == "no strand"
    thick_start: Optional[int] = None,
    thick_end: Optional[int] = None,
    item_rgb: Optional[Tuple[int, int, int]] = None,
    block_count: Optional[int] = None,
    block_sizes: Optional[List[int]] = None,
    block_starts: Optional[List[int]] = None,
)

Validators (asserted at construction; AssertionError, not ValueError):

  • end > start
  • thick_start and thick_end are both set or both None
  • block fields all set together; lengths match block_count; first block_start is 0; last block ends at end

BedRecord.from_interval(interval) round-trips an Interval back to a BedRecord. BedStrand has only Positive / Negative — unstranded round-trips through Interval collapse to Positive.

BedStrand

class BedStrand(Enum):
    Positive = "+"
    Negative = "-"
    @property
    def opposite(self) -> BedStrand: ...

BedStrand("+") is BedStrand.Positive    # True
BedStrand.Positive.opposite             # BedStrand.Negative

There is no unstranded variant; use BedRecord.strand = None for that.

BedSource — streaming reader

pybedlite.BedSource (also exposed as pybed.reader(path)). ContextManager[BedSource], Iterable[BedRecord].

import pybedlite as pybed
from pathlib import Path

with pybed.reader(Path("regions.bed")) as src:
    for rec in src:
        print(rec.chrom, rec.start, rec.end, rec.name)
    src.num_fields                      # set after first parse — e.g. 6 for BED6

Behavior:

  • Skips blank lines and lines starting with #, browser, track.
  • Asserts at least 3 columns per row.
  • Field-aware parsing: . becomes None; item_rgb parsed as (r, g, b); block_sizes / block_starts parsed as list[int].
  • Iterating without with auto-opens / auto-closes on exhaustion.

BedPath accepts Path | str | TextIOWrapper.

BedWriter — streaming writer

pybedlite.BedWriter (also pybed.writer(path, num_fields=None)).

with pybed.writer(Path("out.bed"), num_fields=6) as w:
    w.write_all(records)
    w.write(extra_record)
class BedWriter(ContextManager):
    def __init__(self, path: BedPath, num_fields: Optional[int] = None) -> None: ...
    def write(self, record: BedRecord, truncate: bool = False, add_missing: bool = False) -> None: ...
    def write_all(self, records: Iterable[BedRecord], truncate: bool = False, add_missing: bool = False) -> None: ...
  • Leave num_fields=None to inherit from the first record's bed_field_num.
  • Mismatched widths raise ValueError unless truncate=True (record has more fields) or add_missing=True (record has fewer; pads with .).

Interval — lightweight stranded interval

pybedlite.overlap_detector.Interval. Use this when you don't need the BED12 surface area but want something hashable that satisfies StrandedSpan.

from pybedlite.overlap_detector import Interval

iv = Interval("chr1", 100, 200)
iv2 = Interval("chr1", 150, 250, negative=True, name="r2")
iv.overlap(iv2)                                    # 50  (0 if different refname or no overlap)
iv.length()                                        # 100

# UCSC 1-based-closed → 0-based half-open
Interval.from_ucsc("chr1:127140001-127140050(-)")
# Interval(refname='chr1', start=127140000, end=127140050, negative=True)

Interval.from_bedrecord(rec)

Validators (raise ValueError): start >= 0, end > start. Round-tripping through BedRecord.from_interval collapses unstranded to BedStrand.Positive.

OverlapDetector — interval-tree overlap queries

pybedlite.overlap_detector.OverlapDetector. Generic over any object that satisfies the Span (refname/start/end) or StrandedSpan (also negative) protocol — BedRecord and Interval both qualify out of the box.

from pathlib import Path
from pybedlite.overlap_detector import OverlapDetector, Interval

# From a list of typed intervals or BedRecords
detector = OverlapDetector([
    Interval("chr1",   0, 100),
    Interval("chr1", 150, 250, negative=True),
    Interval("chr2", 500, 600),
])
detector.add(Interval("chr1", 200, 400))

# Or directly from a BED file (returns OverlapDetector[BedRecord])
detector = OverlapDetector.from_bed(Path("targets.bed"))

# Queries — `interval` is anything satisfying Span (refname/start/end)
query = Interval("chr1", 90, 210)
detector.overlaps_any(query)              # True
detector.get_overlaps(query)              # list of stored items overlapping the query
detector.get_enclosing_intervals(query)   # items that fully contain the query
detector.get_enclosed(query)              # items fully contained inside the query

len(detector)                             # 4
list(detector)                            # all stored items

Key signatures:

class OverlapDetector(Generic[SpanType], Iterable[SpanType], Sized):
    def __init__(self, intervals: Iterable[SpanType] | None = None) -> None: ...
    def add(self, interval: SpanType) -> None: ...
    def add_all(self, intervals: Iterable[SpanType]) -> None: ...
    def overlaps_any(self, interval: Span) -> bool: ...
    def get_overlaps(self, interval: Span) -> list[SpanType]: ...
    def get_enclosing_intervals(self, interval: Span) -> list[SpanType]: ...
    def get_enclosed(self, interval: Span) -> list[SpanType]: ...
    @classmethod
    def from_bed(cls, path: Path) -> OverlapDetector[BedRecord]: ...

Performance:

  • Indexes are built lazily on first query and invalidated on every add / add_all. Load all intervals first, then query.
  • Backed by superintervals.IntervalSet (one tree per refname).
  • Query results are deduplicated and sorted by (start, end, negative, refname).

Span / StrandedSpan — structural protocols

class Span(Hashable, Protocol):
    @property
    def refname(self) -> str: ...
    @property
    def start(self) -> int: ...
    @property
    def end(self) -> int: ...

class StrandedSpan(Span, Protocol):
    @property
    def negative(self) -> bool: ...

Any hashable user-defined attrs class with refname, start, end (and optionally negative) can be stored in an OverlapDetector directly — no inheritance required. BedRecord.refname is an alias for chrom, so it satisfies Span despite using chrom as its underlying field name.

Common mistakes

  • Importing OverlapDetector from top-level pybedlite. It lives in pybedlite.overlap_detector.
  • Calling the reader BedReader. It's BedSource (or pybedlite.reader(path)).
  • Adding intervals to an OverlapDetector between queries. Indexes rebuild on every add — load all intervals first, then query in a tight loop.
  • Treating BedRecord.score as 0..1000 only. It's a plain int | None; pybedlite does not clamp.
  • Expecting Interval to track unstranded. It only has negative: bool; round-tripping unstranded BedRecordInterval lossy-converts to BedStrand.Positive.
  • Using BedRecord.strand == "+". It's a BedStrand enum (or None); compare with BedStrand.Positive or use rec.negative.
  • Catching ValueError from a malformed BedRecord. Construction validates with assert — that's AssertionError (and is disabled by python -O). Validate inputs before constructing if you need controlled error handling.

fgpyo.sam reference

Module: fgpyo.sam (and submodules fgpyo.sam.builder, fgpyo.sam.clipping).

reader / writer — file-type-aware open

from fgpyo.sam import reader, writer, SamFileType

# Auto-detects from extension (.sam / .bam / .cram); supports stdin/stdout.
with reader(Path("aln.bam")) as bam:
    for rec in bam:
        ...

with writer(Path("out.bam"), header=bam.header) as out:
    out.write(rec)

Signatures:

reader(path: SamPath, file_type: SamFileType | None = None,
       unmapped: bool = False) -> AlignmentFile

writer(path: SamPath, header: dict | AlignmentHeader | str,
       file_type: SamFileType | None = None) -> AlignmentFile

SamFileType enum: SAM, BAM, CRAM. SamFileType.from_path(p) for explicit detection.

SamBuilder — fluent test record construction

fgpyo.sam.builder.SamBuilder. Default: 100bp paired reads, qual 30, mapq 60, HG19 chr1–chrM sequence dictionary, deterministic seed=42. Construct, call .add_pair() / .add_single(), then materialize.

from fgpyo.sam.builder import SamBuilder

builder = SamBuilder(r1_len=100, r2_len=100, base_quality=30,
                     mapping_quality=60, seed=42)

# Same-contig pair
r1, r2 = builder.add_pair(chrom="chr1", start1=1000, start2=1200)

# Different-contig pair (chimeric / translocation)
builder.add_pair(chrom1="chr1", start1=1000, chrom2="chr2", start2=2000)

# Specific cigar / mapq / strand
builder.add_pair(chrom="chr1", start1=1000, start2=1200,
                 cigar1="50M5I45M", cigar2="100M",
                 mapq1=40, mapq2=60, strand1="+", strand2="-")

# Single read (supplementary or secondary go through add_single)
sup = builder.add_single(name="A", chrom="chr1", start=1500, cigar="40M",
                         read_num=1, supplementary=True, mapq=20)
sec = builder.add_single(name="A", chrom="chr3", start=99,
                         read_num=1, secondary=True)

# Materialize
records = builder.to_unsorted_list()                  # in creation order
sorted_recs = builder.to_sorted_list()                 # coordinate sorted
bam_path: Path = builder.to_path()                     # temp coord-sorted BAM (+ index)
bam_path = builder.to_path(path=tmp_path / "in.bam")   # explicit path

SamBuilder constructor

SamBuilder(
    r1_len: int | None = None,
    r2_len: int | None = None,
    base_quality: int = 30,
    mapping_quality: int = 60,
    sd: list[dict] | None = None,        # SAM @SQ records; default = HG19 chr1..chrM
    rg: dict | None = None,              # @RG; default ID=1 SM=1_AAAAAA LB=default PL=ILLUMINA
    seed: int = 42,
    sort_order: SamOrder = SamOrder.Coordinate,
)

SamBuilder.default_sd(), SamBuilder.default_rg(), and builder.header give you the active dict / RG / pysam AlignmentHeader.

add_pair — keyword-only

builder.add_pair(
    *,
    name: str | None = None,                   # auto-named if None
    bases1: str | None = _UNSET,               # _UNSET = random; pass `None` to omit sequence
    bases2: str | None = _UNSET,
    quals1: list[int] | None = _UNSET,         # _UNSET = base_quality; `None` to omit
    quals2: list[int] | None = _UNSET,
    chrom: str | None = None,                  # both reads same contig
    chrom1: str | None = None,                 # OR: separate per-read contigs
    chrom2: str | None = None,
    start1: int = sam.NO_REF_POS,
    start2: int = sam.NO_REF_POS,
    cigar1: str | None = None,                 # default = all-M of read length
    cigar2: str | None = None,
    mapq1: int | None = None,
    mapq2: int | None = None,
    strand1: str = "+",                        # default FR orientation
    strand2: str = "-",
    attrs: dict[str, Any] | None = None,       # extra SAM tags applied to both records
) -> tuple[AlignedSegment, AlignedSegment]

add_pair does not take secondary / supplementary kwargs — use add_single(read_num=1, supplementary=True, ...) / read_num=2 for those records and reuse the same name= to attach them to a template.

add_single — keyword-only

builder.add_single(
    *,
    name: str | None = None,
    read_num: int | None = None,               # None = unpaired; 1 or 2 for paired
    bases: str | None = _UNSET,
    quals: list[int] | None = _UNSET,
    chrom: str = sam.NO_REF_NAME,              # "*" = unmapped
    start: int = sam.NO_REF_POS,               # -1 = unmapped
    cigar: str | None = None,
    mapq: int | None = None,
    strand: str = "+",
    secondary: bool = False,
    supplementary: bool = False,
    attrs: dict[str, Any] | None = None,
) -> AlignedSegment

to_path

builder.to_path(
    path: Path | None = None,                  # tempfile if None
    index: bool = True,                        # auto-index BAM/CRAM
    pred: Callable[[AlignedSegment], bool] = lambda _: True,
    tmp_file_type: SamFileType | None = None,
) -> Path

Template / TemplateIterator — group records by template

fgpyo.sam.Template. Use this whenever you need r1+r2+supplementaries+secondaries together. Replaces hand-rolled query-name buffering.

from fgpyo.sam import Template, reader

with reader(Path("queryname-grouped.bam")) as bam:
    for tpl in Template.iterator(bam):
        tpl.name                       # query name
        tpl.r1                         # AlignedSegment | None — primary R1
        tpl.r2                         # AlignedSegment | None — primary R2
        tpl.r1_supplementals           # list[AlignedSegment]
        tpl.r2_supplementals
        tpl.r1_secondaries
        tpl.r2_secondaries

        for rec in tpl.primary_recs():  # r1 then r2 (skipping None)
            ...
        for rec in tpl.all_recs():      # primaries + supps + secondaries
            ...
        for rec in tpl.all_r1s():       # r1 + r1 supps + r1 secondaries
            ...

Other useful Template ops:

tpl.set_mate_info()                    # reset mate flags + MC/MQ across the template
tpl.set_tag("RX", umi)                 # set/remove a tag on every record
tpl.write_to(out_writer, primary_only=False)
tpl.validate()                         # asserts r1/r2 flags / supps consistent

# Build from an explicit list (e.g. inside a test)
tpl = Template.build(records, validate=True)

Template.iterator(...) requires query-name-grouped input (sorted by qname or sorted-by-name groupby). Coord-sorted input will produce broken templates.

Cigar — typed CIGAR ops

fgpyo.sam.Cigar. Replaces hand-counted CIGAR length math.

from fgpyo.sam import Cigar, CigarOp, CigarElement

cig = Cigar.from_cigarstring("5S50M2D5M10S")
cig.length_on_query()                   # 70  (5+50+5+10) — METHOD, needs parens
cig.length_on_target()                  # 57  (50+2+5)    — METHOD, needs parens
cig.elements                            # tuple[CigarElement, ...]
cig.reversed()                          # Cigar with element order reversed
cig.coalesce()                          # merge adjacent same-op runs
cig.truncate_to_query_length(50)        # clip from the right by query bases
cig.truncate_to_target_length(50)
str(cig)                                # back to "5S50M2D5M10S"

Cigar.from_cigartuples(rec.cigartuples) # from pysam tuples

# Per-element — these are @property, no parens
elem = cig.elements[1]                  # CigarElement(length=50, operator=CigarOp.M)
elem.length_on_query                    # 50 if op consumes query else 0  (PROPERTY)
elem.length_on_target                   # 50 if op consumes reference else 0  (PROPERTY)

CigarOp.M.consumes_query                # True
CigarOp.S.consumes_reference            # False
CigarOp.I.is_indel                      # True
CigarOp.S.is_clipping                   # True

CigarOp enum: M, I, D, N, S, H, P, EQ, X. Each carries .consumes_query, .consumes_reference, .is_indel, .is_clipping.

SupplementaryAlignment — SA-tag parser

from fgpyo.sam import SupplementaryAlignment

# From a single comma-string
sa = SupplementaryAlignment.parse("chr1,123,-,100M,60,0")
sa.reference_name, sa.start, sa.is_forward, sa.cigar, sa.mapq, sa.nm

# From the full SA tag value
sas = SupplementaryAlignment.parse_sa_tag(rec.get_tag("SA"))

# Or directly from a record (handles missing tag → empty list)
sas = SupplementaryAlignment.from_read(rec)

Pair / template helpers

from fgpyo.sam import (
    PairOrientation, isize, is_proper_pair, set_mate_info,
    set_mate_info_on_secondary, set_mate_info_on_supplementary,
    sum_of_base_qualities, calculate_edit_info,
)

PairOrientation.from_recs(r1, r2)         # FR | RF | TANDEM | None
isize(r1, r2)                              # template length / TLEN
is_proper_pair(r1, r2,
               max_insert_size=1000,
               orientations={PairOrientation.FR})

# Reset all mate-pair fields and MC/MQ tags (replaces manual flag fiddling)
set_mate_info(r1, r2)
set_mate_info_on_secondary(secondary=secondary_r2, mate_primary=primary_r1)
set_mate_info_on_supplementary(supp=supp_r2, mate_primary=primary_r1)

sum_of_base_qualities(rec, min_quality_score=15)
calculate_edit_info(rec, reference_sequence=ref_seq)   # ReadEditInfo

Soft-clipping — fgpyo.sam.clipping

from fgpyo.sam.clipping import (
    softclip_start_of_alignment_by_query,
    softclip_end_of_alignment_by_query,
    softclip_start_of_alignment_by_ref,
    softclip_end_of_alignment_by_ref,
    ClippingInfo,
)

info: ClippingInfo = softclip_start_of_alignment_by_query(
    rec, bases_to_clip=10, clipped_base_quality=None,
    tags_to_invalidate=["MD", "NM", "UQ"],
)
info.query_bases_clipped
info.ref_bases_clipped

The _by_query variants clip a fixed number of query bases off the start/end of the alignment; _by_ref variants clip until N reference bases are dropped. All four mutate the record in place and invalidate alignment-dependent tags.

Constants

from fgpyo.sam import NO_REF_POS, NO_REF_NAME, NO_QUERY_BASES
# NO_REF_POS = -1   NO_REF_NAME = "*"   NO_QUERY_BASES = "*"

Common mistakes

  • Passing a coord-sorted BAM to Template.iterator. It groups by adjacency only — input must be query-name-grouped.
  • Calling add_pair(secondary=True). Not a valid kwarg; add_pair only emits a primary pair. Use add_single(read_num=..., secondary=True, ...) and reuse the same name= to attach to a template.
  • Passing bases=None expecting random bases. That produces a record with no sequence. Omit the kwarg (use the _UNSET default) for random.
  • Writing CIGAR length math by hand. Cigar.length_on_query() / .length_on_target() already exists.
  • Calling cig.length_on_query without parens. Methods on Cigar, properties on CigarElement — easy to mix up. cig.length_on_query returns a bound method object, not an int.
  • Using pysam.AlignmentFile(str(path), "rb") and switching on .endswith(".cram"). fgpyo.sam.reader(path) already does this.

fgpyo util / sequence / read structure / illumina reference

Smaller-surface helpers across fgpyo.read_structure, fgpyo.collections, fgpyo.sequence, and fgpyo.platform.illumina.

ReadStructure — typed sequencing read layouts

fgpyo.read_structure.ReadStructure parses Illumina-style read structure strings like "75T8B75T" (75 template bases, 8 sample-barcode bases, 75 template bases).

from fgpyo.read_structure import ReadStructure, ReadSegment, SegmentType

rs = ReadStructure.from_string("75T8B8M+T")        # last "+T" = variable-length template
rs.segments                                         # tuple[ReadSegment, ...]
[str(s) for s in rs]                                # ["75T", "8B", "8M", "+T"]

# Slice raw bases / qualities by segment
subreads = rs.extract("A" * 158 + "G" * 4)
subreads_q = rs.extract_with_quals(bases, quals)
# Each subread has .bases (and .quals for the with_quals variant), plus the originating ReadSegment.

# Filter segments by kind
templates    = rs.template_segments()
sample_bcs   = rs.sample_barcode_segments()
mol_bcs      = rs.molecular_barcode_segments()
cell_bcs     = rs.cell_barcode_segments()
skips        = rs.skip_segments()

# Helper to make an existing structure end with a variable-length last segment
rs.with_variable_last_segment()

SegmentType enum:

Code Member Meaning
T Template Template / insert bases
B SampleBarcode Sample barcode
M MolecularBarcode UMI / molecular barcode
C CellBarcode Cell barcode
S Skip Bases to discard

ReadSegment exposes .offset, .length (None for variable +), .kind, .has_fixed_length, and .extract(bases) / .extract_with_quals(bases, quals).

PeekableIterator / is_sorted — fgpyo.collections

from fgpyo.collections import PeekableIterator, is_sorted

p = PeekableIterator(iter([1, 2, 3, 4]))
p.peek()                    # 1   (no consumption)
p.can_peek()                # True
next(p)                     # 1
p.takewhile(lambda x: x < 4)  # [2, 3]   (consumes through the predicate)
p.dropwhile(lambda x: x < 5)  # PeekableIterator (returns self after skipping)

is_sorted([1, 2, 2, 3])     # True
is_sorted([3, 1])           # False

PeekableIterator is the standard "I need to look ahead to decide whether to consume" tool — use it instead of manually buffering one item.

fgpyo.sequence — DNA/RNA helpers

from fgpyo.sequence import reverse_complement, complement, gc_content, hamming, levenshtein

reverse_complement("ACGT")          # "ACGT"
reverse_complement("AAGCN")         # "NGCTT"
complement("A")                      # "T"
gc_content("ACGT")                   # 0.5
hamming("ACGT", "ACCT")              # 1     (raises if lengths differ)
levenshtein("ACGT", "AGGT")          # 1

Use these instead of writing custom complement tables or edit-distance loops.

fgpyo.platform.illumina — read-name UMI helpers

from fgpyo.platform.illumina import (
    extract_umis_from_read_name,
    copy_umi_from_read_name,
    SAM_UMI_DELIMITER,                 # "-"
)

# Extract from a read name
umi = extract_umis_from_read_name(
    "M00111:1:1:1:1:1:1:ACGT+CCGG",
    read_name_delimiter=":",            # split read name fields here
    umi_delimiter="+",                  # split multi-UMIs here (joined back with "-")
    strict=False,                       # True = raise on malformed
)
# umi == "ACGT-CCGG"

# Copy onto a record's RX tag
ok = copy_umi_from_read_name(rec, strict=False, remove_umi=False)

Use these for Illumina UMI extraction instead of writing your own read_name.split(":")[-1] parser — they understand the multi-UMI convention and validate UMI characters.

Common mistakes

  • Calling ReadStructure.extract with the wrong total length. Fixed-length structures must match exactly; if the read has variable length, use a structure ending in +T / +B / etc. or call rs.with_variable_last_segment().
  • Using next(it) and a try/except StopIteration to peek. PeekableIterator(it).peek() is the standard idiom and composes with takewhile / dropwhile.
  • Hand-rolling Counter to compute GC content. fgpyo.sequence.gc_content exists.
  • Splitting Illumina read names manually for UMIs. extract_umis_from_read_name handles the multi-UMI +- translation and validates characters.

fgpyo.vcf reference

VariantBuilder — fluent VCF test fixture

fgpyo.vcf.builder.VariantBuilder. Use this in tests instead of constructing pysam.VariantFile and VariantHeader manually. The builder synthesizes a header (with the supplied sequence dictionary and sample list), tracks VariantRecord instances, and can produce a real on-disk VCF.

from pathlib import Path
from fgpyo.vcf.builder import VariantBuilder

builder = VariantBuilder(sample_ids=["S1", "S2"])

# Defaults: contig="chr1", pos=1, ref="A", alts=["T"], qual=60, filter=["PASS"]
v1 = builder.add()
v2 = builder.add(
    contig="chr2", pos=1001, ref="C", alts=["T", "G"], qual=40,
    samples={
        "S1": {"GT": "0|1"},
        "S2": {"GT": "0|0"},
    },
    info={"AF": 0.5},
    id="rs123",
)

# Materialize
unsorted = builder.to_unsorted_list()             # list[VariantRecord]
sorted_records = builder.to_sorted_list()         # contig, then position
vcf_path = builder.to_path()                       # temp .vcf.gz on disk

Constructor

VariantBuilder(
    sample_ids: list[str] | None = None,           # default: a single "sample" sample
    sd: dict[str, dict] | None = None,             # default: HG19 chr1..chrM
)

builder.sample_ids, builder.sd, builder.header (pysam VariantHeader), and builder.records are accessible if you need lower-level control.

add(...) — keyword-only

builder.add(
    contig: str = "chr1",
    pos: int = 1,
    ref: str = "A",
    alts: list[str] = ["T"],
    qual: int = 60,
    filter: list[str] = ["PASS"],
    samples: dict[str, dict[str, Any]] | None = None,   # per-sample FORMAT fields
    id: str | None = None,
    info: dict[str, Any] | None = None,
) -> VariantRecord

Returns the underlying pysam.VariantRecord so you can mutate it further if needed.

to_path / list views

builder.to_unsorted_list() -> list[VariantRecord]
builder.to_sorted_list() -> list[VariantRecord]
builder.to_path(path: Path | None = None) -> Path        # bgzipped + indexed

Common mistakes

  • Forgetting to provide sample IDs. Default is a single sample called "sample"; for genotype tests, pass sample_ids= explicitly.
  • Passing samples={"S1": "0|1"} (string instead of FORMAT dict). Genotypes go in a per-sample sub-dict: samples={"S1": {"GT": "0|1"}}.
  • Using a custom contig without registering it. Contigs come from sd=; if you need a non-default contig, supply your own sd.
  • Hand-building a VariantHeader and VariantFile to test variant code. VariantBuilder produces a valid bgzipped+indexed VCF in two lines.
name bioinformatics-python
description Use when writing Python that handles SAM/BAM/CRAM, FASTA/FASTQ, VCF, or BED files; building test fixtures with bioinformatics records; parsing CIGARs or sequencing read structures; doing genomic interval overlap queries; or modeling typed delimited (TSV/CSV) records. Surfaces fgpyo, pybedlite, and fgmetric features so common primitives (SamBuilder, Template, Cigar, Metric, MetricWriter, BedSource, OverlapDetector) are used instead of hand-rolled equivalents. fgmetric.Metric is preferred over fgpyo.util.Metric.

Using fgpyo, pybedlite, and fgmetric

Three Fulcrum Genomics libraries (fgpyo, pybedlite, fgmetric — all on PyPI) cover most Python bioinformatics primitives. Default to them. Do not hand-roll the equivalents.

If they aren't installed in the project, add them. If a project bans them, do not push.

Metric preference: For typed delimited records, use fgmetric.Metric (Pydantic-backed), not fgpyo.util.Metric (attrs-backed). The fgpyo version is the legacy fallback when fgmetric is unavailable. See references/fgmetric.md.

Built against (verify before trusting specifics)

This skill — including all reference files — was authored against:

  • fgpyo 1.5.1
  • pybedlite 1.1.0
  • fgmetric 0.3.0 (requires Python ≥ 3.12; Pydantic ≥ 2.11.4)
  • Authored 2026-05-08; fgmetric guidance added 2026-05-12; verified against fgmetric 0.3.0 (fgpyo/pybedlite unchanged) 2026-06-05.

Before relying on this skill in a project, check the project's installed versions of these packages (uv pip show fgpyo pybedlite fgmetric, or read pyproject.toml / uv.lock). If any is newer than the version above:

  • Briefly tell the user the skill is older than their pinned version and that you're going to verify any non-trivial API surface you reach for.
  • For each specific class / method / kwarg you use from this skill on the new version, sanity-check it against the project's installed source (e.g. python -c "import fgpyo.sam; help(fgpyo.sam.SamBuilder.add_pair)" or read the upstream changelog) before writing code. Don't assume the table below still matches.
  • If you hit an unexpected AttributeError, TypeError, or ValueError on something this skill recommended, that's the version drift biting — verify against the installed source rather than retrying.

If the project's versions are equal or older, the skill is authoritative.

Reach for, do not re-implement

Symptom in your code Use instead Reference
Hand-built pysam.AlignedSegment ceremony for tests SamBuilder().add_pair(chrom=..., start1=..., start2=...) for primary pairs; .add_single(name=<same>, read_num=1, supplementary=True, ...) to attach supps/secondaries; .to_unsorted_list() to keep records in memory or .to_path(path=tmp_path / "in.bam") for an on-disk BAM sam.md
Buffered loop grouping records by query name Template.iterator(reader) (yields Template with .r1, .r2, .r1_supplementals, etc.) sam.md
Computing CIGAR query/target lengths by hand Cigar.from_cigarstring(s).length_on_query() / .length_on_target() (methods, not properties) sam.md
pysam.AlignmentFile(str(p), "rb") switched on extension fgpyo.sam.reader(path) sam.md
Parsing SA tag fields by hand SupplementaryAlignment.from_read(rec) sam.md
Manually setting mate-info flags / MC / MQ tags set_mate_info(r1, r2) or Template.set_mate_info() sam.md
Custom CSV/TSV writer for tabular metric records Subclass fgmetric.Metric (Pydantic — class M(Metric): name: str; n: int); write with with MetricWriter(M, path) as w: w.write(row) per record (or w.writeall(iterable)); read via M.read(path) (streaming iterator) fgmetric.md
gzip.open(...) / open(...) switched on suffix fgpyo.io.to_reader(path) / to_writer(path) (auto-detects .gz/.bgz) io.md
Raising your own "file not found / not readable" errors assert_path_is_readable(p) / assert_path_is_writable(p) io.md
Custom BED parser splitting tabs and counting fields pybedlite.BedSource(path) (yields BedRecord); pybedlite.writer(path) to write pybedlite.md
Nested loop overlapping every region against every record OverlapDetector.from_bed(p).get_overlaps(query) (interval-tree) pybedlite.md
Hand-built pysam.VariantFile for tests VariantBuilder(sample_ids=...).add(...).to_path() vcf.md
Hand-built pysam.FastaFile / random sequences for tests FastaBuilder().add("chr1").add("ACGT", 100); .to_file(p) fastx.md
zip() over two pysam.FastxFiles with FastxZipped(r1, r2) as z: ... fastx.md
Loading a sequence dictionary out of a .dict / SAM header SequenceDictionary.from_sam(header) / .from_fasta(p) fastx.md
Parser for 75T8B75T-style read layout ReadStructure.from_string(s).extract(bases) util.md
next() + buffer to peek at iterator front PeekableIterator(it).peek() / .takewhile(pred) util.md
Handwritten reverse complement / GC / edit distance fgpyo.sequence.{reverse_complement, gc_content, hamming, levenshtein} util.md
Splitting Illumina read names to extract a UMI extract_umis_from_read_name(name) / copy_umi_from_read_name(rec) util.md

For anything in these domains not in the table — load the matching reference file before writing code and check whether fgpyo/pybedlite already covers it.

Reference files (load on demand)

  • references/sam.mdfgpyo.sam: SamBuilder, Template / TemplateIterator, Cigar ops, clipping, reader/writer, SupplementaryAlignment, PairOrientation, set_mate_info, isize, is_proper_pair, calculate_edit_info.
  • references/fgmetric.mdfgmetric.Metric (Pydantic-backed typed delimited records) and MetricWriter (streaming context-manager writer). Covers list[T] auto-split, Counter[StrEnum] pivot fields, headerless fieldnames= reads, custom delimiters, and the gotchas that differ from fgpyo.util.Metric (no class-level write, no gzip auto-detect, no threads, no header()). Preferred over fgpyo.util.Metric.
  • references/io.mdfgpyo.io: auto-gzip readers/writers, path assertions, line helpers. (Tabular typed records belong in fgmetric.md.)
  • references/fastx.mdFastaBuilder, SequenceDictionary, FastxZipped.
  • references/vcf.mdVariantBuilder.
  • references/util.mdReadStructure / SegmentType, PeekableIterator / is_sorted, fgpyo.sequence, Illumina UMI helpers.
  • references/pybedlite.mdBedSource, BedWriter, BedRecord, BedStrand, Interval, OverlapDetector (+ Span / StrandedSpan protocols).

SamBuilder gotchas (inline — no need to load sam.md)

  • No add_supplementary / add_secondary method. Attach a supplementary or secondary to a pair via add_single(name=<same name as the pair>, read_num=1 or 2, supplementary=True, ...).
  • add_pair does NOT accept supplementary / secondary kwargs. Its only kwargs are name, bases1/2, quals1/2, chrom / chrom1 / chrom2, start1/2, cigar1/2, mapq1/2, strand1/2, attrs. Use add_single for any non-primary record.
  • sort_order takes the SamOrder enum, not a string. SamBuilder(sort_order=SamOrder.QueryName) (import from fgpyo.sam import SamOrder). Default is SamOrder.Coordinate. Passing a string raises ValueError.
  • In-memory vs on-disk: if the function under test takes records, just use the tuples returned from add_pair / add_single (or builder.to_unsorted_list()). Only call .to_path(path=tmp_path / "in.bam") when the function under test actually opens a BAM by path.
  • add_pair(bases1=None) produces a record with no sequence, not a random sequence. Omit the kwarg (let it default to the _UNSET sentinel) for random bases.

fgmetric gotchas (inline — no need to load fgmetric.md)

  • fgmetric.Metric has NO class-level write method. MyMetric.write(path, *rows) is a fgpyo.util.Metric pattern; on fgmetric you must use with MetricWriter(MyMetric, path) as w: w.write(row) (or w.writeall(iterable)).
  • No @attr.s / @dataclass on the subclass. Metric is a Pydantic BaseModel; declare fields with plain type annotations (and optionally Field(...)). Adding attrs decorators breaks Pydantic's field discovery.
  • No gzip auto-detect. Passing Path("out.tsv.gz") to MetricWriter writes uncompressed text into a .gz filename. Write uncompressed and gzip after, or fall back to fgpyo.util.Metric (which does auto-detect) when gzip output is a hard requirement.
  • Metric already inherits CounterPivotTable. Subclass as class M(Metric): — do not write class M(CounterPivotTable, Metric):.
  • MetricWriter must be used as a context manager. The header is written on __enter__; bare construction without with (or without close()) can leave the file open or unflushed.
  • Headerless TSV/CSV: pass fieldnames=[...] to Metric.read(path, fieldnames=...). Order must match the file's column order.

Cross-cutting conventions

  • Frozen attrs is the project style: @attr.s(frozen=True, kw_only=True, auto_attribs=True, slots=True). Match this in code consuming these libraries.
  • 0-based half-open coordinates everywhere in pybedlite. Interval.from_ucsc(...) is the only converter from UCSC chr:start-end 1-based-closed strings.
  • Streaming first. BedSource, Template.iterator, Metric.read, fgpyo.io.to_reader, FastxZipped are line/record streaming. Don't materialize whole files unless required.
  • Builders return real on-disk paths. SamBuilder.to_path(), VariantBuilder.to_path(), FastaBuilder.to_file(path) produce files for tests; to_unsorted_list() / to_sorted_list() keep records in memory.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment