smithery/gptomics

bio-read-sequences

Read biological sequence files (FASTA, FASTQ, GenBank, EMBL, ABI, SFF) with Biopython Bio.SeqIO, choosing between streaming, in-memory, and on-disk-indexed access.

Installation

$ npx skills add smithery/gptomics --skill bio-read-sequences

Summary

  • Read biological sequence files (FASTA, FASTQ, GenBank, EMBL, ABI, SFF) with Biopython Bio.SeqIO, choosing between streaming, in-memory, and on-disk-indexed access.
  • Use when parsing sequence files, iterating multi-record files, randomly accessing records by ID in large files, or maximizing parse throughput.

Similar popular skills

Related neighbors and high-traction skills in the same topics — useful to compare before installing.

Also in this package

Other skills from smithery/gptomics · top by installs.

npx skills add smithery/gptomics

Browse all from smithery/gptomics

More details

Agent compatibility

Declared targets from SKILL.md / docs. Unmarked agents are not listed — the skill may still install via the CLI.

Claude Code Not declared
Cursor Not declared
Codex Not declared
GitHub Copilot Not declared
Windsurf Not declared
Gemini CLI Not declared
Cline Not declared
OpenCode Not declared

Package contents

Files included with this skill beyond the listing page.

  • skill md SKILL.md 14,090 B
  • docs SUMMARY.md 247 B

History

  1. First recorded snapshot · 0 installs

SKILL.md

Version Compatibility

Reference examples tested with: BioPython 1.83+

Before using code patterns, verify installed versions match. If versions differ:

  • Python: pip show biopython then help(module.function) to check signatures

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

Read Sequences

Read biological sequence data from files using Biopython's Bio.SeqIO module.

"Read sequences from a file" -> Parse a file into SeqRecord objects exposing id, sequence, and annotations.

  • Python: SeqIO.parse() / SeqIO.read() (BioPython)
  • R: readDNAStringSet() / readAAStringSet() (Biostrings)

The Governing Principle

Stream by default. SeqIO.parse() yields one record at a time and never holds the whole file in RAM, so it scales to any size. Reach for an in-memory or indexed structure only when the access pattern demands it: load all records (todict) only for small files needing random access; build an index (index / indexdb) for random access into large files. Never list() a huge file or to_dict() it - that defeats streaming and can exhaust memory.

Which Function to Use

Method Returns Memory model Random access Persists Multi-file
parse(handle, format) generator of SeqRecord one record at a time no no no
read(handle, format) one SeqRecord one record n/a no no
to_dict(records) real dict ALL records in RAM yes no feed combined iterators
index(filename, format) dict-like (read-only) byte offsets only, re-parses on access yes no no
indexdb(idxfile, files, format) dict-like (read-only) on-disk SQLite index yes yes yes

Decision rule: parse for streaming; read for a known single-record file; todict when the file is small and random access by ID is needed; index for random access into one large file; indexdb for files larger than RAM, many files indexed together, or an index reused across runs.

Behavioral traps these methods hide:

  • parse() is a one-pass generator. It is NOT subscriptable (parse(...)[3] raises TypeError), and it EXHAUSTS SILENTLY: a second for loop over the same generator object yields nothing with no error. Re-call parse() for each pass, or list() it once if the file is small.
  • read() fails LOUDLY: zero records raise ValueError: No records found in handle; more than one raises ValueError: More than one record found in handle. Use it as an assertion that the file holds exactly one sequence.
  • todict(), index(), and indexdb() all raise ValueError on a DUPLICATE id (Duplicate key '...'). Supply a key_function to derive unique keys when ids collide.
  • index() needs a FILENAME, not a handle (it must seek). It stores only byte offsets and re-parses the record from disk on every access, so it returns a fresh object each time and mutations do not persist. It is read-only (setitem raises NotImplementedError).
  • index_db() stores the offset index in an on-disk SQLite file. It PERSISTS across sessions (reopen later with just the index filename), and scales beyond RAM and across multiple files (pass a list of filenames). This is the right answer for data larger than memory.

The alphabet= argument still appears in some signatures for back-compatibility but is a no-op since BioPython 1.78; leave it None.

Required Import

from Bio import SeqIO

Reading Records

SeqIO.parse() - Stream Multiple Records

Returns a one-pass iterator of SeqRecord objects. Always pass the format explicitly as the second argument.

for record in SeqIO.parse('sequences.fasta', 'fasta'):
    print(record.id, len(record.seq))

SeqIO.read() - Exactly One Record

Use when the file must contain a single sequence; raises on zero or multiple records.

