gwas-catalog
dataRetrieve variant-trait associations from the NHGRI-EBI GWAS Catalog REST API — by trait, study accession, variant or gene — together with the discovery and replication cohort ancestry, sample sizes and effect estimates that decide whether an association transfers to anyone else. Handles the HAL response shape, the 20-row default page, and p-values that underflow a float to 0.0.
GWAS Catalog — variant-trait associations with study context
The NHGRI-EBI GWAS Catalog is the curated record of what published genome-wide association studies actually reported: on 2026-08-27, 229,952 studies and 1,191,948 variant-trait associations, each carrying the paper it came from, the cohort it was measured in, and the ontology term a curator assigned to the trait.
This skill is about getting associations out of the REST API and onto disk with the study record attached, because an association without its cohort is not a result. It is a number estimated in some particular set of people, at a threshold the Catalog sets lower than genome-wide significance, and both of those facts live in the study record rather than in the association.
Nothing here needs an account, a key or a click-through. The REST API and the FTP release files are open over HTTPS. EMBL-EBI's Terms of Use state that it "places no additional restrictions on the use or redistribution of the data available via its Data Resources and Tools other than those provided by the original data owners", and the Catalog's own summary-statistics page says "The majority of data are made available either through 'CC0' or 'EMBL-EBI's standard terms of use' with a small number of exceptions" — those exceptions are per-study, flagged as a Usage License on the download tables, and matter when you redistribute a study's summary statistics rather than when you query the catalogue. Cite the Catalog papers under Sources.
Four mechanical traps cause most wrong GWAS Catalog results, and all four return a well-formed HTTP 200 rather than an error. They are the next four sections; two more — what counts as one association, and what a "mapped gene" actually is — follow after them.
Trap 1 — this is HAL, and there are three different response shapes
Every response is application/hal+json. Records are never at the top level; they are
nested under _embedded.<collectionName>, and the collection name is the plural
resource name, not results or data. Code written against a plain JSON list finds
nothing, raises nothing, and reports no hits.
The part that is easy to miss is that the same API answers in three shapes, and a single unwrap helper written against the first one breaks on the other two:
| Shape | Looks like | Returned by |
|---|---|---|
| Paged collection | _embedded + _links + page |
/studies, /associations, /efoTraits, and the searches whose URI template lists {?…,page,size,…} |
| Unpaged collection | _embedded + _links, and no page key at all |
/studies/{acc}/associations, /efoTraits/{id}/studies, and the searches whose template lists only {?…,projection} |
| Single resource | a bare object with _links, and no _embedded |
/studies/{acc}, /singleNucleotidePolymorphisms/{rsId}, and findByAccessionId |
curl -s "https://www.ebi.ac.uk/gwas/rest/api/studies?size=1" \
| python3 -c "import json,sys
d = json.load(sys.stdin)
print('top-level keys :', sorted(d))
print('collection name:', list(d['_embedded']))
print('page :', d['page'])
print('link rels :', sorted(d['_links']))"
top-level keys : ['_embedded', '_links', 'page']
collection name: ['studies']
page : {'size': 1, 'totalElements': 229952, 'totalPages': 229952, 'number': 0}
link rels : ['first', 'last', 'next', 'profile', 'search', 'self']
profile and search are always there; first, last, next and prev appear only
as the position in the result set allows, which is what makes next a usable
end-of-pages signal.
Write one reader that copes with all three, and make it say which shape it got rather than guessing:
import json, urllib.error, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
def fetch(path, **params):
"""One GET against the API. Returns the parsed HAL document."""
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
with urllib.request.urlopen(url, timeout=300) as fh:
return json.loads(fh.read())
def shape(doc):
"""Which of the three shapes came back."""
if "_embedded" not in doc:
return "single"
return "paged" if "page" in doc else "unpaged"
def records(doc, name):
"""Rows out of a collection document, or an exception naming what arrived.
An empty result still carries `_embedded`, so a missing key means the shape
is wrong -- usually a single resource where a collection was expected."""
if "_embedded" not in doc:
raise TypeError(f"single resource, not a {name} collection: {sorted(doc)}")
if name not in doc["_embedded"]:
raise KeyError(f"no '{name}' in _embedded: {list(doc['_embedded'])}")
return doc["_embedded"][name]
# 1. Paged collection.
d = fetch("studies", size=1)
print(f"{'/studies':52} {shape(d):8} {len(records(d, 'studies'))} row(s), "
f"{d['page']['totalElements']} total")
# 2. Unpaged collection -- no `page`, so `d['page']` is a KeyError, and any
# loop keyed on totalElements never runs.
d = fetch("studies/GCST012219/associations")
print(f"{'/studies/GCST012219/associations':52} {shape(d):8} "
f"{len(records(d, 'associations'))} row(s), page key present: {'page' in d}")
# 3. Single resource -- no `_embedded` at all. This is the one that silently
# returns zero rows if you reach for _embedded with .get().
d = fetch("studies/search/findByAccessionId", accessionId="GCST012219")
print(f"{'findByAccessionId':52} {shape(d):8} accessionId={d['accessionId']}")
try:
records(d, "studies")
except TypeError as e:
print(" reaching for a collection here:", e)
/studies paged 1 row(s), 229952 total
/studies/GCST012219/associations unpaged 26 row(s), page key present: False
findByAccessionId single accessionId=GCST012219
reaching for a collection here: single resource, not a studies collection: ['_links', 'accessionId', 'ancestries', 'diseaseTrait', 'fullPvalueSet', 'genotypingTechnologies', 'gxe', 'gxg', 'imputed', 'initialSampleSize', 'platforms', 'pooled', 'publicationInfo', 'qualifier', 'replicationSampleSize', 'snpCount', 'studyDesignComment', 'userRequested']
Read the URI template on /{resource}/search before writing a query. The root
document and each /search document advertise exactly which parameters an endpoint
takes, and whether page,size is among them. That is the difference between shape one
and shape two, and it is the only reliable way to know:
curl -s "https://www.ebi.ac.uk/gwas/rest/api/studies/search" \
| python3 -c "import json,sys
for rel, link in sorted(json.load(sys.stdin)['_links'].items()):
print(f\"{rel:30} {link['href'].split('/')[-1]}\")"
findByAccessionId findByAccessionId{?accessionId,projection}
findByDiseaseTrait findByDiseaseTrait{?diseaseTrait,page,size,sort,projection}
findByEfoTrait findByEfoTrait{?efoTrait,page,size,sort,projection}
findByEfoUri findByEfoUri{?uri,page,size,sort,projection}
findByFullPvalueSet findByFullPvalueSet{?fullPvalueSet,page,size,sort,projection}
findByPublicationIdPubmedId findByPublicationIdPubmedId{?pubmedId,page,size,sort,projection}
findByUserRequested findByUserRequested{?userRequested,page,size,sort,projection}
self search
Trap 2 — the default page is 20 rows, and elsewhere size does nothing at all
Two opposite failures, and which one you get depends on the endpoint you happened to pick.
On a paged endpoint the default is 20 and there is no warning. Ask for the studies
of a common trait and you get twenty of them, in a complete-looking document, with the
evidence of truncation only in page.totalElements and in whether _links carries a
next.
On an unpaged endpoint size is accepted and ignored — the whole result set comes
back regardless. Probing such an endpoint with size=1 to "look at one record" fetches
everything; associations/search/findByRsId?rsId=rs7412&size=1 returned all 1,182
associations for APOE, and took minutes.
And size on a paged endpoint is capped at 1000: ?size=5000 answers HTTP 200 with
page.size quietly rewritten to 1000, and page.totalPages recomputed to match. Code
that asks for 5000 and then loops total // 5000 times gets 1000 rows a page and stops
after a fifth of them — with no error, and with every row it did read being correct.
import json, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
def fetch_url(url):
with urllib.request.urlopen(url, timeout=300) as fh:
return json.loads(fh.read())
def fetch(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
return fetch_url(url)
def walk(path, name, per_page=100, cap=100_000, **params):
"""Follow _links.next to exhaustion, then prove nothing was left behind.
Paging is driven by the `next` link, not by incrementing an offset: the
server owns the page size and rewrites it if you ask for too much. On the
last page `next` is simply absent, which is what ends the loop."""
doc = fetch(path, size=per_page, **params)
if "page" not in doc:
return doc["_embedded"][name] # unpaged endpoint: this is all of it
total = doc["page"]["totalElements"]
if total > cap:
raise RuntimeError(f"{path}: {total} rows exceeds cap {cap}")
rows = list(doc["_embedded"][name])
while "next" in doc["_links"]:
doc = fetch_url(doc["_links"]["next"]["href"])
rows += doc["_embedded"][name]
assert len(rows) == total, f"{path}: collected {len(rows)} of {total}"
return rows
# a. The 20-row default, on a trait with far more than 20 studies.
d = fetch("studies/search/findByEfoTrait", efoTrait="type 2 diabetes mellitus")
print(f"default page : {len(d['_embedded']['studies'])} studies returned, "
f"{d['page']['totalElements']} exist, next link: {'next' in d['_links']}")
# b. `size` is honoured up to 1000 and silently rewritten above it.
for asked in (100, 1000, 5000):
p = fetch("efoTraits", size=asked)["page"]
print(f"asked size={asked:<5} -> page.size={p['size']:<5} "
f"totalPages={p['totalPages']}")
# c. Following `next` collects the lot. Small trait, so this is four quick pages.
rows = walk("studies/search/findByEfoTrait", "studies", per_page=5,
efoTrait="eye color")
print(f"walked : {len(rows)} eye-colour studies, "
f"{len({r['accessionId'] for r in rows})} distinct accessions")
# d. An unpaged endpoint ignores `size` completely -- both calls cost the same.
one = fetch("studies/GCST012219/associations", size=1)["_embedded"]["associations"]
all_ = fetch("studies/GCST012219/associations")["_embedded"]["associations"]
print(f"size=1 ignored : {len(one)} rows with size=1, {len(all_)} rows without")
default page : 20 studies returned, 343 exist, next link: True
asked size=100 -> page.size=100 totalPages=247
asked size=1000 -> page.size=1000 totalPages=25
asked size=5000 -> page.size=1000 totalPages=25
walked : 16 eye-colour studies, 16 distinct accessions
size=1 ignored : 26 rows with size=1, 26 rows without
Rules that fall out of this:
- Never read
_embeddedwithout readingpagein the same breath.len(rows) == page.totalElementsis the assertion; anextlink is the alarm. - Page by following
_links.next, not by incrementing your own offset. The server decides the page size, and it will rewrite yours without saying so. - Do not probe an unpaged endpoint with a small
size. It costs the full result set. Check the URI template first, or expect the whole thing. - The API is not fast: the global
/studiescollection returns roughly 100 study records per 40 seconds, so a 343-study trait is about three minutes and the full 229,952 studies is not a REST job at all. See Get the files for the bulk route.
Trap 3 — p-values are a mantissa and an exponent, and the float underflows to zero
Every association carries three p-value fields, and they do not agree at the extremes:
| Field | Type | At p = 2e-28539 |
|---|---|---|
pvalueMantissa |
integer, 1-9 | 2 |
pvalueExponent |
negative integer | -28539 |
pvalue |
IEEE-754 double | 0.0 |
(The same association: rs1129038-T in GCST012219.)
An IEEE-754 double cannot represent anything below about 5e-324, and the Catalog's
strongest associations are far below that — meta-analyses of hundreds of thousands of
people routinely report exponents in the hundreds or thousands. So pvalue is a
convenience field that silently becomes 0.0 exactly for the associations you most
wanted to find.
The damage is rarely a crash. 0.0 < 5e-8 is True, so a significance filter still
passes. What breaks is everything that treats zero as a value: -log10(pvalue) raises,
ranking ties every strong hit at the same place, and a "drop rows with no p-value" step
keyed on falsiness throws away the top of the table.
Reconstruct from the mantissa and exponent, and keep them as the sort key.
import json, math, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
def fetch(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
with urllib.request.urlopen(url, timeout=300) as fh:
return json.loads(fh.read())
def pkey(assoc):
"""Sortable, printable p-value that survives underflow.
Returns (exponent, mantissa) -- ascending puts the strongest first, and it
is exact, because both fields are integers. Never round-trip through the
`pvalue` double to compare two associations."""
return (assoc["pvalueExponent"], assoc["pvalueMantissa"])
def neglog10p(assoc):
"""-log10(p) from the integer fields. math.log10(assoc['pvalue']) raises
ValueError on exactly the rows this exists for."""
return -(math.log10(assoc["pvalueMantissa"]) + assoc["pvalueExponent"])
rows = fetch("studies/GCST012219/associations")["_embedded"]["associations"]
rows.sort(key=pkey)
underflowed = [a for a in rows if a["pvalue"] == 0.0]
print(f"GCST012219 (eye colour, Simcoe 2021): {len(rows)} associations, "
f"{len(underflowed)} with pvalue == 0.0")
print(f"{'risk allele':22} {'mantissa':>9} {'exponent':>9} "
f"{'pvalue (double)':>16} {'-log10(p)':>10}")
for a in rows[:3] + rows[-2:]:
allele = a["loci"][0]["strongestRiskAlleles"][0]["riskAlleleName"]
print(f"{allele:22} {a['pvalueMantissa']:>9} {a['pvalueExponent']:>9} "
f"{a['pvalue']:>16} {neglog10p(a):>10.1f}")
try:
math.log10(rows[0]["pvalue"])
except ValueError as e:
print("\nmath.log10(pvalue) on the strongest association:", e)
GCST012219 (eye colour, Simcoe 2021): 26 associations, 2 with pvalue == 0.0
risk allele mantissa exponent pvalue (double) -log10(p)
rs1129038-T 2 -28539 0.0 28538.7
rs4778218-G 1 -582 0.0 582.0
rs12203592-T 2 -321 1.996e-321 320.7
rs2351061-A 4 -8 4e-08 7.4
rs80308281-T 5 -8 5e-08 7.3
math.log10(pvalue) on the strongest association: math domain error
Note the third row. 2e-321 is representable, but only as a subnormal double, so
pvalue comes back as 1.996e-321 — the nearest double, three significant figures
short. Between roughly 1e-308 and 5e-324 the float is present and lossy rather than
absent, which is worse than zero because nothing about it looks wrong.
Sorting /associations server-side is fragile — sort locally instead. Two separate
failures, both reproduced on 2026-08-27:
sort=…,descon a p-value field returns HTTP 500./associations?sort=pvalue,desc,sort=pvalueExponent,descandsort=pvalueMantissa,descall fail with a server error. The same fields ascending are fine, as isdescon other fields (snpType,standardError) and on/studiesand/efoTraits.- A sorted page can arrive truncated, at HTTP 200. The last page of
/associationsundersort=pvalueExponent,asccame back as a200with the body cut off mid-record and unparseable — three times out of three atsize=100&page=11919, and twice out of twice atsize=1000&page=1191. The same page withoutsortparsed cleanly, 48 rows. So the status code is not the check:json.loadsraising is the only signal, and a client that swallows it loses the tail of the collection silently.
Sort the rows you have in hand instead — pkey above is exact and costs nothing.
Effect sizes come in two mutually exclusive forms, and both can be null at once.
orPerCopyNum is an odds ratio, for case-control traits; betaNum with betaUnit and
betaDirection is a linear effect, for quantitative traits, and the direction is a word
("increase" / "decrease") rather than a sign, so a beta is unsigned until you read
two fields. range is the reported confidence interval as a string, standardError the
SE. Never coerce a beta into an OR — and do not assume every row has either: GCST012219
carries a beta on 24 of its 26 associations and neither an OR nor a beta on the other
two, because the paper did not report one.
Trap 4 — traits are ontology terms, matched exactly, and they are not all EFO
The trait endpoints are named efoTraits, and the name is now the smaller half of the
story. Counting the distinct terms in gwas-efo-trait-mappings.tsv from the 2026-08-24
release — 20,047 of them, over 110,338 curated trait strings — gives:
| Prefix | Distinct terms | Ontology |
|---|---|---|
OBA |
10,027 | Ontology of Biological Attributes — measurements |
EFO |
9,032 | Experimental Factor Ontology |
MONDO |
491 | Mondo Disease Ontology |
HP |
397 | Human Phenotype Ontology |
GO, Orphanet, PATO, NCIT, UBERON, GSSO, MP, HANCESTRO |
100 between them | — |
Most terms are no longer EFO, and the uri switches namespace with the prefix:
http://www.ebi.ac.uk/efo/EFO_0003949,
http://purl.obolibrary.org/obo/MONDO_0005148,
http://www.orpha.net/ORDO/Orphanet_654. Three namespaces, twelve prefixes. Any code
that recovers an id by stripping http://www.ebi.ac.uk/efo/ hands back a full URI for
roughly a third of the mappings and every one of the 10,027 OBA terms.
The label search is an exact match on the whole string — case is ignored, and
nothing else is. Not a prefix, not a substring, no synonyms, no spelling variants, and a
trailing space is fatal. And the two "not found" answers differ: findByEfoTrait
returns HTTP 200 with an empty list, findByShortForm returns HTTP 404.
import json, urllib.error, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
def fetch(path, **params):
"""Returns (status, document). 404 is a real answer here, not a failure."""
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
try:
with urllib.request.urlopen(url, timeout=300) as fh:
return fh.status, json.loads(fh.read())
except urllib.error.HTTPError as e:
return e.code, None
def resolve_trait(label):
"""Exact trait label -> (shortForm, uri). Raises rather than guessing."""
status, doc = fetch("efoTraits/search/findByEfoTrait", trait=label)
hits = doc["_embedded"]["efoTraits"] if doc else []
if not hits:
raise LookupError(f"{label!r}: no trait with that exact label "
f"(HTTP {status} with an empty list, not an error)")
if len(hits) > 1:
raise LookupError(f"{label!r}: {[h['shortForm'] for h in hits]}")
return hits[0]["shortForm"], hits[0]["uri"]
for label in ("type 2 diabetes", # a substring is not a match
"type II diabetes mellitus", # the wording of the retired EFO term
"eye colour", # British spelling
"type 2 diabetes mellitus ", # one trailing space
"TYPE 2 DIABETES MELLITUS", # case, and only case, is ignored
"type 2 diabetes mellitus",
"body mass index",
"eye color"):
try:
sf, uri = resolve_trait(label)
print(f"{label!r:29} -> {sf:16} {uri}")
except LookupError:
print(f"{label!r:29} -> no match (HTTP 200, empty list)")
# The retired EFO id for type 2 diabetes 404s; the live term is MONDO.
for sf in ("EFO_0001360", "MONDO_0005148", "EFO_0003949"):
status, doc = fetch("efoTraits/search/findByShortForm", shortForm=sf)
label = doc["trait"] if doc else "-"
print(f"shortForm {sf:16} HTTP {status} {label}")
'type 2 diabetes' -> no match (HTTP 200, empty list)
'type II diabetes mellitus' -> no match (HTTP 200, empty list)
'eye colour' -> no match (HTTP 200, empty list)
'type 2 diabetes mellitus ' -> no match (HTTP 200, empty list)
'TYPE 2 DIABETES MELLITUS' -> MONDO_0005148 http://purl.obolibrary.org/obo/MONDO_0005148
'type 2 diabetes mellitus' -> MONDO_0005148 http://purl.obolibrary.org/obo/MONDO_0005148
'body mass index' -> EFO_0004340 http://www.ebi.ac.uk/efo/EFO_0004340
'eye color' -> EFO_0003949 http://www.ebi.ac.uk/efo/EFO_0003949
shortForm EFO_0001360 HTTP 404 -
shortForm MONDO_0005148 HTTP 200 type 2 diabetes mellitus
shortForm EFO_0003949 HTTP 200 eye color
EFO_0001360 — the id every pre-2023 pipeline used for type 2 diabetes — is gone. So:
- Resolve the label to a short form once, then key on the short form, and record which one you resolved alongside the results. A hardcoded ontology id has a shelf life.
- There is no fuzzy trait search in this API. When you do not already know the exact
label, get the vocabulary rather than guessing at it — page
/efoTraits(24,618 terms, ~1.3 s per 100-row page) and match locally, or readgwas-efo-trait-mappings.tsvfrom the release, below. - Do not construct short forms by string-surgery on
uri. Read theshortFormfield, which the API supplies for both namespaces. - The Catalog also keeps a curator-written free-text trait on the study
(
diseaseTrait.trait, e.g."Eye color", capitalised) distinct from the ontology term on the association ("eye color"). They are different fields with different vocabularies; the filegwas-efo-trait-mappings.tsvin each release maps one to the other offline, for all 110,338 curated trait strings. findByDiseaseTraitsearches that free-text study trait,findByEfoTraitthe ontology label. Picking the wrong one is a silent empty result.
An association is not a result
Everything above is plumbing. This is the part that decides whether the number you return means anything, and it is why the study record has to travel with the association rather than being fetched only when someone asks.
Three facts change the interpretation, and all three are in the study record:
- The Catalog's curation threshold is 1×10⁻⁵, not genome-wide significance.
Associations weaker than 5×10⁻⁸ are in the table by design — the weakest association
in the whole catalogue is exactly
1×10⁻⁵, and nothing beyond it is curated. Filtering at 5×10⁻⁸ is the reader's job, and the Catalog does not do it. - Ancestry of the discovery cohort.
ancestriesis a list withtypeof"initial"(discovery) or"replication", each withnumberOfIndividuals,ancestralGroupsandcountryOfRecruitment. An effect estimated in Europeans has no guaranteed transfer — allele frequency, LD structure and effect size all differ — and the great majority of the Catalog is European. - Whether anything replicated. A study with no
"replication"ancestry entry andreplicationSampleSizeofNAreported a discovery result and stopped there.
projection=associationByStudy and projection=associationByEfoTrait inline the whole
study record — ancestries included — into each association, which turns N+1 requests
into one:
import json, urllib.parse, urllib.request
from collections import Counter
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
GWAS_SIG = 5e-8
def fetch(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
with urllib.request.urlopen(url, timeout=300) as fh:
return json.loads(fh.read())
def p_of(a):
"""Exact p-value as a float where representable, else the true magnitude
as a (mantissa, exponent) pair. Only ever used for a threshold test."""
return a["pvalueMantissa"] * 10.0 ** a["pvalueExponent"]
def cohort(study, stage):
"""(individuals, ancestral groups) for one stage of a study."""
entries = [a for a in study["ancestries"] if a["type"] == stage]
n = sum(a["numberOfIndividuals"] or 0 for a in entries)
groups = {g["ancestralGroup"] for a in entries for g in a["ancestralGroups"]}
return n, groups
# One request: every association for the trait, each with its study inlined.
rows = fetch("efoTraits/EFO_0003949/associations",
projection="associationByEfoTrait")["_embedded"]["associations"]
studies = {a["study"]["accessionId"]: a["study"] for a in rows}
print(f"eye color (EFO_0003949): {len(rows)} associations across "
f"{len(studies)} studies, in ONE request")
weak = [a for a in rows if p_of(a) > GWAS_SIG]
print(f"weaker than 5e-8 : {len(weak)} of {len(rows)} "
f"({100 * len(weak) / len(rows):.0f}%) -- curated at p < 1e-5, not 5e-8")
disc = Counter()
for s in studies.values():
n, groups = cohort(s, "initial")
for g in groups:
disc[g] += n
total = sum(disc.values())
print("\ndiscovery individuals by ancestral group")
for g, n in disc.most_common():
print(f" {g:24} {n:>9,} {100 * n / total:5.1f}%")
no_rep = [acc for acc, s in sorted(studies.items()) if not cohort(s, "replication")[0]]
print(f"\nstudies with no replication stage: {len(no_rep)} of {len(studies)}")
print(" " + ", ".join(no_rep))
# The sub-phenotype hides in a free-text field on the ASSOCIATION, not the trait.
descs = Counter(repr(a["pvalueDescription"]) for a in rows)
print("\npvalueDescription, top 6 of "
f"{len(descs)} distinct values on one ontology trait")
for d, n in descs.most_common(6):
print(f" {d:36} {n:>4}")
eye color (EFO_0003949): 201 associations across 16 studies, in ONE request
weaker than 5e-8 : 41 of 201 (20%) -- curated at p < 1e-5, not 5e-8
discovery individuals by ancestral group
European 199,091 98.7%
East Asian 2,616 1.3%
studies with no replication stage: 8 of 16
GCST000710, GCST005093, GCST005094, GCST005095, GCST005096, GCST007489, GCST90026489, GCST90255684
pvalueDescription, top 6 of 13 distinct values on one ontology trait
None 138
'(EA)' 22
'(Eye color)' 8
' ' 6
'(eye color)' 5
'(Eye color phototype score)' 5
What that output is actually telling you, and what to carry into a write-up:
- A fifth of the rows are below genome-wide significance. Report the threshold you
applied, every time.
p_of()above is the only safe way to test it, becausea["pvalue"]is 0.0 for the strongest rows. - 98.7% of the discovery sample is European, and the entire non-European remainder
is one East Asian cohort. Any statement of the form "variant X raises risk of Y"
inherits that, and a polygenic score built on it is known to lose accuracy out of
ancestry. Where an
ancestralGroupis"NR"it means not reported — treat it as missing, not as a category, and do not let it dilute a percentage silently. - Eight of sixteen studies have no replication stage. That is not a defect; it is what a discovery paper looks like. It is a defect to report it as though it had one.
- One ontology trait is several phenotypes.
pvalueDescriptionis free text on the association, and thirteen distinct values sit under the single term "eye color" —'(EA)','(Eye color phototype score)', and elsewhere in the same 201 rows'(Blue eye color)','(Brown eye color)'and'(Intermediate/Green eye color)', which take opposite effect directions at the same locus. It is inconsistently cased ('(Eye color)'and'(eye color)'are both present),Noneon 138 rows, and on six rows it is a single space rather than empty. Strip it, read it, and never pool effect sizes across distinct values.
Two further study-record fields worth carrying: fullPvalueSet (whether full summary
statistics were deposited, which decides if fine-mapping or Mendelian randomisation is
even possible) and snpCount with imputed (how many variants were tested, which sets
the real multiple-testing burden).
Variants, mapped genes and genomic context
/singleNucleotidePolymorphisms/{rsId} gives the variant's mapped position, cytogenetic
region, functional class and the genes around it. Two things about genomicContexts are
easy to get wrong:
- It is mapped twice, once from Ensembl and once from NCBI, so every gene appears more than once and the two sources disagree on distances. Counting rows counts sources.
isClosestGeneis not unique, and it is not the gene the variant is in. It is flagged per source and per direction — nearest upstream and nearest downstream both get it — so a single variant can have four rows marked closest. The gene the variant actually sits inside is the one withisIntergenic: false, and it is flaggedisClosestGene: false.
import json, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
def fetch(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
with urllib.request.urlopen(url, timeout=300) as fh:
return json.loads(fh.read())
snp = fetch("singleNucleotidePolymorphisms/rs1129038")
loc = snp["locations"][0]
print(f"{snp['rsId']} chr{loc['chromosomeName']}:{loc['chromosomePosition']} "
f"{loc['region']['name']} {snp['functionalClass']} merged={snp['merged']}")
ctx = snp["genomicContexts"]
print(f"\n{len(ctx)} genomicContexts rows, "
f"{len({c['gene']['geneName'] for c in ctx})} distinct genes, "
f"sources {sorted({c['source'] for c in ctx})}")
print(f"{'gene':10} {'source':9} {'distance':>9} {'intergenic':>11} "
f"{'upstream':>9} {'closest':>8}")
for c in sorted(ctx, key=lambda c: (c["gene"]["geneName"], c["source"])):
print(f"{c['gene']['geneName']:10} {c['source']:9} {c['distance']:>9} "
f"{str(c['isIntergenic']):>11} {str(c['isUpstream']):>9} "
f"{str(c['isClosestGene']):>8}")
# Genes the variant is actually inside, deduplicated across the two sources.
inside = sorted({c["gene"]["geneName"] for c in ctx if not c["isIntergenic"]})
closest = sorted({c["gene"]["geneName"] for c in ctx if c["isClosestGene"]})
print(f"\nvariant lies within : {inside}")
print(f"flagged closest : {closest}")
rs1129038 chr15:28111713 15q13.1 3_prime_UTR_variant merged=0
6 genomicContexts rows, 3 distinct genes, sources ['Ensembl', 'NCBI']
gene source distance intergenic upstream closest
HERC2 Ensembl 0 False False False
HERC2 NCBI 0 False False False
OCA2 Ensembl 12049 True True True
OCA2 NCBI 12398 True True True
RPL41P2 Ensembl 203938 True False True
RPL41P2 NCBI 203940 True False True
variant lies within : ['HERC2']
flagged closest : ['OCA2', 'RPL41P2']
That is the textbook eye-colour locus, and every part of the trap is visible in six
rows. The variant is in the 3′ UTR of HERC2; the gene it regulates is OCA2, about
12 kb upstream; and RPL41P2 — a ribosomal protein pseudogene 204 kb downstream — is
flagged isClosestGene just as OCA2 is, because each direction gets its own nearest
gene. Ensembl and NCBI put OCA2 at 12,049 and 12,398 bp respectively, so even the
distance depends on which row you read. Taking genomicContexts[0], or filtering on
isClosestGene alone, returns a pseudogene as the answer to "which gene is this?"
Lookups by variant, gene and region all hang off the same resource, and they split
across both collection shapes — which is why the template, not the name, tells you what
you are getting. On /singleNucleotidePolymorphisms/search, findByRsId,
findByBpLocation, findByEfoTrait and findByDiseaseTrait are unpaged and return
the whole result set; findByGene, findByChromBpLocationRange,
findIdsByLocations… and findByPubmedId are paged and default to 20.
And the two "not found" answers differ again: an rsID the Catalog does not hold gives
HTTP 200 with an empty list from associations/search/findByRsId?rsId=…, but 404
from /singleNucleotidePolymorphisms/{rsId}. Same fact, two signals, and only one of
them raises.
Not every association is one variant
loci[0]["strongestRiskAlleles"][0] is the obvious way to read a risk allele, and it is
right for almost every row and silently wrong for two whole classes of them. They are
rare — of 1,000 consecutive /associations rows sampled on 2026-08-27, 999 were single
variants and one was a 3-SNP haplotype — which is exactly why nobody writes the code
that handles them, and why the error survives review. The Catalog records three forms,
distinguished by two booleans on the association:
| Form | multiSnpHaplotype |
snpInteraction |
Shape |
|---|---|---|---|
| Single variant | false |
false |
one locus, one risk allele |
| Haplotype | true |
false |
one locus, N risk alleles, haplotypeSnpCount = N |
| SNP × SNP interaction | false |
true |
N loci, one risk allele each |
So the two failure modes point in opposite directions: the haplotype hides its extra
variants inside loci[0], and the interaction hides them in the loci you never
looked at. Both come back as a perfectly ordinary association row.
import json, textwrap, urllib.parse, urllib.request
from collections import Counter
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
def fetch(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
with urllib.request.urlopen(url, timeout=300) as fh:
return json.loads(fh.read())
def risk_alleles(assoc):
"""Every risk allele in an association, flattened across loci, in order."""
return [a["riskAlleleName"] for locus in assoc["loci"]
for a in locus["strongestRiskAlleles"]]
def form(assoc):
if assoc["snpInteraction"]:
return "interaction"
return "haplotype" if assoc["multiSnpHaplotype"] else "single"
for acc in ("GCST012219", "GCST000071", "GCST006494"):
rows = fetch(f"studies/{acc}/associations")["_embedded"]["associations"]
counts = Counter(form(r) for r in rows)
# Ties broken on the allele list itself, so this picks the same row every
# run -- the API does not promise an order and max() keeps the first max.
widest = max(rows, key=lambda r: (len(risk_alleles(r)), risk_alleles(r)))
alleles = risk_alleles(widest)
locus = widest["loci"][0]
print(f"{acc}: {len(rows)} associations {dict(counts)}")
print(f" widest row : {len(widest['loci'])} loci, {len(alleles)} risk "
f"alleles, description {locus['description']!r}, "
f"haplotypeSnpCount {locus['haplotypeSnpCount']}")
print(f" the naive read : {alleles[0]}")
print(" what is there :")
print(textwrap.fill(", ".join(alleles), width=70,
initial_indent=" ", subsequent_indent=" "))
GCST012219: 26 associations {'single': 26}
widest row : 1 loci, 1 risk alleles, description 'Single variant', haplotypeSnpCount None
the naive read : rs9971729-C
what is there :
rs9971729-C
GCST000071: 3 associations {'haplotype': 2, 'single': 1}
widest row : 1 loci, 16 risk alleles, description '17 marker haplotype-1', haplotypeSnpCount 16
the naive read : rs11209003-?
what is there :
rs11209003-?, rs11209002-?, rs2064689-?, rs1004819-?, rs2902440-?,
rs11465802-?, rs2201841-?, rs11465804-?, rs11209026-?,
rs1343151-?, rs10889676-?, rs10889677-?, rs9988642-?,
rs12567232-?, rs6669582-?, rs10789230-?
GCST006494: 7 associations {'interaction': 7}
widest row : 2 loci, 2 risk alleles, description 'SNP x SNP interaction', haplotypeSnpCount None
the naive read : rs10830963-G
what is there :
rs10830963-G, rs73659517-G
Reporting rs11209003-? as the variant for that row turns a 16-marker haplotype into
a single-SNP claim, and reporting rs10830963-G alone turns an interaction into a main
effect. Three further details visible in that output:
-?means the risk allele was not reported. The suffix onriskAlleleNameis the allele, and?is a real value in the data — do not parse it as a base.descriptionandhaplotypeSnpCountcan disagree. That row is labelled'17 marker haplotype-1'and carries 16 alleles withhaplotypeSnpCount: 16. Trust the array you can count, not the curator's free-text label.authorReportedGeneslives on the locus, not the association, so an interaction has one gene list per locus. Flatten across loci, asrisk_allelesdoes above.
Get the files
Everything above prints. This writes, because a target dossier, a figure or a meta-analysis needs files — and because the catalogue grows every fortnight, so a directory of TSVs with no release stamp cannot be compared against a later pull.
Route A — the REST API, for one trait. One request for the associations with the study record inlined, then flatten to two TSVs plus a manifest recording what was resolved and when.
import csv, datetime, json, os, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
TRAIT_LABEL = "eye color"
OUT = "Data/gwas-catalog"
GWAS_SIG = 5e-8
def fetch(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
with urllib.request.urlopen(url, timeout=300) as fh:
return url, json.loads(fh.read())
def write_tsv(path, rows, cols):
with open(path, "w", newline="") as fh:
w = csv.DictWriter(fh, fieldnames=cols, delimiter="\t",
extrasaction="ignore")
w.writeheader()
w.writerows(rows)
return {"path": path, "rows": len(rows), "bytes": os.path.getsize(path)}
os.makedirs(OUT, exist_ok=True)
sources, files = {}, []
# 1. Exact label -> ontology short form. Never hardcode the id.
url, doc = fetch("efoTraits/search/findByEfoTrait", trait=TRAIT_LABEL)
hits = doc["_embedded"]["efoTraits"]
if len(hits) != 1:
raise LookupError(f"{TRAIT_LABEL!r}: {[h['shortForm'] for h in hits]}")
trait = hits[0]
sources["trait"] = url
print("resolved:", trait["shortForm"], "|", trait["uri"])
# 2. One request for every association, study record inlined.
url, doc = fetch(f"efoTraits/{trait['shortForm']}/associations",
projection="associationByEfoTrait")
assocs = doc["_embedded"]["associations"]
sources["associations"] = url
# 3. Associations, flat. Mantissa and exponent are the p-value of record; the
# double is carried too, so a reader can see where it underflowed.
arows = []
for a in assocs:
# Flatten across every locus: a haplotype is one locus with N alleles, an
# interaction is N loci with one each, and taking loci[0][0] loses both.
alleles = [al["riskAlleleName"] for L in a["loci"]
for al in L["strongestRiskAlleles"]]
genes = sorted({g["geneName"] for L in a["loci"]
for g in L["authorReportedGenes"]})
p = a["pvalueMantissa"] * 10.0 ** a["pvalueExponent"]
arows.append({
"studyAccession": a["study"]["accessionId"],
"associationForm": ("interaction" if a["snpInteraction"] else
"haplotype" if a["multiSnpHaplotype"] else "single"),
"lociCount": len(a["loci"]),
"riskAlleleCount": len(alleles),
"riskAlleles": ";".join(alleles),
"riskFrequency": a["riskFrequency"],
"authorReportedGenes": ",".join(genes),
"pvalueMantissa": a["pvalueMantissa"],
"pvalueExponent": a["pvalueExponent"],
"pvalueDouble": a["pvalue"],
"genomeWideSignificant": p <= GWAS_SIG,
"pvalueDescription": (a["pvalueDescription"] or "").strip(),
"orPerCopyNum": a["orPerCopyNum"],
"betaNum": a["betaNum"],
"betaUnit": a["betaUnit"],
"betaDirection": a["betaDirection"],
"standardError": a["standardError"],
"range": a["range"],
"snpType": a["snpType"],
})
arows.sort(key=lambda r: (r["pvalueExponent"], r["pvalueMantissa"]))
files.append(write_tsv(f"{OUT}/{trait['shortForm']}_associations.tsv", arows,
list(arows[0])))
# 4. Studies, one row each, with the cohort that decides interpretation.
srows = []
for acc, s in sorted({a["study"]["accessionId"]: a["study"]
for a in assocs}.items()):
stage = {}
for st in ("initial", "replication"):
e = [x for x in s["ancestries"] if x["type"] == st]
stage[st] = (sum(x["numberOfIndividuals"] or 0 for x in e),
sorted({g["ancestralGroup"] for x in e
for g in x["ancestralGroups"]}))
srows.append({
"studyAccession": acc,
"pubmedId": s["publicationInfo"]["pubmedId"],
"firstAuthor": s["publicationInfo"]["author"]["fullname"],
"publicationDate": s["publicationInfo"]["publicationDate"],
"journal": s["publicationInfo"]["publication"],
"diseaseTrait": s["diseaseTrait"]["trait"],
"discoveryN": stage["initial"][0],
"discoveryAncestry": "|".join(stage["initial"][1]),
"replicationN": stage["replication"][0],
"replicationAncestry": "|".join(stage["replication"][1]),
"initialSampleSize": s["initialSampleSize"],
"replicationSampleSize": s["replicationSampleSize"],
"snpCount": s["snpCount"],
"imputed": s["imputed"],
"fullPvalueSet": s["fullPvalueSet"],
"genotypingTechnology": "|".join(
t["genotypingTechnology"] for t in s["genotypingTechnologies"]),
})
files.append(write_tsv(f"{OUT}/{trait['shortForm']}_studies.tsv", srows,
list(srows[0])))
# 5. The manifest. Without it these numbers cannot be compared to a later pull.
_, root = fetch("studies", size=1)
manifest = {
"retrieved": datetime.date.today().isoformat(),
"traitLabel": TRAIT_LABEL,
"traitShortForm": trait["shortForm"],
"traitUri": trait["uri"],
"catalogueStudiesTotal": root["page"]["totalElements"],
"significanceThreshold": GWAS_SIG,
"associationsRetrieved": len(arows),
"associationsAtThreshold": sum(r["genomeWideSignificant"] for r in arows),
"studiesRetrieved": len(srows),
"sources": sources,
"files": files,
}
with open(f"{OUT}/manifest.json", "w") as fh:
json.dump(manifest, fh, indent=2)
for f in files:
print(f"{f['path']:52} {f['rows']:>4} rows {f['bytes']:>7,} bytes")
print(f"{OUT + '/manifest.json':52} "
f"{manifest['associationsAtThreshold']} of "
f"{manifest['associationsRetrieved']} at p <= 5e-8")
resolved: EFO_0003949 | http://www.ebi.ac.uk/efo/EFO_0003949
Data/gwas-catalog/EFO_0003949_associations.tsv 201 rows 24,094 bytes
Data/gwas-catalog/EFO_0003949_studies.tsv 16 rows 3,492 bytes
Data/gwas-catalog/manifest.json 160 of 201 at p <= 5e-8
Route B — the release files, for anything catalogue-wide. At roughly half a second
per study record the REST API cannot answer questions about all 229,952 studies. Each
fortnightly release publishes the whole catalogue as flat files over anonymous HTTPS at
https://ftp.ebi.ac.uk/pub/databases/gwas/releases/latest/:
| File | Size | What it carries |
|---|---|---|
gwas-catalog-download-studies-v1.0.3.1.txt |
123 MB | every study, with STUDY ACCESSION, MAPPED_TRAIT_URI, COHORT, FULL SUMMARY STATISTICS, SUMMARY STATS LOCATION |
gwas-catalog-ancestry.tsv |
49 MB | one row per study per ancestry stage — the cohort table |
gwas-catalog-associations_ontology-annotated-full.zip |
70 MB | every association, ontology-annotated |
gwas-efo-trait-mappings.tsv |
19 MB | free-text study trait → ontology term and parent |
Stream it rather than loading it, and join on the accessions Route A already collected. This block downloads about 49 MB and takes a few minutes on a slow link:
import csv, io, json, os, urllib.request
from collections import Counter
FTP = "https://ftp.ebi.ac.uk/pub/databases/gwas/releases/latest"
OUT = "Data/gwas-catalog"
wanted = {r["studyAccession"] for r in
csv.DictReader(open(f"{OUT}/EFO_0003949_studies.tsv"), delimiter="\t")}
print(f"{len(wanted)} accessions from Route A")
# The header is misspelled upstream -- NUMBER OF INDIVDUALS, ADDITONAL ANCESTRY
# DESCRIPTION. Read the header row; never hand-type these names.
url = f"{FTP}/gwas-catalog-ancestry.tsv"
rows, seen, stages = [], 0, Counter()
with urllib.request.urlopen(url, timeout=600) as fh:
reader = csv.DictReader(io.TextIOWrapper(fh, encoding="utf-8"), delimiter="\t")
header = reader.fieldnames
for rec in reader:
seen += 1
stages[rec["STAGE"]] += 1
if rec["STUDY ACCESSION"] in wanted:
rows.append(rec)
print(f"streamed {seen:,} ancestry rows for the whole catalogue")
print(f"stages: {dict(stages.most_common())}")
print(f"kept {len(rows)} rows for {len({r['STUDY ACCESSION'] for r in rows})} studies")
with open(f"{OUT}/EFO_0003949_ancestry.tsv", "w", newline="") as fh:
w = csv.DictWriter(fh, fieldnames=header, delimiter="\t")
w.writeheader()
w.writerows(sorted(rows, key=lambda r: (r["STUDY ACCESSION"], r["STAGE"])))
with open(f"{OUT}/manifest.json") as fh:
manifest = json.load(fh)
manifest["bulkRelease"] = {"url": url, "ancestryRowsInRelease": seen,
"ancestryRowsKept": len(rows), "columns": header}
with open(f"{OUT}/manifest.json", "w") as fh:
json.dump(manifest, fh, indent=2)
path = f"{OUT}/EFO_0003949_ancestry.tsv"
print(f"{path:52} {len(rows):>4} rows {os.path.getsize(path):>7,} bytes")
16 accessions from Route A
streamed 298,269 ancestry rows for the whole catalogue
stages: {'initial': 260878, 'replication': 37391}
kept 26 rows for 16 studies
Data/gwas-catalog/EFO_0003949_ancestry.tsv 26 rows 4,452 bytes
That output settles something the per-trait view only hints at. Across all 229,952
studies the release carries 260,878 initial rows and 37,391 replication rows —
and since a study with a replication stage contributes at least one row, at most one
study in six replicated anything at all. The eight-of-sixteen seen above for eye
colour is the ordinary case, not an unlucky draw, and "was it replicated?" is a question
worth asking of every association you plan to report.
Full summary statistics, when a study deposited them. Studies with
fullPvalueSet: true have a directory under
https://ftp.ebi.ac.uk/pub/databases/gwas/summary_statistics/, bucketed a thousand
accessions at a time — GCST012001-GCST013000/GCST012219/. Do not build the filename
from the accession: submissions since the GWAS-SSF standard are {ACCESSION}.tsv.gz
with a -meta.yaml beside them, but legacy studies are not
(GCST000028/Saxena-17463246.txt, and GCST000510 holds fifteen separate .assoc.gz
files, one per blood trait). List the directory, read the -meta.yaml for the genome
assembly and imputation panel, and check md5sum.txt. A study with
fullPvalueSet: false has no directory at all, and the URL 404s.
Try it
A self-contained check that this skill still holds. Public data, no account, no key,
Python standard library only — no pip install step, and nothing written to disk.
Data — the GWAS Catalog REST API, reached through nine endpoints:
https://www.ebi.ac.uk/gwas/rest/api/studies?size=1
https://www.ebi.ac.uk/gwas/rest/api/studies/search/findByAccessionId?accessionId=GCST012219
https://www.ebi.ac.uk/gwas/rest/api/studies/GCST012219/associations
https://www.ebi.ac.uk/gwas/rest/api/studies/GCST000071/associations
https://www.ebi.ac.uk/gwas/rest/api/studies/GCST006494/associations
https://www.ebi.ac.uk/gwas/rest/api/studies/search/findByEfoTrait?efoTrait=eye%20color&size=5
https://www.ebi.ac.uk/gwas/rest/api/efoTraits/search/findByEfoTrait?trait=type%202%20diabetes%20mellitus
https://www.ebi.ac.uk/gwas/rest/api/efoTraits/search/findByShortForm?shortForm=EFO_0001360
https://www.ebi.ac.uk/gwas/rest/api/efoTraits/EFO_0003949/associations?projection=associationByEfoTrait
The Catalog is open EMBL-EBI infrastructure; no account or licence acceptance is needed
to read it. Each probe is chosen for a trap rather than for being typical.
GCST012219 (Simcoe et al. 2021, PMID 33692100, eye colour in 192,986 people)
exercises three at once: it is reached through an unpaged sub-resource, its strongest
association is 2×10⁻²⁸⁵³⁹ and therefore underflows the pvalue double to 0.0, and
it has no replication stage. EFO_0003949 and MONDO_0005148 probe the trait
vocabulary, EFO_0001360 is the retired id that now 404s, and GCST000071 (a
16-marker haplotype) and GCST006494 (SNP × SNP interactions) are the two
association shapes that are not one variant. Last confirmed reachable 2026-08-27.
Budget about two minutes. Most calls answer in a second or two; the last one — every association for a trait, with the study record inlined — took roughly 70 seconds, which is the price of the one request that replaces seventeen.
Run — this block asserts the shape of the API, not only its values, so a renamed field or restructured pagination fails it loudly:
import json, math, urllib.error, urllib.parse, urllib.request
BASE = "https://www.ebi.ac.uk/gwas/rest/api"
GWAS_SIG = 5e-8
DOUBLE_FLOOR = -324 # below this an IEEE-754 double is exactly 0.0
def fetch(url):
"""Returns (status, document). 404 is an answer here, not a failure."""
try:
with urllib.request.urlopen(url, timeout=300) as fh:
return fh.status, json.loads(fh.read())
except urllib.error.HTTPError as e:
return e.code, None
def get(path, **params):
url = f"{BASE}/{path}"
if params:
url += "?" + urllib.parse.urlencode(params, doseq=True)
return fetch(url)
# 1. SHAPE: a paged collection is _embedded + _links + page, and nothing else.
status, doc = get("studies", size=1)
assert status == 200, status
assert sorted(doc) == ["_embedded", "_links", "page"], sorted(doc)
assert list(doc["_embedded"]) == ["studies"], list(doc["_embedded"])
assert sorted(doc["page"]) == ["number", "size", "totalElements", "totalPages"], \
sorted(doc["page"])
assert doc["page"]["size"] == 1 and doc["page"]["number"] == 0, doc["page"]
assert "next" in doc["_links"], sorted(doc["_links"])
n_studies = doc["page"]["totalElements"]
print(f"1 paged shape : _embedded.studies + page, {n_studies:,} studies")
# 2. SHAPE: findByAccessionId is a SINGLE resource -- no _embedded at all.
status, doc = get("studies/search/findByAccessionId", accessionId="GCST012219")
assert status == 200 and "_embedded" not in doc, sorted(doc)
assert doc["accessionId"] == "GCST012219", doc["accessionId"]
assert doc["fullPvalueSet"] is False and doc["imputed"] is True
print(f"2 single resource : GCST012219, {doc['snpCount']:,} variants tested, "
f"fullPvalueSet={doc['fullPvalueSet']}")
# 3. SHAPE: a sub-resource is UNPAGED -- no `page` key, and `size` is ignored.
status, one = get("studies/GCST012219/associations", size=1)
status2, all_ = get("studies/GCST012219/associations")
assert status == status2 == 200
assert "page" not in one and "page" not in all_, "sub-resource grew a page key"
rows = all_["_embedded"]["associations"]
assert len(one["_embedded"]["associations"]) == len(rows), "size is now honoured"
print(f"3 unpaged shape : size=1 returned all {len(rows)} associations")
# 4. INVARIANT: mantissa x 10^exponent is the p-value of record, and the
# `pvalue` double is exactly 0.0 for everything below the double floor.
for a in rows:
assert 1 <= a["pvalueMantissa"] <= 9, a["pvalueMantissa"]
assert a["pvalueExponent"] < 0, a["pvalueExponent"]
if a["pvalueExponent"] <= DOUBLE_FLOOR:
assert a["pvalue"] == 0.0, (a["pvalueMantissa"], a["pvalueExponent"])
else:
assert math.isclose(a["pvalue"],
a["pvalueMantissa"] * 10.0 ** a["pvalueExponent"],
rel_tol=1e-9), a["pvalue"]
best = min(rows, key=lambda a: (a["pvalueExponent"], a["pvalueMantissa"]))
under = [a for a in rows if a["pvalue"] == 0.0]
assert best["pvalue"] == 0.0 and best in under
print(f"4 p-value underflow: strongest {best['pvalueMantissa']}e"
f"{best['pvalueExponent']} -> pvalue {best['pvalue']}, "
f"{len(under)} of {len(rows)} underflowed")
# 5. INVARIANT: paging is driven by _links.next, and the walk is exhaustive.
status, doc = get("studies/search/findByEfoTrait", efoTrait="eye color", size=5)
total, seen, pages = doc["page"]["totalElements"], [], 0
while True:
pages += 1
seen += doc["_embedded"]["studies"]
if "next" not in doc["_links"]:
break
status, doc = fetch(doc["_links"]["next"]["href"])
assert status == 200, status
assert len(seen) == total, f"collected {len(seen)} of {total}"
assert pages == math.ceil(total / 5), pages
print(f"5 next-link paging : {pages} pages of 5 -> {len(seen)} of {total} studies")
# 6. INVARIANT: the trait label search is exact, and the two misses differ.
status, doc = get("efoTraits/search/findByEfoTrait", trait="type 2 diabetes")
assert status == 200 and doc["_embedded"]["efoTraits"] == [], "prefix match now?"
status, doc = get("efoTraits/search/findByEfoTrait",
trait="type 2 diabetes mellitus")
hit = doc["_embedded"]["efoTraits"][0]
assert hit["shortForm"] == "MONDO_0005148", hit["shortForm"]
assert hit["uri"].endswith("MONDO_0005148") and "purl.obolibrary.org" in hit["uri"]
retired, _ = get("efoTraits/search/findByShortForm", shortForm="EFO_0001360")
live, _ = get("efoTraits/search/findByShortForm", shortForm="EFO_0003949")
assert (retired, live) == (404, 200), (retired, live)
print(f"6 trait vocabulary : 'type 2 diabetes' -> 0 hits; exact label -> "
f"{hit['shortForm']}; EFO_0001360 -> HTTP {retired}")
# 7. INVARIANT: the projection inlines the study, so ancestry travels with the
# association -- and the Catalog curates below genome-wide significance.
status, doc = get("efoTraits/EFO_0003949/associations",
projection="associationByEfoTrait")
assoc = doc["_embedded"]["associations"]
assert all("study" in a and "ancestries" in a["study"] for a in assoc)
studies = {a["study"]["accessionId"]: a["study"] for a in assoc}
weak = [a for a in assoc
if a["pvalueMantissa"] * 10.0 ** a["pvalueExponent"] > GWAS_SIG]
assert weak, "no sub-threshold associations -- curation threshold changed?"
euro = sum(x["numberOfIndividuals"] or 0
for s in studies.values() for x in s["ancestries"]
if x["type"] == "initial"
and any(g["ancestralGroup"] == "European" for g in x["ancestralGroups"]))
disc = sum(x["numberOfIndividuals"] or 0
for s in studies.values() for x in s["ancestries"]
if x["type"] == "initial")
print(f"7 study context : {len(assoc)} associations / {len(studies)} studies, "
f"{len(weak)} weaker than 5e-8, {100 * euro / disc:.0f}% European discovery")
# 8. INVARIANT: an association is not always one variant, and the two
# multi-variant forms hide their extra alleles in different places.
status, doc = get("studies/GCST000071/associations")
hap = [a for a in doc["_embedded"]["associations"] if a["multiSnpHaplotype"]]
assert hap, "GCST000071 no longer carries a haplotype association"
wide = max(hap, key=lambda a: len(a["loci"][0]["strongestRiskAlleles"]))
n_alleles = len(wide["loci"][0]["strongestRiskAlleles"])
assert len(wide["loci"]) == 1 and n_alleles > 1, (len(wide["loci"]), n_alleles)
assert wide["loci"][0]["haplotypeSnpCount"] == n_alleles, wide["loci"][0]
status, doc = get("studies/GCST006494/associations")
inter = [a for a in doc["_embedded"]["associations"] if a["snpInteraction"]]
assert inter, "GCST006494 no longer carries an interaction association"
assert all(len(a["loci"]) > 1 for a in inter), "interaction is no longer multi-locus"
assert all(len(L["strongestRiskAlleles"]) == 1
for a in inter for L in a["loci"]), "interaction locus grew alleles"
print(f"8 multi-variant : haplotype = 1 locus x {n_alleles} alleles; "
f"interaction = {len(inter[0]['loci'])} loci x 1 allele")
Expect — two kinds of line, and they fail for different reasons.
Invariants. These hold whatever the release, so a failure means the skill is
wrong, not stale: three response shapes distinguished by the presence of page and
_embedded; page carrying exactly number, size, totalElements, totalPages;
size ignored on a sub-resource; pvalueMantissa in 1–9 with a negative exponent;
pvalue exactly 0.0 at or below exponent −324 and matching the reconstruction above
it; a _links.next walk collecting exactly totalElements; exact-match trait labels;
HTTP 404 from findByShortForm versus HTTP 200 with an empty list from
findByEfoTrait; EFO_0001360 retired and type 2 diabetes living at
MONDO_0005148; every associationByEfoTrait row carrying an inlined study with
ancestries; at least one association weaker than 5×10⁻⁸; and a haplotype association
being one locus with many risk alleles while an interaction is many loci with one each.
Observed values, 2026-08-27. These move as the Catalog grows — a mismatch is drift to investigate, not a bug. Study and association counts rise every fortnight; the eye colour trait had 201 associations across 16 studies; GCST012219 had 26 associations, 2 of them underflowed; the widest haplotype in GCST000071 carried 16 risk alleles.
1 paged shape : _embedded.studies + page, 229,952 studies
2 single resource : GCST012219, 11,532,091 variants tested, fullPvalueSet=False
3 unpaged shape : size=1 returned all 26 associations
4 p-value underflow: strongest 2e-28539 -> pvalue 0.0, 2 of 26 underflowed
5 next-link paging : 4 pages of 5 -> 16 of 16 studies
6 trait vocabulary : 'type 2 diabetes' -> 0 hits; exact label -> MONDO_0005148; EFO_0001360 -> HTTP 404
7 study context : 201 associations / 16 studies, 41 weaker than 5e-8, 99% European discovery
8 multi-variant : haplotype = 1 locus x 16 alleles; interaction = 2 loci x 1 allele
Limits worth stating in a write-up
- The Catalog records what was published, not what is true. It is a curated literature index. An association appears because a paper reported it at 1×10⁻⁵ or stronger, not because anyone re-analysed the data.
- It is overwhelmingly European. Say so with a number from the study records you actually retrieved, not as a generic caveat.
- Curated associations are the reported lead variant, not a fine-mapped causal one.
The risk allele is in linkage disequilibrium with whatever is causal, and
authorReportedGenesandgenomicContextsare proximity, not mechanism. - Effect sizes are not comparable across studies without care — different
covariates, different phenotype definitions, different reference alleles, and
pvalueDescriptionsub-phenotypes pooled under one ontology term. - Sample sizes overlap. Many studies draw on the same biobanks, so summing
numberOfIndividualsacross studies of a trait counts people more than once. The 99% figure above is a composition estimate, not a headcount. - State the release. A new release lands roughly every fortnight and the totals move with it — 229,952 studies and 1,191,948 associations on 2026-08-27, against a release dated 2026-08-24. Any count you report is a dated observation, not a constant, which is what the manifest in Get the files exists to record.
Sources
- GWAS Catalog — https://www.ebi.ac.uk/gwas/
- REST API root, self-describing — https://www.ebi.ac.uk/gwas/rest/api/
- API documentation — https://www.ebi.ac.uk/gwas/rest/docs/api
- Release files over HTTPS — https://ftp.ebi.ac.uk/pub/databases/gwas/releases/latest/
- Summary statistics, and the per-study Usage License column — https://www.ebi.ac.uk/gwas/downloads/summary-statistics
- EMBL-EBI Terms of Use — https://www.ebi.ac.uk/about/terms-of-use
- Hindorff et al. (2009) PNAS 106, 9362-9367 — https://doi.org/10.1073/pnas.0903103106
- Welter et al. (2014) Nucleic Acids Res 42, D1001-D1006 — https://doi.org/10.1093/nar/gkt1229
- Buniello et al. (2019) Nucleic Acids Res 47, D1005-D1012 — https://doi.org/10.1093/nar/gky1120
- Sollis et al. (2023) Nucleic Acids Res 51, D977-D985 — https://doi.org/10.1093/nar/gkac1010
- Cerezo et al. (2025) Nucleic Acids Res 53, D998-D1005 — https://doi.org/10.1093/nar/gkae1070
The GWAS Catalog is a collaboration between the National Human Genome Research Institute and EMBL-EBI. Cite Cerezo et al. (2025) for the resource, and cite the original study behind any association you report — the PubMed id travels with every study record for exactly that reason.