# Biopython Databases

> Fetch sequences from NCBI Entrez/UniProt/PDB/Ensembl REST APIs via BioPython SeqIO; read/write FASTA/GenBank/FASTQ. Use when downloading by accession/gene name, batch-fetching records, or converting sequence file formats.

- Skill: `pavel-kravchenko/biopython-databases` (Agent Skill)
- Install (CLI): `npx skillmds@latest add pavel-kravchenko/biopython-databases`
- Raw SKILL.md: https://api.skillmd.com/api/skills/pavel-kravchenko/biopython-databases/raw
- Safety review: pending
- Works with: Claude Code, Claude.ai, OpenAI Codex
- Category: Integrations & APIs
- Author: pavel-kravchenko (https://skillmd.com/u/pavel-kravchenko)
- Updated: 2026-09-17
- Page: https://skillmd.com/skills/pavel-kravchenko/biopython-databases

---


# BioPython — Biological Databases

## When to Use
- Fetching sequences/annotations from NCBI (Entrez), UniProt, PDB, or Ensembl by accession or gene name
- Reading/writing FASTA, GenBank, EMBL, or FASTQ files with `Bio.SeqIO`
- Manipulating sequences: reverse complement, transcription, translation, GC%, ORF finding
- Batch-downloading many records programmatically instead of one-by-one browser downloads
- Converting between sequence file formats (e.g. FASTQ → FASTA, GenBank → FASTA)

## Version Compatibility
- Biopython ≥ 1.81, Python ≥ 3.9
- NCBI E-utilities (Entrez) REST API — no version string, but requires `Entrez.email`
- UniProt REST API at `rest.uniprot.org` (current since 2022; the old `www.uniprot.org/uniprot/` endpoint is deprecated)
- RCSB PDB Data API (`data.rcsb.org`) and Search API v2 (`search.rcsb.org/rcsbsearch/v2`)
- Ensembl REST API at `rest.ensembl.org` (release-independent JSON endpoints)

## Prerequisites
- `pip install biopython`
- stdlib `urllib.request` / `json` for the UniProt, PDB, and Ensembl REST calls (no extra dependencies needed)
- A registered `Entrez.email`; a free NCBI API key raises the rate limit from 3 to 10 requests/sec
- Familiarity with FASTA/GenBank/FASTQ file formats

## Quick Reference

### NCBI Accession Prefixes
| Prefix | Type | DB |
|--------|------|----|
| NM_ | curated mRNA | RefSeq |
| NR_ | curated ncRNA | RefSeq |
| NP_ | curated protein | RefSeq |
| NC_ | complete chromosome | RefSeq |
| XM_/XP_ | predicted mRNA/protein | RefSeq |
| No prefix (e.g. AY123456) | submitted sequence | GenBank |
| SRP/SRS/SRX/SRR | study/sample/experiment/run | SRA |
| GSE/GSM/GPL | series/sample/platform | GEO |

### SeqIO Format Strings
| Format | Read | Write | Notes |
|--------|------|-------|-------|
| `"fasta"` | yes | yes | no quality or annotation |
| `"genbank"` / `"gb"` | yes | yes | full annotations + features |
| `"fastq"` | yes | yes | per-base quality scores |
| `"embl"` | yes | yes | EBI equivalent of GenBank |

## Key Patterns

### Setup
```python
from Bio.Seq import Seq, MutableSeq
from Bio.SeqRecord import SeqRecord
from Bio.SeqFeature import SeqFeature, FeatureLocation
from Bio import SeqIO, Entrez, AlignIO, Align
from Bio.SeqUtils import gc_fraction, molecular_weight

Entrez.email = "your.email@example.com"  # required before any E-utility call
```

### Seq object
```python
dna = Seq("ATGCGATCGATCGTAA")
dna.complement()             # 3'->5' complement
dna.reverse_complement()     # 5'->3' complement of the other strand
dna.transcribe()             # DNA -> RNA (T -> U)
dna.translate()               # returns protein, * = stop codon
dna.translate(to_stop=True)  # stops before stop codon (most common)
dna.translate(table=2)       # mitochondrial genetic code

gc_fraction(dna) * 100       # GC percentage
molecular_weight(dna)        # DNA MW in Da
molecular_weight(protein, seq_type="protein")

m = MutableSeq("ATGCGATCG")
m[3] = "T"  # point mutation
```

### SeqRecord object
```python
record = SeqRecord(
    Seq("ATGCGATCG"),
    id="GENE_001",
    name="GENE",
    description="Description string"
)
record.annotations["organism"] = "Homo sapiens"
record.annotations["molecule_type"] = "DNA"  # required for GenBank output

feat = SeqFeature(FeatureLocation(0, 27), type="CDS",
                  qualifiers={"gene": ["BRCA1"], "product": ["BRCA1 protein"]})
record.features.append(feat)

record.letter_annotations["phred_quality"] = [30, 28, 35]
cds_seq = feat.location.extract(record.seq)
```

**Goal:** fetch a GenBank record by accession/gene name and translate its CDS.

**Approach:** `Entrez.esearch` to resolve a gene/query to accessions, `Entrez.efetch` with `rettype="gb"` to pull the full annotated record, then walk `record.features` for the `CDS` feature and use `feature.location.extract()` + `.translate()`.

```python
from Bio import Entrez, SeqIO

Entrez.email = "your.email@example.com"

def search_and_fetch_gene(term, db="nucleotide", retmax=5):
    """Resolve a free-text query to accessions, then fetch the top GenBank record."""
    handle = Entrez.esearch(db=db, term=term, retmax=retmax)
    results = Entrez.read(handle)
    handle.close()
    ids = results["IdList"]
    if not ids:
        return None
    handle = Entrez.efetch(db=db, id=ids[0], rettype="gb", retmode="text")
    record = SeqIO.read(handle, "genbank")
    handle.close()
    return record

record = search_and_fetch_gene(
    "insulin[Gene] AND Homo sapiens[Organism] AND mRNA[Filter] AND RefSeq[Filter]"
)
for feature in record.features:
    if feature.type == "CDS":
        protein = feature.location.extract(record.seq).translate(to_stop=True)
        print(feature.qualifiers.get("product", ["?"])[0], protein)
```

### Batch fetch by accession list
```python
accessions = ["NM_000518.5", "NM_000207.3", "NM_000546.6"]  # HBB, INS, TP53
handle = Entrez.efetch(db="nucleotide", id=",".join(accessions), rettype="fasta", retmode="text")
records = list(SeqIO.parse(handle, "fasta"))
handle.close()
SeqIO.write(records, "batch.fasta", "fasta")
```

### elink: cross-database navigation
```python
handle = Entrez.elink(dbfrom="nucleotide", db="protein", id="NM_000207.3")
link_results = Entrez.read(handle)
handle.close()
protein_ids = [link["Id"] for linkset in link_results
               for linkdb in linkset["LinkSetDb"]
               for link in linkdb["Link"]]

handle = Entrez.elink(dbfrom="gene", db="pubmed", id="7157")  # TP53
```

**Goal:** query UniProt, PDB, and Ensembl over their REST APIs.

**Approach:** these databases don't need Biopython — plain `urllib.request` + `json` against their public REST endpoints is simpler and has no extra dependency.

```python
import json
import urllib.parse
import urllib.request


def fetch_uniprot(accession, fmt="json"):
    """Fetch a single UniProtKB entry by accession, e.g. fetch_uniprot('P01308') for human insulin."""
    url = f"https://rest.uniprot.org/uniprotkb/{accession}.{fmt}"
    with urllib.request.urlopen(url) as r:
        return json.loads(r.read()) if fmt == "json" else r.read().decode()


def search_uniprot(query, limit=5):
    """Full-text/field search against UniProtKB, e.g. 'insulin AND reviewed:true AND organism_id:9606'."""
    q = urllib.parse.quote(query)
    url = f"https://rest.uniprot.org/uniprotkb/search?query={q}&size={limit}&format=json"
    with urllib.request.urlopen(url) as r:
        return json.loads(r.read())


def get_pdb_info(pdb_id):
    """Fetch entry-level metadata (title, resolution, method) for a PDB ID, e.g. '1ATP'."""
    url = f"https://data.rcsb.org/rest/v1/core/entry/{pdb_id}"
    with urllib.request.urlopen(url) as r:
        return json.loads(r.read())


def download_pdb(pdb_id, fmt="pdb"):
    """Download a structure file (fmt='pdb' or 'cif') to the current directory."""
    url = f"https://files.rcsb.org/download/{pdb_id}.{fmt}"
    urllib.request.urlretrieve(url, f"{pdb_id}.{fmt}")


def ensembl_lookup(ensembl_id, expand=1):
    """Look up gene/transcript metadata by Ensembl stable ID, e.g. 'ENSG00000141510' (TP53)."""
    url = f"https://rest.ensembl.org/lookup/id/{ensembl_id}?expand={expand}"
    req = urllib.request.Request(url, headers={"Content-Type": "application/json"})
    with urllib.request.urlopen(req) as r:
        return json.loads(r.read())


def ensembl_sequence(ensembl_id, seq_type="cds"):
    """Fetch sequence (seq_type='genomic'|'cds'|'cdna'|'protein') for an Ensembl stable ID."""
    url = f"https://rest.ensembl.org/sequence/id/{ensembl_id}?type={seq_type}"
    req = urllib.request.Request(url, headers={"Content-Type": "application/json"})
    with urllib.request.urlopen(req) as r:
        return json.loads(r.read())
```

### SeqIO: reading, writing, converting
```python
records = list(SeqIO.parse("file.fasta", "fasta"))   # iterator -> wrap in list() to reuse
record = SeqIO.read("single.gb", "genbank")           # exactly one record

SeqIO.write(records, "out.fasta", "fasta")

d = SeqIO.to_dict(SeqIO.parse("file.fasta", "fasta"))
d["INS_HUMAN"].seq

count = SeqIO.convert("reads.fastq", "fastq", "reads.fasta", "fasta")  # quality info is lost
```

### Pairwise alignment
```python
aligner = Align.PairwiseAligner()
aligner.mode = "global"   # Needleman-Wunsch; "local" = Smith-Waterman
alignments = aligner.align(seq1, seq2)
print(alignments[0].score)
```

### ORF finder (all 6 frames)
```python
def find_orfs(seq, min_length=100):
    """Return (strand, frame, start, end, orf_seq) for every ORF >= min_length nt on both strands."""
    orfs = []
    for strand, s in [('+', seq), ('-', seq.reverse_complement())]:
        for frame in range(3):
            for i in range(frame, len(s) - 2, 3):
                if s[i:i+3] == 'ATG':
                    for j in range(i, len(s) - 2, 3):
                        if str(s[j:j+3]) in ('TAA', 'TAG', 'TGA'):
                            if j - i >= min_length:
                                orfs.append((strand, frame, i, j + 3, s[i:j + 3]))
                            break
    return orfs
```

## Pitfalls
- **Always set `Entrez.email`** before any E-utility call; NCBI blocks unidentified clients.
- `SeqIO.parse()` returns an **iterator** — iterating twice yields nothing. Wrap in `list()` to reuse.
- `SeqIO.read()` raises if the file has 0 or >1 records — use `parse()` for multi-record files.
- `record.annotations["molecule_type"]` is **required** for GenBank output; omitting it raises `ValueError`.
- Use `translate(to_stop=True)` for clean protein sequences (omits the trailing `*`).
- `efetch` with `rettype="gb"` + SeqIO: use `"genbank"` (not `"gb"`) as the SeqIO format string.
- Rate-limit NCBI: max 3 requests/sec without an API key, 10/sec with one. Add `time.sleep(0.4)` in loops.
- UniProt REST base URL changed in 2022: use `rest.uniprot.org`, not the deprecated `www.uniprot.org/uniprot/`.

## See Also
- `bio-database-access-entrez-search` — building complex Entrez search queries
- `bio-database-access-entrez-fetch` — fetching records at scale with retries/rate limiting
- `pdb-database` — deeper PDB structure search and download workflows
- `ensembl-database` — Ensembl-specific comparative genomics and homology queries