record = SeqIO.read('single.fasta', 'fasta')

Random Access

SeqIO.to_dict() - Small Files

Loads every record into a dictionary keyed by id. Fast random access, but holds all records in RAM.

records = SeqIO.to_dict(SeqIO.parse('sequences.fasta', 'fasta'))
seq = records['sequence_id'].seq

SeqIO.index() - One Large File

Goal: Random access by id into a large file without loading every record into memory.

Approach: Build an in-memory map of byte offsets keyed by id; each lookup re-parses one record from disk.

Reference (BioPython 1.83+):

records = SeqIO.index('large.fasta', 'fasta')
seq = records['sequence_id'].seq
records.close()

A keyfunction maps the id STRING to a custom key (note: todict's key_function receives the whole record instead):

def get_accession(identifier):
    return identifier.split('.')[0]  # drop the version suffix

records = SeqIO.index('sequences.fasta', 'fasta', key_function=get_accession)

SeqIO.index_db() - Huge / Multiple Files

Goal: Random access into data larger than RAM, or across many files, with the index reusable across runs.

Approach: Persist the offset index in an on-disk SQLite database; reopen it later without re-parsing.

Reference (BioPython 1.83+):

# First call parses the file(s) and builds the SQLite index
records = SeqIO.index_db('index.sqlite', 'large.fasta', 'fasta')
seq = records['sequence_id'].seq
records.close()

# Later sessions reopen instantly with just the index filename
records = SeqIO.index_db('index.sqlite')

# Index multiple files as one database
records = SeqIO.index_db('combined.sqlite', ['file1.fasta', 'file2.fasta'], 'fasta')

High-Performance Parsing

For maximum throughput on large files, low-level parsers (SimpleFastaParser, FastqGeneralIterator) yield raw tuples and skip SeqRecord construction, so they run substantially faster than SeqIO.parse.

SimpleFastaParser

Goal: Parse large FASTA files at maximum speed without SeqRecord overhead.

Approach: Iterate (title, sequence) string tuples directly from the handle.

Reference (BioPython 1.83+):

from Bio.SeqIO.FastaIO import SimpleFastaParser

with open('large.fasta') as handle:
    for title, sequence in SimpleFastaParser(handle):
        if len(sequence) > 1000:
            seq_id = title.split()[0]  # first whitespace token is the id

FastqGeneralIterator

Goal: Parse large FASTQ files at maximum speed.

Approach: Iterate (title, sequence, quality_string) string tuples; decode quality manually if needed.

Reference (BioPython 1.83+):

from Bio.SeqIO.QualityIO import FastqGeneralIterator

with open('reads.fastq') as handle:
    for title, sequence, quality in FastqGeneralIterator(handle):
        avg_qual = sum(ord(c) - 33 for c in quality) / len(quality)  # Phred+33

SeqRecord Attributes

After parsing, each record exposes:

record.id          # first whitespace token of the header (string)
record.name        # same first token (for FASTA, name == id)
record.description # the ENTIRE header after '>', including the id token
record.seq         # sequence data (Seq object; case-preserving)
record.features    # list of SeqFeature objects (GenBank/EMBL)
record.annotations # dict of annotations (organism, molecule_type, ...)
record.letter_annotations  # per-letter dict (e.g. 'phred_quality' list)
record.dbxrefs     # database cross-references

id vs name vs description - the first-space split

A FASTA header >FIRST rest of the line parses to: id = FIRST (the first whitespace token), name = FIRST (same token), description = FIRST rest of the line (the WHOLE header after >, including the id). So >seq1 some desc gives id seq1, name seq1, description seq1 some desc. The id is therefore the leading word of the description, not a separate field - relevant when writing records back out.

Common Formats

Format String Typical Extension Notes
FASTA 'fasta' .fasta, .fa, .fna, .faa Most common
FASTA 2-line 'fasta-2line' .fasta One line per sequence (no wrapping)
FASTQ 'fastq' .fastq, .fq Alias of fastq-sanger (Phred+33)
FASTQ Solexa 'fastq-solexa' .fastq Old Solexa (Solexa+64, scores -5..62)
FASTQ Illumina 'fastq-illumina' .fastq Illumina 1.3-1.7 (Phred+64)
GenBank 'genbank' or 'gb' .gb, .gbk With features/annotations
EMBL 'embl' .embl European format with features
Swiss-Prot 'swiss' .dat UniProt format

FASTQ quality encoding cannot be auto-detected reliably: the same quality line can be valid Phred+33 and Phred+64. Picking the wrong string can silently shift every score by 31. Confirm the encoding before parsing; see fastq-quality for the full encoding decision.

Specialized Formats

Format String Use Case
ABI 'abi' Sanger sequencing trace files (.ab1)
ABI Trimmed 'abi-trim' ABI with low-quality ends trimmed
SFF 'sff' 454/Ion Torrent flowgram data
SFF Trimmed 'sff-trim' SFF with adapter/quality trimming
QUAL 'qual' Quality scores file (pairs with FASTA)
PDB SEQRES 'pdb-seqres' Protein sequences from PDB SEQRES records
PDB ATOM 'pdb-atom' Sequences from ATOM records in PDB
SnapGene 'snapgene' SnapGene .dna files

Reading ABI Trace Files

record = SeqIO.read('sample.ab1', 'abi')
qualities = record.letter_annotations['phred_quality']
record_trimmed = SeqIO.read('sample.ab1', 'abi-trim')  # low-quality ends removed

Reading 454/Ion Torrent SFF

for record in SeqIO.parse('reads.sff', 'sff'):
    print(record.id, len(record.seq))

Reading PDB Sequences

for record in SeqIO.parse('structure.pdb', 'pdb-seqres'):
    print(record.id, record.seq)

Alignment Formats (Read-Only)

Format String Notes
PHYLIP 'phylip' Interleaved; 'phylip-relaxed' allows longer names
Clustal 'clustal' ClustalW output
Stockholm 'stockholm' Rfam/Pfam alignments
NEXUS 'nexus' PAUP/MrBayes format
MAF 'maf' Multiple Alignment Format

Code Patterns

Count Records Without Loading All

count = sum(1 for _ in SeqIO.parse('sequences.fasta', 'fasta'))

Read GenBank with Features

for record in SeqIO.parse('sequence.gb', 'genbank'):
    for feature in record.features:
        if feature.type == 'CDS':
            product = feature.qualifiers.get('product', ['Unknown'])[0]
            cds_seq = feature.extract(record.seq)  # spliced feature sequence

Access FASTQ Quality Scores

for record in SeqIO.parse('reads.fastq', 'fastq'):
    qualities = record.letter_annotations['phred_quality']
    avg_quality = sum(qualities) / len(qualities)

Read From a File Handle

with open('sequences.fasta') as handle:
    for record in SeqIO.parse(handle, 'fasta'):
        print(record.id)

Common Errors

Symptom Cause Fix
Second loop over a parser yields nothing, no error parse() generator exhausted after the first pass Re-call parse() per pass, or list() once for small files
TypeError: 'generator' object is not subscriptable Indexed/sliced a parse() result Wrap in list(), or use to_dict/index for keyed access
ValueError: More than one record found in handle read() on a multi-record file Use parse()
ValueError: No records found in handle read() on an empty/zero-record file Check the file and format string; use parse() if multi-record
ValueError: Duplicate key '...' todict/index/indexdb hit a repeated id Pass a key_function that derives unique keys
Random access by id silently slow / re-reads disk index() re-parses each access; mutations don't persist Expected; cache needed records, or use to_dict for small files
MemoryError / process killed on a huge file list() or to_dict() loaded everything into RAM Stream with parse(); use index_db() for random access
ValueError: unknown format Misspelled format string Use a lowercase string from the format tables
ValueError/AssertionError naming the LOCUS line GenBank parser reads fixed LOCUS columns (molecule type ~44-54, topology ~55-63); ICE/SnapGene/Ensembl/assembler LOCUS lines violate the spec Biologically valid content can still fail the strict column parse; fix the LOCUS columns or re-export from a spec-compliant writer
FASTQ scores all off by ~31 with no error Wrong FASTQ variant string (Phred+33 vs +64 overlap) Confirm encoding; see fastq-quality
AttributeError referencing .alphabet Code assumes pre-1.78 alphabet API Drop alphabet usage; molecule type lives in annotations['molecule_type']

Related Skills

  • write-sequences - Write parsed sequences to new files
  • filter-sequences - Filter sequences by criteria after reading
  • format-conversion - Convert between formats (GenBank->FASTA silently drops annotations)
  • compressed-files - Read gzip/bzip2/BGZF compressed files; only BGZF supports indexed random access
  • fastq-quality - FASTQ encoding (Phred vs Solexa) and offset selection
  • sequence-manipulation/seq-objects - Work with parsed SeqRecord and Seq objects
  • database-access/entrez-fetch - Fetch sequences from NCBI instead of local files
  • alignment-files/sam-bam-basics - For SAM/BAM/CRAM alignment files, use samtools/pysam