Metadata-Version: 2.4
Name: afdb-query
Version: 0.4.0
Summary: Sequence-based programmatic access to the AlphaFold Protein Structure Database
License: MIT
Requires-Python: >=3.10
Requires-Dist: httpx>=0.27
Provides-Extra: dev
Requires-Dist: pytest-cov>=5; extra == 'dev'
Requires-Dist: pytest>=8; extra == 'dev'
Requires-Dist: respx>=0.21; extra == 'dev'
Requires-Dist: ruff<0.16,>=0.15; extra == 'dev'
Description-Content-Type: text/markdown

# afdb-query

Sequence-based programmatic access to the [AlphaFold Protein Structure Database](https://alphafold.ebi.ac.uk/) (AFDB). Query a protein by its amino-acid sequence, then pull per-residue pLDDT without hand-rolling URL derivation and JSON fetching.

## Install

```bash
pip install afdb-query
```

## Quickstart

```python
from afdb_query import AlphaFold

with AlphaFold() as af:
    hits = af.search(sequence)        # Tier 1: list[Structure], in AFDB's returned order
    s = hits[0]

    s.global_plddt        # mean pLDDT AFDB reports for the model (from the summary)
    s.sequence_identity   # 1.0 == exact match
    s.oligomeric_state    # "MONOMER", "HOMODIMER", ...
    s.uniprot_accession   # e.g. "P12345", or None

    p = s.plddt()         # Tier 2: per-residue pLDDT (fetched once, then cached)
    p.scores              # full per-residue list[float]
    p.residue_numbers     # parallel residue numbers — do not assume 1..N
```

`search` raises `InvalidSequenceError` for sequences that cannot be queried
(internal stop `*`, shorter than 20 residues, or non-standard amino acids), and
returns `[]` when AFDB has no entry for a valid sequence. Every HTTP or transport
failure is raised as `AFDBHTTPError`, so catching `AFDBError` is enough — you never
need to import `httpx`.

## Choosing a structure

A sequence query can match several structures, and **they are not interchangeable**:

- `hits[0]` is *not* guaranteed to be the canonical `AF-<accession>-F1` model. For
  some sequences a multi-chain or numeric AB-INITIO model ranks first.
- `global_plddt` (`confidence_avg_local_score`) is averaged over the **whole
  deposited model**. For a HOMODIMER or HETERODIMER that spans every chain, not just
  the one matching your query, so it is not comparable with a monomer's.

`select_group` applies three preference tiers — caller's accession, then monomers over
complexes, then canonical `-F1` over numeric ids — and returns **everything still tied**:

```python
from afdb_query import AlphaFold, select_group, mean_global_plddt, is_monomer

with AlphaFold() as af:
    doc = af.fetch_summary(sequence)           # None when AFDB has no entry
    group = select_group(doc["structures"]) if doc else []

    if not all(is_monomer(s) for s in group):
        ...                                   # your call: skip, or use it knowingly

    plddt = mean_global_plddt(group)           # average across the group, not one member
```

**No tier reads the confidence scores.** Ranking candidates by pLDDT and returning the
winner is not merely arbitrary, it is *biased*: the expected maximum of N draws rises
with N, so a comparison whose two arms match different numbers of candidates gets
different inflation on each side.

**There is no fourth tier, and that is deliberate.** After the three tiers a tie is
common — 25.5% of records on a real cache. Breaking it by any rule at all throws away
N−1 predictions of the same sequence for nothing. `select_group` keeps them and
`mean_global_plddt` averages across them: unbiased for the same reason an arbitrary
pick is (neither consults the value), with strictly lower variance because it uses
every prediction.

Tied candidates are near-always one protein reached through different UniProt entries —
identical-sequence orthologs, isoform duplicates, TrEMBL redundancy. On a real cache
4,288 of 4,289 tied sets were entirely full-length `-F1` models, so their per-residue
arrays are the same length and average elementwise too.

Do **not** take `group[0]`. The ordering is for reproducible output only; treating it
as "the" structure reintroduces exactly the arbitrary choice this API exists to avoid.

### Length is not visible to the tiers

`coverage == 1.0` means *the query is fully covered by the model*, **not** that the
model is the query's size — that is how an 860-residue multi-chain entry reports
coverage 1.0 against a 430-residue query. Filtering to monomers removes that case, but
an ortholog carrying the query sequence plus a few extra terminal residues passes every
visible tier, and then breaks two things: positional slicing (an offset from the
query's amino-acid length indexes the wrong residues) and the group average (its mean
spans residues your query does not contain).

The summary carries no residue count, so `filter_by_length` takes the lengths you
learned from fetching per-residue confidence:

```python
from afdb_query import filter_by_length

lengths = {s["model_identifier"]: len(fetched[s["model_identifier"]]) for s in group}
kept, dropped = filter_by_length(group, lengths, expected_length=len(sequence))
if dropped:
    log(f"{len(dropped)} of {len(group)} matched entries were the wrong length")
```

Unknown lengths are dropped, never assumed to conform, and dropped members are returned
rather than discarded so the loss is reportable.

`select_group` does not filter, either. "No usable structure" and "a group whose
average spans a complex" need different handling, and only the caller knows which its
analysis can tolerate — so test `is_monomer` on the members and decide.

## Batch lookups

`search_many` runs many sequences concurrently with resumable on-disk caching:

```python
report = af.search_many(
    [{"id": "rec1", "sequence": seq1}, {"id": "rec2", "sequence": seq2}],
    out_dir="afdb_cache",
    concurrency=6,
)
```

```python
{
  "total":    2,
  "skipped":  0,                                                   # already cached
  "filtered": {"internal_stop": 0, "too_short": 0, "nonstandard_aa": 0, "total": 0},
  "queried":  {"hits": 2, "misses": 0, "errors": 0, "total": 2},
}
```

`total == skipped + filtered["total"] + queried["total"]`; no count appears twice.

- You supply a generic `id` per sequence; it keys the cache file and maps back to
  your own records.
- `out_dir/summaries/{id}.json` stores each hit (a 404 miss stores
  `{"structures": []}`); existing files are left untouched, so re-runs resume.
- Real HTTP errors are counted but not saved, so they retry on the next run.

`search_many` fetches **summaries only**. It does not choose a structure and it does
not fetch per-residue confidence: which of several exact-sequence matches answers
your question is a property of your analysis, not of AFDB, and making that choice
inside a batch runner would hide it. Run `select_group` over the cached summaries and
pass its members to `fetch_plddt_many` when you need per-residue data.

## Per-residue pLDDT in bulk

`fetch_plddt_many` caches the **full** per-residue array for each structure you chose:

```python
from afdb_query import (
    AlphaFold, select_group, fetch_plddt_many, load_plddt, mean_per_residue,
)
import json, pathlib

cache = pathlib.Path("afdb_cache")

# summaries were cached earlier by search_many. One record per GROUP MEMBER: the id
# keys the cache file, so members of the same group must not collide on it.
records, members = [], {}
for f in (cache / "summaries").glob("*.json"):
    group = select_group(json.loads(f.read_text()).get("structures") or [])
    members[f.stem] = []
    for s in group:
        rid = f"{f.stem}::{s['model_identifier']}"
        members[f.stem].append(rid)
        records.append({"id": rid, "model_url": s["model_url"]})

with AlphaFold() as af:
    fetch_plddt_many(af, records, cache)

# Read one record's group back -- with the SAME composite ids that were written -- and
# average it. mean_per_residue raises if the members disagree in length.
plddts = [p for rid in members["rec1"] if (p := load_plddt(cache, rid))]
consensus = mean_per_residue([p.scores for p in plddts])

plddts[0].mean()             # one definition of "mean pLDDT"
plddts[0].mean(start=42)     # ...and of a region of it
```

Resumability keys on `afdb_cache/plddt/{id}.json`, **not** on the summary — so running
summaries first and residues later back-fills every already-cached record instead of
skipping it.

## Mean pLDDT, and why it takes bounds

`mean_plddt(scores, start=None, stop=None)` is the single definition. AFDB reports a
pre-computed `confidence_avg_local_score`; other predictors give you only the array. A
codebase that reads the field in one place and averages in another has two definitions
of its headline number and no way to tell them apart.

Nothing rounds — AFDB's field arrives rounded server-side, so `mean_plddt` over the same
array can differ in the last digit. Round when you format.

**Bounds exist because global means are not comparable between a protein and a
sub-range of itself.** Drop a disordered N-terminal tail and the mean rises purely
because low scorers left the set. `shared_suffix_means(long, short)` does the
length-controlled version:

```python
from afdb_query import shared_suffix_means

shared_suffix_means(canonical.scores, truncated.scores)
# {"offset": 71,
#  "shared_long": ...,   # mean over the residues both have — same count each side
#  "shared_short": ...,
#  "displaced": ...}     # mean over the residues only the longer one has
```

`residue_index` and `is_contiguous` are there because per-residue arrays are parallel to
`residueNumber`, and treating array position as residue number is an assumption, not a
guarantee.

### Averaging a group's arrays

`mean_per_residue(arrays, expected_length=None)` builds the consensus array for a group
`select_group` returned. It **raises on ragged input** rather than `zip`-truncating:
unequal lengths make the index set non-rectangular, so position `i` stops denoting the
same residue in every member and every downstream region mean silently describes a
comparison that does not exist.

Pass `expected_length=len(sequence)` to additionally require the arrays be the query's
size, and use `filter_by_length` to drop non-conforming members first.

The commutation this rests on — `mean over region of (elementwise mean)` equals
`mean over group of (region mean)` — holds because both are unweighted means over the
same rectangular index set, and is asserted in the test suite.

## Not (yet) supported

- UniProt-accession lookup (sequence-only for now)
- PAE (Predicted Aligned Error)

## Development

```bash
pip install -e ".[dev]"
ruff check . && ruff format --check .
pytest                    # unit suite; 100% branch coverage is enforced
pytest -m integration     # live AFDB tests (network required)
```

The integration suite includes a ground-truth check that the per-residue pLDDT in
AFDB's confidence JSON matches the B-factor column of the deposited mmCIF — two
independently generated files — so upstream API drift surfaces here rather than in
your results.
