references/citation.bib
@ARTICLE{Harrison2024,
ABSTRACT = {Ensembl (https://www.ensembl.org) is a freely available genomic resource that has produced high-quality annotations, tools, and services for vertebrates and model organisms for more than two decades. In recent years, there has been a dramatic shift in the genomic landscape, with a large increase in the number and phylogenetic breadth of high-quality reference genomes, alongside major advances in the pan-genome representations of higher species. In order to support these efforts and accelerate downstream research, Ensembl continues to focus on scaling for the rapid annotation of new genome assemblies, developing new methods for comparative analysis, and expanding the depth and quality of our genome annotations. This year we have continued our expansion to support global biodiversity research, doubling the number of annotated genomes we support on our Rapid Release site to over 1700, driven by our close collaboration with biodiversity projects such as Darwin Tree of Life. We have also strengthened support for key agricultural species, including the first regulatory builds for farmed animals, and have updated key tools and resources that support the global scientific community, notably the Ensembl Variant Effect Predictor. Ensembl data, software, and tools are freely available.},
AUTHOR = {Harrison, Peter W and Amode, M Ridwan and Austine-Orimoloye, Olanrewaju and Azov, Andrey G and Barba, Matthieu and Barnes, If and Becker, Arne and Bennett, Ruth and Berry, Andrew and Bhai, Jyothish and Bhurji, Simarpreet Kaur and Boddu, Sanjay and Branco Lins, Paulo R and Brooks, Lucy and Ramaraju, Shashank Budhanuru and Campbell, Lahcen I and Martinez, Manuel Carbajo and Charkhchi, Mehrnaz and Chougule, Kapeel and Cockburn, Alexander and Davidson, Claire and De Silva, Nishadi H and Dodiya, Kamalkumar and Donaldson, Sarah and El Houdaigui, Bilal and Naboulsi, Tamara El and Fatima, Reham and Giron, Carlos Garcia and Genez, Thiago and Grigoriadis, Dionysios and Ghattaoraya, Gurpreet S and Martinez, Jose Gonzalez and Gurbich, Tatiana A and Hardy, Matthew and Hollis, Zoe and Hourlier, Thibaut and Hunt, Toby and Kay, Mike and Kaykala, Vinay and Le, Tuan and Lemos, Diana and Lodha, Disha and Marques-Coelho, Diego and Maslen, Gareth and Merino, Gabriela Alejandra and Mirabueno, Louisse Paola and Mushtaq, Aleena and Hossain, Syed Nakib and Ogeh, Denye N and Sakthivel, Manoj Pandian and Parker, Anne and Perry, Malcolm and Piližota, Ivana and Poppleton, Daniel and Prosovetskaia, Irina and Raj, Shriya and Pérez-Silva, José G and Salam, Ahamed Imran Abdul and Saraf, Shradha and Saraiva-Agostinho, Nuno and Sheppard, Dan and Sinha, Swati and Sipos, Botond and Sitnik, Vasily and Stark, William and Steed, Emily and Suner, Marie-Marthe and Surapaneni, Likhitha and Sutinen, Kyösti and Tricomi, Francesca Floriana and Urbina-Gómez, David and Veidenberg, Andres and Walsh, Thomas A and Ware, Doreen and Wass, Elizabeth and Willhoft, Natalie L and Allen, Jamie and Alvarez-Jarreta, Jorge and Chakiachvili, Marc and Flint, Bethany and Giorgetti, Stefano and Haggerty, Leanne and Ilsley, Garth R and Keatley, Jon and Loveland, Jane E and Moore, Benjamin and Mudge, Jonathan M and Naamati, Guy and Tate, John and Trevanion, Stephen J and Winterbottom, Andrea and Frankish, Adam and Hunt, Sarah E and Cunningham, Fiona and Dyer, Sarah and Finn, Robert D and Martin, Fergal J and Yates, Andrew D},
URL = {https://doi.org/10.1093/nar/gkad1049},
DATE = {2024-01},
DOI = {10.1093/nar/gkad1049},
EPRINT = {https://academic.oup.com/nar/article-pdf/52/D1/D891/55040594/gkad1049.pdf},
ISSN = {0305-1048},
JOURNALTITLE = {Nucleic Acids Research},
NUMBER = {D1},
PAGES = {D891--D899},
TITLE = {Ensembl 2024},
VOLUME = {52},
}
references/ensembl_rest_api_reference.md
# Ensembl REST API Reference
This document provides a concise reference for the Ensembl REST API
(`https://rest.ensembl.org`). Use it to build custom queries when the
`ensembl_api.py` script does not cover a specific use case.
## General Conventions
- **Base URL:** `https://rest.ensembl.org` (GRCh38). For GRCh37:
`https://grch37.rest.ensembl.org`
- **Content Negotiation:** Set the `Content-Type` header to control the
response format:
- `application/json` — structured JSON (default for most endpoints)
- `text/plain` — raw sequence string
- `text/x-fasta` — FASTA-formatted sequence
- `text/x-gff3` — GFF3 annotation output
- **Rate Limit:** Max 15 requests/second. On HTTP 429, honour the
`Retry-After` header.
- **Region Format:** `CHR:START..END:STRAND` where STRAND is `1` (forward)
or `-1` (reverse). A hyphen (`START-END`) also works for most endpoints.
- **Query Parameter Separator:** Ensembl uses `;` (semicolon) to separate
query parameters, e.g. `?expand=1;mane=1`. Standard `&` also works.
---
## Lookup Endpoints
### `GET /lookup/id/{id}`
Look up any Ensembl stable ID (gene, transcript, protein) and retrieve
metadata.
- **`expand`** (0/1): Include child objects (Transcript array for genes, Exon
array for transcripts, Translation for coding transcripts)
- **`mane`** (0/1): Include MANE Select/Plus Clinical annotations on transcripts
- **`db_type`** (string): Database (default: `core`). Options: `core`,
`otherfeatures`
- **`format`** (string): `full` (default) or `condensed`
- **`species`** (string): Override species if the ID is ambiguous
**Key response fields (Gene):**
`id`, `display_name` (symbol), `biotype`, `description`,
`seq_region_name` (chromosome), `start`, `end`, `strand`,
`assembly_name`, `Transcript[]` (when expanded).
**Key response fields (Transcript, expanded):**
`id`, `biotype`, `display_name`, `is_canonical` (0 or 1), `length`,
`MANE[]` (array with `type`: `MANE_Select` or `MANE_Plus_Clinical`),
`TSL` (Transcript Support Level object with `value`),
`Exon[]`, `Translation` (with `id`, `start`, `end`, `length`).
### `GET /lookup/symbol/{species}/{symbol}`
Resolve a gene symbol to its Ensembl stable ID.
- **`expand`** (0/1): Include child objects
Returns the same structure as `/lookup/id/`.
### `POST /lookup/id`
Batch lookup: send `{"ids": ["ENSG...", "ENST..."]}` as JSON body.
Returns a dict keyed by ID.
---
## Cross-Reference Endpoints
### `GET /xrefs/id/{id}`
Retrieve external database references for an Ensembl ID.
- **`external_db`** (string): Filter by database name (e.g. `UniProt`, `HGNC`,
`RefSeq_mRNA`, `UCSC`, `EntrezGene`)
- **`all_levels`** (0/1): Include xrefs from parent/child features
**Response:** Array of objects with `primary_id`, `display_id`,
`db_display_name`, `dbname`, `description`, `info_type`.
### `GET /xrefs/symbol/{species}/{symbol}`
Find Ensembl IDs matching an external symbol.
### `GET /xrefs/name/{species}/{name}`
Broader search — looks up any name across all external databases.
---
## Sequence Endpoints
### `GET /sequence/id/{id}`
Fetch sequence for an Ensembl feature by stable ID.
- **`type`** (string): `genomic` (default), `cdna`, `cds`, `protein`
- **`expand_5prime`** (int): Extend N bases upstream
- **`expand_3prime`** (int): Extend N bases downstream
- **`mask`** (string): Masking: `hard` or `soft`
Set `Accept: text/x-fasta` for FASTA output, `text/plain` for raw string.
### `GET /sequence/region/{species}/{region}`
Fetch genomic DNA for a coordinate window.
- **`coord_system_version`** (string): Assembly version (e.g. `GRCh38`)
- **`expand_5prime`** (int): Extend N bases upstream
- **`expand_3prime`** (int): Extend N bases downstream
- **`mask`** (string): `hard` or `soft` repeat masking
- **`mask_feature`** (0/1): Apply feature-level masking
Region format: `CHR:START..END:STRAND` (e.g., `X:1000000..1000100:1`).
### `POST /sequence/region/{species}`
Batch: send `{"regions": ["X:1000..2000", "7:100..200"]}`.
---
## Overlap Endpoints
### `GET /overlap/region/{species}/{region}`
Find features overlapping a genomic region. This is useful for finding
genes at a locus, variants in a window, or regulatory features.
- **`feature`** (string): Feature types to return. Repeat for multiple. Values:
`gene`, `transcript`, `cds`, `exon`, `repeat`, `simple`, `misc`, `variation`,
`somatic_variation`, `structural_variation`, `somatic_structural_variation`,
`constrained`, `regulatory`, `motif`, `chipseq`, `array_probe`
- **`biotype`** (string): Filter by biotype (e.g. `protein_coding`)
- **`variant_set`** (string): Short set name for variant filtering
**Example:** Find all genes and transcripts in a region:
```
/overlap/region/human/7:140424943-140624564?feature=gene;feature=transcript
```
### `GET /overlap/id/{id}`
Features that overlap with an Ensembl feature (gene, transcript, etc.).
Same `feature` parameter options as above.
### `GET /overlap/translation/{id}`
Protein features overlapping a translation. Used for domain annotations.
- **`feature`** (string): `protein_feature`, `residue_overlap`,
`translation_exon`
- **`type`** (string): Filter by source database name (e.g. `Pfam`, `Gene3D`,
`CDD`, `Smart`, `SuperFamily`, `PANTHER`, `Prosite_patterns`, `PRINTS`,
`MobiDBLite`)
**Example:** Get Pfam domain annotations for a protein:
```
/overlap/translation/ENSP00000269305?feature=protein_feature;type=Pfam
```
**Response fields (protein_feature):** `type` (source DB), `id` (domain
accession), `description`, `start` (amino acid start), `end`.
---
## Comparative Genomics Endpoints
### `GET /homology/id/{species}/{id}`
Retrieve homologues (orthologs/paralogs) for a gene.
- **`type`** (string): `orthologues`, `paralogues`, `projections`, `all`
- **`target_species`** (string): Restrict to a specific target species
- **`target_taxon`** (int): Restrict by NCBI taxon ID
- **`sequence`** (string): `none`, `cdna`, `protein` — include aligned sequences
### `GET /homology/symbol/{species}/{symbol}`
Same as above, by gene symbol instead of ID.
### `GET /genetree/id/{id}`
Fetch a full gene tree by its Ensembl Compara tree ID.
### `GET /genetree/member/id/{species}/{id}`
Fetch the gene tree containing a specific gene.
---
## Variation Endpoints
### `GET /variation/{species}/{id}`
Retrieve details for a known variant (by rsID or Ensembl variation ID).
**Response fields:** `name` (rsID), `mappings[]` (with `location`,
`allele_string`, `start`, `end`, `strand`), `ancestral_allele`,
`minor_allele`, `MAF`, `clinical_significance[]`, `source`.
### `GET /variant_recoder/{species}/{id}`
Recode a variant between different formats (HGVS, VCF, SPDI, rsID).
**Response fields:** `spdi[]`, `hgvsg[]`, `hgvsc[]`, `hgvsp[]`,
`vcf_string[]`, `id[]` (rsIDs).
---
## VEP (Variant Effect Predictor) Endpoints
### `GET /vep/{species}/region/{region}/{allele}`
Predict consequences for a genomic variant.
Region format: `CHR:START-END:STRAND` (e.g., `9:21971147-21971147:1`).
Allele is the alternate allele string.
### `GET /vep/{species}/id/{id}`
Predict consequences by rsID.
### `GET /vep/{species}/hgvs/{hgvs_notation}`
Predict consequences by HGVS notation.
**Plugin parameters (append as query params):**
- `AlphaMissense=1` — AlphaMissense pathogenicity predictions
- `Conservation=1` — PhyloP conservation scores
- `DosageSensitivity=1` — Haploinsufficiency / triplosensitivity
- `LoF=loftee` — LOFTEE loss-of-function assessment
- `LOEUF=1` — Loss-of-function observed/expected upper bound fraction
- `NMD=1` — Nonsense-mediated decay prediction
- `UTRAnnotator=1` — 5'/3' UTR variant annotation
- `mutfunc=1` — Functional impact prediction
- `IntAct=1` — Protein interaction impact
- `MaveDB=1` — Multiplex assay scores
- `OpenTargets=1` — Open Targets genetics data
**Key response fields:** `most_severe_consequence`,
`transcript_consequences[]` (with `gene_symbol`, `transcript_id`,
`consequence_terms[]`, `amino_acids`, `sift_prediction`, `sift_score`,
`polyphen_prediction`, `polyphen_score`, `am_class`,
`am_pathogenicity`, `conservation`, `lof`, `loeuf`).
### `POST /vep/{species}/region`
Batch VEP: send `{"variants": ["1 100 . A T . . ."]}` in VCF-like format.
---
## Mapping Endpoints
### `GET /map/{species}/{asm_one}/{region}/{asm_two}`
Convert coordinates between assemblies (e.g., GRCh37 → GRCh38).
**Example:**
```
/map/human/GRCh37/17:43044295-43125370/GRCh38
```
**Response:** `mappings[]` with `original` and `mapped` coordinate blocks.
### `GET /map/cdna/{id}/{region}`
Map cDNA coordinates to genomic coordinates for a transcript.
### `GET /map/cds/{id}/{region}`
Map CDS coordinates to genomic coordinates.
### `GET /map/translation/{id}/{region}`
Map protein (amino acid) positions to genomic coordinates.
---
## Phenotype Endpoints
### `GET /phenotype/gene/{species}/{gene}`
Phenotype annotations for a gene (by symbol or Ensembl ID).
### `GET /phenotype/region/{species}/{region}`
Phenotype-associated variants in a genomic region.
### `GET /phenotype/term/{species}/{term}`
Find variants/genes associated with a phenotype term (ontology ID or
description string, e.g., `coffee consumption`).
---
## Linkage Disequilibrium Endpoints
### `GET /ld/{species}/{id}/{population_name}`
Compute LD for variants in a window around a variant.
- **`window_size`** (int): Window size in kb (default: 500)
- **`r2`** (float): Minimum r² threshold
- **`d_prime`** (float): Minimum D' threshold
### `GET /ld/{species}/pairwise/{id1}/{id2}`
Pairwise LD between two specific variants.
### `GET /ld/{species}/region/{region}/{population_name}`
LD for all variant pairs in a region.
Population names: e.g., `1000GENOMES:phase_3:CEU`, `1000GENOMES:phase_3:YRI`.
---
## Information Endpoints
### `GET /info/assembly/{species}`
Assembly metadata (karyotype, top-level regions, coordinate systems).
### `GET /info/assembly/{species}/{region_name}`
Details for a specific chromosome/region (length, bands).
### `GET /info/species`
List all species available in the Ensembl REST API.
### `GET /info/external_dbs/{species}`
List all external database names available for cross-references.
Useful for finding the correct `external_db` parameter value.
### `GET /info/biotypes/{species}`
List all biotype classifications for a species.
---
## Best Practices for Custom Queries
1. **Always set `Content-Type`** to `application/json` for JSON responses
or `text/plain` / `text/x-fasta` for sequence endpoints.
2. **Use `expand=1`** on lookup endpoints to get child features in a single
call instead of making separate requests for each transcript/exon.
3. **Prefer batch endpoints** (`POST /lookup/id`, `POST /sequence/region`)
when querying multiple IDs or regions — this reduces the number of HTTP
round-trips.
4. **Check `info/external_dbs`** before using the `external_db` filter on
xrefs — the exact database name strings are case-sensitive and not
always obvious (e.g., `UniProt_gn` not `UniProt`).
5. **Region size limits**: The `/overlap/region/` endpoint has a maximum
region size of 5 Mb for most feature types. Split larger regions.
6. **Use `grch37.rest.ensembl.org`** as the base URL when your coordinates
are on the GRCh37 (hg19) assembly. Most endpoints support the same
paths on this older server.
7. **Save responses to temp files** — do not try to read large JSON
responses into context. Use `jq` or Python one-liners to extract
specific fields.
scripts/ensembl_api.py
# Copyright 2026 Google LLC
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
"""A command-line tool to query the Ensembl REST API.
This script provides subcommands for gene/transcript/protein lookup, ID
resolution and cross-referencing, sequence retrieval, and variant effect
prediction (VEP). All rich data is saved to a temporary JSON file; a concise
human-readable summary is printed to stdout.
"""
# /// script
# requires-python = ">=3.10"
# dependencies = [
# "polite-http",
# ]
# ///
import argparse
import json
import os
import sys
import tempfile
from polite_http import http_client
BASE_URL = "https://rest.ensembl.org"
GRCH37_URL = "https://grch37.rest.ensembl.org"
VEP_PLUGINS = (
"?AlphaMissense=1&Conservation=1&DosageSensitivity=1&IntAct=1"
"&MaveDB=1&OpenTargets=1&LoF=loftee&NMD=1&UTRAnnotator=1"
"&mutfunc=1&LOEUF=1"
)
_CLIENT_REGULAR = http_client.HttpClient(BASE_URL, qps=15)
_CLIENT_GRCH37 = http_client.HttpClient(GRCH37_URL, qps=15)
def _get_client(assembly=None):
"""Return the correct client for the requested assembly."""
if assembly and assembly.upper() == "GRCH37":
print("[*] Using GRCh37 assembly.")
return _CLIENT_GRCH37
return _CLIENT_REGULAR
def _get_species(args):
"""Return the species, defaulting to 'human' if not specified."""
species = args.species
if not species:
species = "human"
print("[*] No species specified. Defaulting to 'human'.")
return species
def _save_json(data, prefix, output_path=None):
"""Write *data* as pretty-printed JSON to a temp file and return the path."""
if output_path:
path = output_path
with open(path, "w") as fh:
json.dump(data, fh, indent=2)
else:
fd, path = tempfile.mkstemp(
prefix=f"ensembl_{prefix}_", suffix=".json", text=True
)
with os.fdopen(fd, "w") as fh:
json.dump(data, fh, indent=2)
return path
def _try_fallback(url: str, query: str, client=None) -> list[dict[str, str]]:
"""Tries to resolve gene symbol via cross-references if initial lookup fails.
Args:
url: The URL for the cross-reference lookup.
query: The original gene symbol query.
client: The HttpClient to use. Defaults to _CLIENT_REGULAR.
Returns:
A list of dictionaries, where each dictionary contains gene information
resolved from the cross-reference.
"""
if client is None:
client = _CLIENT_REGULAR
fallback_data = client.fetch_json(url)
if not fallback_data:
return []
matches = []
for item in fallback_data:
if item.get("type") != "gene":
continue
ens_id = item.get("id")
url_lookup = f"/lookup/id/{ens_id}"
try:
gene_data = client.fetch_json(url_lookup)
matches.append(gene_data)
except http_client.HttpError:
print(f"[!] Failed to fetch details for resolved ID {ens_id}")
matches.append({
"id": ens_id,
"biotype": "N/A",
"seq_region_name": "?",
"start": "?",
"end": "?",
"strand": "?",
})
if matches:
print(f"[*] Resolved via synonym '{query}' to {len(matches)} gene(s):")
return matches
def cmd_resolve_gene(args):
"""Resolve a symbol / alias / RefSeq ID to one or more ENSG IDs."""
query = args.query
species = _get_species(args)
client = _get_client(args.assembly)
url = f"/lookup/symbol/{species}/{query}?expand=0"
try:
data = client.fetch_json(url)
matches = data if isinstance(data, list) else [data]
print(f"[*] Symbol '{query}' resolved ({len(matches)} match(es)):")
except http_client.HttpError as e:
if e.status_code == 404 or (
e.status_code == 400
and e.body
and b"No valid lookup found for symbol" in e.body
):
print(f"[*] Symbol '{query}' not found. Trying synonym resolution...")
url_fallback = f"/xrefs/symbol/{species}/{query}"
try:
data = _try_fallback(url_fallback, query, client=client)
except http_client.HttpError as exc:
raise e from exc
if not data:
raise e
matches = data
else:
raise e
for m in matches:
eid = m.get("id", "N/A")
biotype = m.get("biotype", "N/A")
chrom = m.get("seq_region_name", "?")
start = m.get("start", "?")
end = m.get("end", "?")
strand = m.get("strand", "?")
print(
f" ENSG: {eid} | biotype: {biotype} "
f"| location: chr{chrom}:{start}-{end} (strand {strand})"
)
path = _save_json(data, f"resolve_{query}", output_path=args.output)
print(f"[*] Full JSON saved to {path}")
def cmd_map_id(args):
"""Cross-reference an Ensembl ID to an external database."""
eid = args.id
ext_db = args.external_db
client = _get_client(args.assembly)
params = f"?external_db={ext_db}" if ext_db else ""
url = f"/xrefs/id/{eid}{params}"
data = client.fetch_json(url)
if not data:
print(f"[*] No cross-references found for {eid}.")
else:
db_label = f" in {ext_db}" if ext_db else ""
print(f"[*] {len(data)} cross-reference(s) for {eid}{db_label}:")
for entry in data[:10]:
primary = entry.get("primary_id", "N/A")
display = entry.get("display_id", "N/A")
dbname = entry.get("db_display_name", entry.get("dbname", "N/A"))
print(f" {dbname}: {primary} ({display})")
if len(data) > 10:
print(f" … and {len(data) - 10} more (see JSON).")
path = _save_json(data, f"mapid_{eid}", output_path=args.output)
print(f"[*] Full JSON saved to {path}")
def cmd_get_sequence(args):
"""Fetch raw genomic DNA for a coordinate window."""
coords = args.coords
species = _get_species(args)
assembly = args.assembly
# Normalise the region string: accept chr17:100-200, 17:100-200,
# 17:100..200
region = coords.replace(",", "").lower().removeprefix("chr")
region = region.replace("-", "..")
client = _get_client(assembly)
url = f"/sequence/region/{species}/{region}?"
if assembly:
url += f"coord_system_version={assembly}&"
headers = {"Accept": "text/plain"}
seq = client.fetch_text(url, headers=headers)
if args.output:
path = args.output
with open(path, "w") as fh:
fh.write(str(seq))
else:
fd, path = tempfile.mkstemp(prefix="ensembl_seq_", suffix=".txt", text=True)
with os.fdopen(fd, "w") as fh:
fh.write(str(seq))
length = len(str(seq))
print(f"[*] Fetched genomic sequence for {coords} (length: {length} bp).")
print(f"[*] Sequence saved to {path}")
def cmd_gene_summary(args):
"""Return high-level metadata for a gene by its ENSG ID."""
ensg = args.ensg_id
client = _get_client(args.assembly)
url = f"/lookup/id/{ensg}"
data = client.fetch_json(url)
symbol = data.get("display_name", "N/A")
biotype = data.get("biotype", "N/A")
desc = data.get("description", "N/A")
chrom = data.get("seq_region_name", "?")
start = data.get("start", "?")
end = data.get("end", "?")
strand = "+" if data.get("strand", 1) == 1 else "-"
assembly = data.get("assembly_name", "N/A")
print(f"[*] Gene summary for {ensg}:")
print(f" Symbol: {symbol}")
print(f" Biotype: {biotype}")
print(f" Description: {desc}")
print(f" Location: chr{chrom}:{start}-{end} ({strand})")
print(f" Assembly: {assembly}")
path = _save_json(data, f"gene_{ensg}", output_path=args.output)
print(f"[*] Full JSON saved to {path}")
def cmd_transcripts(args):
"""List transcripts for a gene, with optional MANE / canonical filtering."""
ensg = args.ensg_id
client = _get_client(args.assembly)
url = f"/lookup/id/{ensg}?expand=1;mane=1"
data = client.fetch_json(url)
transcripts = data.get("Transcript", [])
if not transcripts:
print(f"[*] No transcripts found for {ensg}.")
path = _save_json(data, f"transcripts_{ensg}", output_path=args.output)
print(f"[*] Full JSON saved to {path}")
return
# Apply filters
filtered = transcripts
if args.only_mane:
filtered = [t for t in transcripts if t.get("MANE")]
elif args.only_canonical:
filtered = [t for t in transcripts if t.get("is_canonical") == 1]
if not filtered and (args.only_mane or args.only_canonical):
label = "MANE Select" if args.only_mane else "Canonical"
print(
f"[*] No {label} transcript found for {ensg}. "
f"Total transcripts: {len(transcripts)}."
)
return
filter_label = ""
if args.only_mane:
filter_label = " (MANE Select only)"
elif args.only_canonical:
filter_label = " (Canonical only)"
print(
f"[*] {len(filtered)} transcript(s) for {ensg}{filter_label} "
f"(from {len(transcripts)} total):"
)
print()
print("| Transcript ID | Biotype | TSL | Length (bp) | Flags |")
print("| --- | --- | --- | --- | --- |")
for t in filtered:
tid = t.get("id", "N/A")
biotype = t.get("biotype", "N/A")
tsl = t.get("TSL", {})
tsl_val = tsl.get("value") if isinstance(tsl, dict) else tsl
if tsl_val is None:
tsl_val = "N/A"
length = t.get("length", "N/A")
flags = []
if t.get("is_canonical") == 1:
flags.append("Canonical")
mane = t.get("MANE")
if mane:
for m in mane:
flags.append(m.get("type", "MANE"))
flag_str = ", ".join(flags) if flags else "-"
print(f"| {tid} | {biotype} | {tsl_val} | {length} | {flag_str} |")
path = _save_json(data, f"transcripts_{ensg}", output_path=args.output)
print(f"\n[*] Full JSON saved to {path}")
def cmd_canonical_tss(args):
"""Retrieve the TSS for the canonical transcript of a gene."""
query = args.gene
species = _get_species(args)
client = _get_client(args.assembly)
# 1. Resolve to ID if symbol
ensg = query
if not query.lower().startswith("ens"):
url = f"/lookup/symbol/{species}/{query}?expand=0"
data = client.fetch_json(url)
if isinstance(data, list) and len(data) > 1:
print(
f"[!] Warning: '{query}' resolved to {len(data)} genes. "
f"Using first match: {data[0].get('id')} "
f"({data[0].get('display_name', 'N/A')}). "
f"Other matches: {', '.join(d.get('id', '?') for d in data[1:])}"
)
ensg = data.get("id") if not isinstance(data, list) else data[0].get("id")
if not ensg:
print(
f"[!] Could not resolve symbol {query} to Ensembl ID.",
file=sys.stderr,
)
sys.exit(1)
# 2. Fetch transcripts
url = f"/lookup/id/{ensg}?expand=1;mane=1"
data = client.fetch_json(url)
transcripts = data.get("Transcript", [])
if not transcripts:
print(f"[!] No transcripts found for {ensg}.", file=sys.stderr)
sys.exit(1)
canonical = None
for t in transcripts:
if t.get("is_canonical") == 1:
canonical = t
break
if not canonical:
print(f"[!] No Canonical transcript found for {ensg}.", file=sys.stderr)
sys.exit(1)
strand = canonical.get("strand", 1)
start = canonical.get("start")
end = canonical.get("end")
tid = canonical.get("id")
if strand == 1:
tss = start
else:
tss = end
chrom = data.get("seq_region_name", "?")
print(f"[*] Gene {ensg} (chr{chrom})")
print(f"[*] Canonical Transcript: {tid}")
print(f"[*] Strand: {'+' if strand == 1 else '-'}")
print(f"[*] TSS Coordinate: {tss}")
path = _save_json(canonical, f"canonical_tss_{ensg}", output_path=args.output)
print(f"[*] Full canonical transcript JSON saved to {path}")
def cmd_transcript_structure(args):
"""Return exon, CDS, and UTR layout for a transcript."""
enst = args.transcript_id
client = _get_client(args.assembly)
url = f"/lookup/id/{enst}?expand=1;mane=1"
data = client.fetch_json(url)
exons = data.get("Exon", [])
strand = data.get("strand", 1)
trans_start = data.get("start")
trans_end = data.get("end")
print(f"[*] Transcript structure for {enst}:")
print(f" Biotype: {data.get('biotype', 'N/A')}")
print(
f" Genomic span: chr{data.get('seq_region_name', '?')}:"
f"{trans_start}-{trans_end} (strand {'+' if strand == 1 else '-'})"
)
print(f" Exons: {len(exons)}")
translation = data.get("Translation")
utr5 = None
utr3 = None
if translation:
cds_start = translation.get("start")
cds_end = translation.get("end")
cds_length_aa = translation.get("length", "N/A")
ensp = translation.get("id", "N/A")
print(f" CDS: {cds_start}-{cds_end} ({cds_length_aa} aa, {ensp})")
# Compute UTRs
utr5 = (
(cds_end + 1, trans_end)
if cds_end and trans_end and cds_end < trans_end
else None
)
utr3 = (
(trans_start, cds_start - 1)
if cds_start and trans_start and cds_start > trans_start
else None
)
if strand == 1:
utr3, utr5 = utr5, utr3
if utr5:
print(f" 5' UTR: {utr5[0]}-{utr5[1]}")
if utr3:
print(f" 3' UTR: {utr3[0]}-{utr3[1]}")
else:
print(" (non-coding – no CDS/UTR)")
if exons:
print()
print("| Exon # | ID | Start | End | Length (bp) |")
print("| --- | --- | --- | --- | --- |")
sorted_exons = sorted(exons, key=lambda e: e.get("start", 0))
for i, ex in enumerate(sorted_exons, 1):
eid = ex.get("id", "N/A")
estart = ex.get("start", "?")
eend = ex.get("end", "?")
elen = (
eend - estart + 1
if isinstance(estart, int) and isinstance(eend, int)
else "?"
)
print(f"| {i} | {eid} | {estart} | {eend} | {elen} |")
# Enrich the saved data with computed UTR info
if translation:
data["_computed_utrs"] = {}
if utr5:
data["_computed_utrs"]["5_prime"] = {"start": utr5[0], "end": utr5[1]}
if utr3:
data["_computed_utrs"]["3_prime"] = {"start": utr3[0], "end": utr3[1]}
path = _save_json(data, f"structure_{enst}", output_path=args.output)
print(f"\n[*] Full JSON saved to {path}")
def cmd_protein_info(args):
"""Fetch ENSP ID and sequence length for a transcript."""
enst = args.transcript_id
client = _get_client(args.assembly)
url = f"/lookup/id/{enst}?expand=1"
data = client.fetch_json(url)
translation = data.get("Translation")
if not translation:
print(f"[*] {enst} has no translation (likely non-coding).")
path = _save_json(data, f"protein_{enst}", output_path=args.output)
print(f"[*] Full JSON saved to {path}")
return
ensp = translation.get("id", "N/A")
length = translation.get("length", "N/A")
print(f"[*] Protein for {enst}:")
print(f" ENSP: {ensp}")
print(f" Length: {length} aa")
path = _save_json(data, f"protein_{enst}", output_path=args.output)
print(f"[*] Full JSON saved to {path}")
def cmd_protein_sequence(args):
"""Fetch the amino acid sequence (FASTA) for a transcript or protein ID."""
target = args.id
client = _get_client(args.assembly)
url = f"/sequence/id/{target}?type=protein"
headers = {"Accept": "text/x-fasta"}
fasta = client.fetch_text(url, headers=headers)
if args.output:
path = args.output
with open(path, "w") as fh:
fh.write(str(fasta))
else:
fd, path = tempfile.mkstemp(
prefix=f"ensembl_protseq_{target}_", suffix=".fasta", text=True
)
with os.fdopen(fd, "w") as fh:
fh.write(str(fasta))
# Count sequence length (exclude header lines)
lines = str(fasta).strip().split("\n")
seq = "".join(l for l in lines if not l.startswith(">"))
print(f"[*] Protein sequence for {target}: {len(seq)} aa")
print(f"[*] FASTA saved to {path}")
def cmd_vep(args):
"""Predict variant consequences using the Ensembl VEP."""
variant = args.variant_str
species = _get_species(args)
client = _get_client(args.assembly)
if variant.startswith("rs"):
url = f"/vep/{species}/id/{variant}{VEP_PLUGINS}"
else:
parts = variant.split(":")
if len(parts) == 4:
chrom, pos, ref, alt = parts
end = int(pos) + len(ref) - 1
region = f"{chrom}:{pos}-{end}:1"
url = f"/vep/{species}/region/{region}/{alt}{VEP_PLUGINS}"
else:
# Fallback: treat as HGVS
url = f"/vep/{species}/hgvs/{variant}{VEP_PLUGINS}"
data = client.fetch_json(url)
if not isinstance(data, list) or not data:
print("[!] No VEP results returned.")
path = _save_json(
data if data else {},
f"vep_{variant.replace(':', '_')}",
output_path=args.output,
)
print(f"[*] JSON saved to {path}")
return
top = data[0]
mcv = top.get("most_severe_consequence", "Unknown")
input_var = top.get("input", variant)
t_conseq = top.get("transcript_consequences", [])
print(f"[*] Variant: {input_var}")
print(f"[*] Most severe consequence: {mcv}")
print(f"[*] Found {len(t_conseq)} transcript consequences.")
# Build the predictions table
open_keys = {
"am_class": "AlphaMissense Class",
"am_pathogenicity": "AlphaMissense Pathogenicity",
"conservation": "Conservation",
"phaplo": "Dosage Sensitivity (Haplo)",
"ptriplo": "Dosage Sensitivity (Triplo)",
"lof": "Loss of Function (Loftee)",
"nmd": "Nonsense-mediated Decay",
"utr_annotator": "UTR Annotator",
"mutfunc": "Mutfunc",
"loeuf": "Loss of Function (LOEUF)",
}
rows = []
for tc in t_conseq:
tid = tc.get("transcript_id", "Unknown")
gene = tc.get("gene_symbol", "Unknown")
terms = tc.get("consequence_terms", [])
if terms:
rows.append((tid, gene, "Consequence", ", ".join(terms)))
aa = tc.get("amino_acids")
if aa:
rows.append((tid, gene, "Amino Acids", aa))
sift = tc.get("sift_prediction")
if sift:
score = tc.get("sift_score")
val = f"{sift} ({score})" if score is not None else sift
rows.append((tid, gene, "SIFT", val))
poly = tc.get("polyphen_prediction")
if poly:
score = tc.get("polyphen_score")
val = f"{poly} ({score})" if score is not None else poly
rows.append((tid, gene, "PolyPhen", val))
for key, label in open_keys.items():
val = tc.get(key)
if val is None and "alphamissense" in tc:
val = tc["alphamissense"].get(key)
if val is not None:
if isinstance(val, (list, dict)):
val = json.dumps(val)
rows.append((tid, gene, label, str(val)))
if rows:
print("\n[*] VEP Predictions Table:")
print("| Transcript ID | Gene | Method/Metric | Value |")
print("| --- | --- | --- | --- |")
for r in rows:
print(f"| {r[0]} | {r[1]} | {r[2]} | {r[3]} |")
else:
print("\n[*] No detailed predictions in transcript consequences.")
safe = variant.replace(":", "_").replace(">", "_")
path = _save_json(data, f"vep_{safe}", output_path=args.output)
print(f"\n[*] Full JSON saved to {path}")
def main():
"""Parse CLI arguments and dispatch to the appropriate subcommand."""
parent_parser = argparse.ArgumentParser(add_help=False)
parent_parser.add_argument(
"--output",
help="Output file path (optional)",
)
parent_parser.add_argument(
"--assembly",
default=None,
help="Assembly (e.g. GRCh38, GRCh37). Default: GRCh38.",
)
parser = argparse.ArgumentParser(
description="Query the Ensembl REST API.",
formatter_class=argparse.RawDescriptionHelpFormatter,
)
sub = parser.add_subparsers(dest="command", required=True)
# ---- resolve-gene ----
p = sub.add_parser(
"resolve-gene",
parents=[parent_parser],
help="Resolve a gene symbol / alias / RefSeq ID to ENSG ID(s).",
)
p.add_argument("query", help="Gene symbol, alias, or RefSeq ID")
p.add_argument(
"--species",
default=None,
help="Species (defaults to 'human' if not specified)",
)
p.set_defaults(func=cmd_resolve_gene)
# ---- map-id ----
p = sub.add_parser(
"map-id",
parents=[parent_parser],
help="Cross-reference an Ensembl ID to external databases.",
)
p.add_argument("id", help="Ensembl ID (ENSG, ENST, ENSP)")
p.add_argument(
"--external-db",
dest="external_db",
default=None,
help="Filter by external DB (e.g., UniProt, HGNC, RefSeq_mRNA)",
)
p.set_defaults(func=cmd_map_id)
# ---- get-sequence ----
p = sub.add_parser(
"get-sequence",
parents=[parent_parser],
help="Fetch raw genomic DNA for a coordinate window.",
)
p.add_argument(
"coords",
help="Genomic region, e.g. 17:7661779-7687550 or chr17:7661779-7687550",
)
p.add_argument(
"--species",
default=None,
help="Species (defaults to 'human' if not specified)",
)
p.set_defaults(func=cmd_get_sequence)
# ---- gene-summary ----
p = sub.add_parser(
"gene-summary",
parents=[parent_parser],
help="Get high-level metadata for a gene (symbol, biotype, location).",
)
p.add_argument("ensg_id", help="Ensembl gene ID (e.g. ENSG00000141510)")
p.set_defaults(func=cmd_gene_summary)
# ---- transcripts ----
p = sub.add_parser(
"transcripts",
parents=[parent_parser],
help="List transcripts for a gene, with optional MANE/canonical filter.",
)
p.add_argument("ensg_id", help="Ensembl gene ID")
grp = p.add_mutually_exclusive_group()
grp.add_argument(
"--only-mane",
action="store_true",
help="Return only the MANE Select transcript (human only).",
)
grp.add_argument(
"--only-canonical",
action="store_true",
help="Return only the Ensembl Canonical transcript.",
)
p.set_defaults(func=cmd_transcripts)
# ---- canonical-tss ----
p = sub.add_parser(
"canonical-tss",
parents=[parent_parser],
help="Get TSS coordinate for the canonical transcript of a gene.",
)
p.add_argument("gene", help="Gene symbol or Ensembl ID")
p.add_argument(
"--species",
default=None,
help="Species (defaults to 'human' if not specified)",
)
p.set_defaults(func=cmd_canonical_tss)
# ---- transcript-structure ----
p = sub.add_parser(
"transcript-structure",
parents=[parent_parser],
help="Get exon, CDS, and UTR layout for a transcript.",
)
p.add_argument(
"transcript_id", help="Ensembl transcript ID (e.g. ENST00000269305)"
)
p.set_defaults(func=cmd_transcript_structure)
# ---- protein-info ----
p = sub.add_parser(
"protein-info",
parents=[parent_parser],
help="Get ENSP ID and sequence length for a transcript.",
)
p.add_argument(
"transcript_id", help="Ensembl transcript ID (e.g. ENST00000269305)"
)
p.set_defaults(func=cmd_protein_info)
# ---- protein-sequence ----
p = sub.add_parser(
"protein-sequence",
parents=[parent_parser],
help="Fetch the amino acid sequence (FASTA) for an ENST or ENSP.",
)
p.add_argument("id", help="Ensembl transcript or protein ID")
p.set_defaults(func=cmd_protein_sequence)
# ---- vep ----
p = sub.add_parser(
"vep",
parents=[parent_parser],
help="Predict variant consequences (VEP) with open-license plugins.",
)
p.add_argument(
"variant_str",
help="Variant as chr:pos:ref:alt (e.g. 9:21971147:T:C) or rsID",
)
p.add_argument(
"--species",
default=None,
help="Species (defaults to 'human' if not specified)",
)
p.set_defaults(func=cmd_vep)
args = parser.parse_args()
args.func(args)
if __name__ == "__main__":
main()
SKILL.md
---
name: ensembl-database
description: >
Query the Ensembl database to resolve gene, transcript, and protein IDs,
fetch genomic or protein sequences, retrieve gene structures (exons), and
get variant consequence and effect predictions (VEP). Use this skill as a
primary ID translator, genomic sequence database and variant effect prediction
tool.
---
# Ensembl Database: ID Mapping and Genomic Features
## Prerequisites
1. **`uv`**: Read the `uv` skill and follow its Setup instructions to ensure
`uv` is installed and on PATH.
2. **User Notification**: If .licenses/ensembl_database_LICENSE.txt does not
already exist in the workspace root directory then (1) prominently notify
the user to check the terms at https://useast.ensembl.org/index.html and
https://github.com/Ensembl/ensembl-rest/wiki, then (2) create the file
recording the notification text and timestamp.
## Overview
The Ensembl database is a resource for genome annotation. This skill allows you
to interact with the Ensembl REST API to resolve ambiguous symbols,
cross-reference IDs (RefSeq, HGNC, UniProt, ENSG), fetch raw sequences, and
retrieve detailed transcript structures.
**Key Concepts:**
- **ENSG (Gene):** Stable identifier for a human gene. Other species will have
different three-letter species codes.
- **ENST (Transcript):** Stable identifier for a transcript (splicing
isoform).
- **ENSP (Protein):** Stable identifier for a translated protein.
- **MANE Select:** The consensus primary transcript agreed upon by Ensembl and
NCBI.
- **Canonical:** Ensembl's representative transcript (used if MANE is not
available or non-human).
## Core Rules
- **Use the Wrapper**: ALWAYS execute the provided helper scripts to query the
database rather than accessing the database directly. The scripts
automatically enforce the required rate limit gracefully.
- **Default Species:** If the species is absent or ambiguous in the prompt,
default to `"human"`. You MUST explicitly flag this default to the user to
ensure they are aware.
- **Primary Transcripts:** When listing transcripts for a gene, only return
the MANE Select transcript (for human) or the Canonical transcript (for
others) unless the user explicitly asks for all alternative isoforms. You
MUST flag to the user when multiple transcripts are available and you are
defaulting to the primary one.
- **Assembly Handling:** The default assembly is GRCh38. For GRCh37 requests,
you MUST use the `--assembly GRCh37` flag. You MUST explicitly flag to the
user when a non-default assembly is being used.
- **Output Location:** The script writes full JSON/FASTA output to temporary
files in `/tmp` by default, or to a user-specified file using the `--output`
flag. It also prints a concise summary to stdout.
- **Notification**: If this skill is used, ensure this is mentioned in the
output.
### Available Commands
**1. Resolve Gene ID** — Resolve a symbol, alias, or RefSeq ID to ENSG ID(s).
Automatically falls back to resolving synonyms if primary symbol is not found.
```bash
uv run scripts/ensembl_api.py resolve-gene TP53 --species human --output tp53.json
uv run scripts/ensembl_api.py resolve-gene PCL2 --output pcl2.json # Falls back to synonym resolution
```
**2. Map ID to External Database** — Cross-reference an Ensembl ID to UniProt,
HGNC, RefSeq, etc.
```bash
uv run scripts/ensembl_api.py map-id ENSG00000141510 --external-db UniProt --output uniprot_map.json
uv run scripts/ensembl_api.py map-id ENST00000269305 --external-db RefSeq_mRNA --output refseq_map.json
```
**3. Get Genomic Sequence** — Fetch raw DNA for a coordinate window. Supports
GRCh37 via `--assembly GRCh37`.
```bash
uv run scripts/ensembl_api.py get-sequence 17:7661779-7687550 --species human --output seq.txt
uv run scripts/ensembl_api.py get-sequence chr9:21971100-21971200 --assembly GRCh37 --output seq_grch37.txt
```
**4. Gene Summary** — High-level metadata: symbol, biotype, description,
chromosomal location.
```bash
uv run scripts/ensembl_api.py gene-summary ENSG00000141510 --output gene_summary.json
```
**5. List Transcripts** — All transcripts for a gene, with optional
`--only-mane` or `--only-canonical` filters. Output includes Transcript Support
Level (TSL).
```bash
uv run scripts/ensembl_api.py transcripts ENSG00000141510 --only-mane --output transcripts_mane.json
uv run scripts/ensembl_api.py transcripts ENSG00000141510 --only-canonical --output transcripts_canonical.json
uv run scripts/ensembl_api.py transcripts ENSG00000141510 --output transcripts_all.json
```
**5b. Canonical TSS** — Get the single coordinate of the Transcription Start
Site (TSS) for the canonical transcript of a gene.
> [!NOTE] Unlike the standard `transcripts` command, `canonical-tss` accepts
> both symbols (e.g., `TP53`) and Ensembl IDs, and automatically resolves them.
> It also does the math for strand orientation (TSS is `Start` for `+` strand
> and `End` for `-` strand), outputting the single integer coordinate directly.
```bash
uv run scripts/ensembl_api.py canonical-tss TP53 --output tp53_tss.json
uv run scripts/ensembl_api.py canonical-tss ENSG00000141510 --output tss.json
```
**6. Transcript Structure** — Exon coordinates, CDS boundaries, and computed
5'/3' UTR regions for a transcript.
```bash
uv run scripts/ensembl_api.py transcript-structure ENST00000269305 --output structure.json
```
**7. Protein Info** — ENSP ID and sequence length for a transcript.
```bash
uv run scripts/ensembl_api.py protein-info ENST00000269305 --output protein_info.json
```
**8. Protein Sequence** — Amino acid FASTA for a transcript (ENST) or protein
(ENSP) ID.
```bash
uv run scripts/ensembl_api.py protein-sequence ENST00000269305 --output protein.fasta
uv run scripts/ensembl_api.py protein-sequence ENSP00000269305 --output protein_ensp.fasta
```
**9. Variant Consequence (VEP)** — Predict molecular consequences for a genomic
variant. Includes open-licensed plugins: AlphaMissense, Conservation,
DosageSensitivity, IntAct, MaveDB, OpenTargets, LoF (Loftee), NMD, UTRAnnotator,
mutfunc, LOEUF.
```bash
uv run scripts/ensembl_api.py vep 9:21971147:T:C --species human --output vep.json
uv run scripts/ensembl_api.py vep rs699 --species human --output vep_rs699.json
```
Example VEP stdout output:
```
[*] Variant: 9:21971147:T>C
[*] Most severe consequence: missense_variant
[*] Found 15 transcript consequences.
[*] VEP Predictions:
- ENST00000304494 (CDKN2A): Consequence = missense_variant
- ENST00000304494 (CDKN2A): Amino Acids = N/S
- ENST00000304494 (CDKN2A): SIFT = deleterious (0.01)
- ENST00000304494 (CDKN2A): AlphaMissense Class = likely_benign
- ENST00000304494 (CDKN2A): AlphaMissense Pathogenicity = 0.2129
- ENST00000304494 (CDKN2A): Conservation = 2.05
- ENST00000304494 (CDKN2A): Dosage Sensitivity (Haplo) = 0.889228328567991
- ENST00000304494 (CDKN2A): Dosage Sensitivity (Triplo) = 0.135514349094646
- ENST00000304494 (CDKN2A): Loss of Function (LOEUF) = 0.791
```
**Presenting VEP Results:** After running the VEP command, you MUST present the
full VEP Predictions list from stdout to the user. This list contains both
standard VEP predictions (Consequence, Amino Acids, SIFT, PolyPhen) and
open-license plugin results (AlphaMissense, Conservation, Dosage Sensitivity,
LOEUF, Loftee LoF, NMD, UTRAnnotator, Mutfunc). Do NOT just summarize — show the
complete list so the user can see all predictions. If the list is very long
(many transcripts), show the MANE Select / canonical transcript rows in full and
note that the complete data is in the JSON output.
## Parsing Outputs
If the user needs detailed, nested structural data (like the precise integer
coordinates of Exon 2 of a transcript) that isn't summarized in stdout:
1. Locate the JSON file (either specified via `--output` or the temporary file
path printed by the script).
2. Use terminal tools like `jq` or write a quick, disposable python snippet to
extract the specific data point requested. Do **not** attempt to read the
entire JSON file into your context if it is very large.
## Custom Queries
If you need to make an API call that the script does not support (e.g., fetching
protein domain annotations, coordinate mapping between assemblies, homology
searches, linkage disequilibrium, or phenotype lookups), read
`references/ensembl_rest_api_reference.md` for a complete reference of available
endpoints, parameters, and response fields.
**CRITICAL:** When writing custom scripts or using alternatives to the provided
scripts, you **MUST** respect the Ensembl REST API rate limits (maximum 15
requests per second) and handle `429 Too Many Requests` errors gracefully (e.g.,
with exponential backoff).