Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
40 changes: 29 additions & 11 deletions protea/api/routers/proteins_stats.py
Original file line number Diff line number Diff line change
Expand Up @@ -85,15 +85,19 @@ def _compute_embeddings_by_plm(factory: sessionmaker[Session]) -> dict[str, Any]

GROUPs ``sequence_embedding`` by ``embedding_config_id``, JOINs back to
``protein`` via ``sequence_id`` so the count reflects PROTEIN coverage
(sequence-dedup is folded back out). Total uses the canonical Swiss-Prot
population (``is_canonical=True AND reviewed=True``) so the percent
matches what the surrounding cards report.
(sequence-dedup is folded back out).

THE POPULATION IS THE CANONICAL CORPUS, not the Swiss-Prot spine. Both
sides carried ``reviewed=True``, making the arithmetic honest about the
wrong thing: 576.450 at 100,0% while 746.014 canonical proteins exist and
all of them are embedded. ``reviewed`` had been added to stop a measured
121,5%, whose cause was an unaligned numerator, so aligning both sides on
``is_canonical`` alone fixes it without hiding 169.564 proteins.
"""
with session_scope(factory) as session:
total = (
session.query(func.count(Protein.accession))
.filter(Protein.is_canonical.is_(True))
.filter(Protein.reviewed.is_(True))
.scalar()
or 0
)
Expand All @@ -108,12 +112,11 @@ def _compute_embeddings_by_plm(factory: sessionmaker[Session]) -> dict[str, Any]
)
.join(SequenceEmbedding, SequenceEmbedding.embedding_config_id == EmbeddingConfig.id)
.join(Protein, Protein.sequence_id == SequenceEmbedding.sequence_id)
# The SAME two filters as the total below. Without them the numerator
# counts isoforms, TrEMBL and merged accessions that the denominator
# excludes, and the percentage goes over 100: measured at 121,5% with
# 700.475 over 576.450 while an embedding pass was running.
# The SAME filter as the total above. Without it the numerator counts
# isoforms that the denominator excludes and the percentage goes over
# 100: measured at 121,5%, 700.475 over 576.450, while an embedding
# pass was running.
.filter(Protein.is_canonical.is_(True))
.filter(Protein.reviewed.is_(True))
.group_by(
EmbeddingConfig.id,
EmbeddingConfig.model_name,
Expand Down Expand Up @@ -204,17 +207,32 @@ def _compute_taxonomy_top(factory: sessionmaker[Session]) -> dict[str, Any]:
Only counts canonical rows so isoforms don't double-weight an organism.
Returns ``taxonomy_id`` alongside the human-readable ``organism`` so the
frontend can deep-link to a UniProt taxonomy page.

GROUPED BY TAXON, NOT BY NAME. The two load paths spell an organism
differently -- the Swiss-Prot path writes ``Homo sapiens`` and the TSV path
``Homo sapiens (Human)`` -- so grouping by the string split one organism
into two bars and ranked the halves. 726 taxa carry more than one name, and
every one of the top six did, so the ordering was wrong from the top down:
taxon 9606 is 41.823 canonical proteins and showed as 21.387 and 20.436.

The taxon is the identity and the name is a label, so the label is picked
with ``min`` to keep it deterministic rather than whichever row the plan
happened to emit. Rows with no taxon are excluded for the same reason they
always were excluded by the ``organism`` filter: the 28.868 deleted
UniProtKB entries carry neither, and collapsing them into one unnamed bar
would invent an organism.
"""
with session_scope(factory) as session:
rows = (
session.query(
Protein.organism,
func.min(Protein.organism),
Protein.taxonomy_id,
func.count(Protein.accession),
)
.filter(Protein.is_canonical.is_(True))
.filter(Protein.organism.is_not(None))
.group_by(Protein.organism, Protein.taxonomy_id)
.filter(Protein.taxonomy_id.is_not(None))
.group_by(Protein.taxonomy_id)
.order_by(func.count(Protein.accession).desc())
.limit(20)
.all()
Expand Down
96 changes: 86 additions & 10 deletions tests/test_proteins_stats_router.py
Original file line number Diff line number Diff line change
Expand Up @@ -510,20 +510,96 @@ def test_the_high_throughput_codes_are_in(self):


class TestCoverageDividesTwoCountsOfOneThing:
def test_the_numerator_carries_the_same_filters_as_the_total(self):
"""A percentage over 100 was the symptom, not the defect.
"""A percentage over 100 was the symptom; two populations was the defect.

Measured while an embedding pass was running: 700.475 proteins with a
vector over a denominator of 576.450 canonical reviewed ones, 121,5%. The
first fix aligned both sides on ``is_canonical AND reviewed``, which made
the arithmetic honest about the wrong thing: it then reported 576.450 at
100,0% while 746.014 canonical proteins exist and every one of them is
embedded. Aligning on ``is_canonical`` alone keeps the percentage sound
and stops the card hiding 169.564 proteins.
"""

@staticmethod
def _filtros(fuente: str) -> set[str]:
"""The ``Protein.<col>`` each ``.filter(...)`` names, as a set."""
import re

return set(re.findall(r"\.filter\((Protein\.\w+)\.", fuente))

Measured while an embedding pass was running: 700.475 proteins with a
vector over a denominator of 576.450 canonical reviewed ones, 121,5%.
The numerator counted isoforms, TrEMBL entries and merged accessions
that the denominator excludes.
def test_both_sides_count_the_same_population(self):
"""The invariant, stated once instead of by counting a literal twice.

The previous version asserted that a particular filter appeared twice
in the source, so it had to be edited whenever the population changed
and said nothing about the two sides agreeing. What matters is that
they agree.
"""
import inspect

from protea.api.routers import proteins_stats as mod

src = inspect.getsource(mod._compute_embeddings_by_plm)
# Two filter CALLS, not two mentions: the docstring names the population
# too, and counting words would pass on a comment alone.
assert src.count(".filter(Protein.is_canonical.is_(True))") == 2
assert src.count(".filter(Protein.reviewed.is_(True))") == 2
total, numerador = src.split("rows = (", 1)
assert self._filtros(total) == self._filtros(numerador), (
"numerator and denominator filter on different populations"
)

def test_the_population_is_the_canonical_corpus(self):
"""Canonical, and NOT narrowed to the Swiss-Prot spine.

``is_canonical`` because an isoform shares its parent's sequence and
would count as separate coverage of the same thing.
"""
import inspect

from protea.api.routers import proteins_stats as mod

src = inspect.getsource(mod._compute_embeddings_by_plm)
total, numerador = src.split("rows = (", 1)
assert self._filtros(total) == {"Protein.is_canonical"}
assert "Protein.reviewed" not in self._filtros(numerador)


class TestOneOrganismIsOneBar:
"""The ranking grouped by NAME, and the two load paths disagree on names.

The Swiss-Prot path writes ``Homo sapiens`` and the TSV path
``Homo sapiens (Human)``, so taxon 9606 appeared as two bars and the
ranking compared halves against wholes. Measured on the live corpus: 726
taxa carry more than one name, and every one of the top six did, so the
order was wrong from the top down -- 41.823 canonical proteins showed as
21.387 and 20.436.
"""

@staticmethod
def _sql() -> str:
import inspect

from protea.api.routers import proteins_stats as mod

return inspect.getsource(mod._compute_taxonomy_top)

def test_it_groups_by_taxon_and_not_by_name(self):
src = self._sql()
assert ".group_by(Protein.taxonomy_id)" in src
assert "Protein.organism, Protein.taxonomy_id" not in src

def test_the_name_is_aggregated_so_two_spellings_cannot_split_a_taxon(self):
"""Selecting the bare column would make the group key include it again."""
src = self._sql()
assert "func.min(Protein.organism)" in src

def test_rows_with_no_taxon_are_excluded(self):
"""Otherwise the 28.868 deleted entries collapse into one unnamed bar.

They carry neither organism nor taxonomy_id, so the organism filter
already removed them; this states it of the new group key too.
"""
src = self._sql()
assert ".filter(Protein.taxonomy_id.is_not(None))" in src

def test_isoforms_still_do_not_double_weight(self):
src = self._sql()
assert ".filter(Protein.is_canonical.is_(True))" in src
Loading