Add the implementation plan for building project taxa lists from a region

docs/claude/analysis/2026-07-02-phase0-regional-coverage-findings.md

@@ -0,0 +1,44 @@
+# Phase 0 findings — regional taxa-list coverage spike (#1364)
+
+**Date:** 2026-07-02
+**Verdict:** GO for Proposal A. The candidate GBIF and iNaturalist endpoints work with the parameters the plan assumed, the A3 reverse-geocode path resolves cleanly, and a region-derived list covers a useful majority of a real classifier's labels.
+**Script:** `docs/claude/analysis/phase0_regional_coverage_spike.py` (run inside the django container; reads DB + live APIs, writes nothing).
+
+## Setup
+
+- Region: Vermont (chosen to match the Quebec & Vermont classifier).
+- Taxon scope: Lepidoptera (GBIF taxonKey 797, iNat taxon_id 47157).
+- Classifier under test: local `Algorithm` pk 10, "Quebec & Vermont Species Classifier - Apr 2024", **2497 labels**.
+- Name match: exact string on the scientific name (same join class masking uses — `Taxon.name` == label).
+
+## Results (measured, this run)
+
+| Metric | Value |
+|---|---|
+| Q&V classifier labels | 2497 |
+| GBIF VT Lepidoptera species | 2164 |
+| iNat VT Lepidoptera species | 1922 |
+| Region union (GBIF ∪ iNat) | 2299 |
+| Q&V ∩ GBIF | 1743 (69.8% of labels) |
+| Q&V ∩ iNat | 1527 (61.2% of labels) |
+| **Q&V ∩ union** (default masking-list size) | **1749 (70.0% of labels)** |
+| Region species not in Q&V (no-coverage bucket) | 550 (23.9% of the union) |
+
+Reverse-geocode check (A3): point `(44.26, -72.58)` → GADM level-1 gid `USA.46_1` ("Vermont"). The `/geocode/reverse` response carries GADM ids at every level, so deriving a region code from a deployment's stored `latitude`/`longitude` is straightforward — the `--all-projects` backfill path is viable.
+
+## Interpretation
+
+- **A region-derived default list is worth building.** For Vermont it keeps 1749 of 2497 classes and masks ~748 (30%) that neither GBIF nor iNat records in the state — a meaningful, defensible reduction, which is exactly what class masking is for.
+- **GBIF is the primary source; iNat is secondary.** GBIF alone reaches 69.8%; adding iNat lifts the *intersection* by only 6 species (to 70.0%). iNat's larger contribution is to the regional *union* (the uncovered bucket), not to classifier coverage. This supports starting Phase 1 with a **GBIF-only** source client and adding iNat later — the wide-union merge still holds, it just starts with one source.
+- **The `include_uncovered` opt-in is a real, non-trivial set** (550 species for Vermont). The default-subset-to-model-known behaviour matters: without it, a naive regional list would add 550 taxa no model can predict.
+
+## Caveats / what still needs verifying (carried into Phase 1)
+
+- **Name-join fragility is the top risk.** 748 Q&V labels (30%) are absent from the region union. Some are genuine regional absences; some are likely name-format or synonym mismatches (authorship, subspecies, GBIF/iNat vs. our backbone naming). Before trusting the default list, audit a sample of those 748 to estimate the true-absence vs. mismatch split, and decide whether to widen matching (e.g. `Taxon.search_names`, synonym resolution) — this is the single measurement that most affects list quality.
+- **GBIF name resolution is slow** (2164 speciesKey→name lookups took ~92s at 16 threads). Production must cache resolved keys and/or resolve lazily; not a blocker, but the service needs a cache layer.
+- **Granularity was GADM level 1 (state).** Level 2 (county) would tighten the list but risk under-including; the region-code field should carry the level, and the right default is an open question.
+- **Single region / single classifier tested.** Vermont + Q&V is the friendliest case (the classifier was built for this region). Coverage for a mismatched region (e.g. a tropical site against Q&V) will be far lower — which is the correct behaviour, but worth confirming the numbers degrade sensibly before enabling auto-masking broadly.
+
+## Recommendation for Phase 1
+
+Proceed. Build the core service with a **GBIF-first** source client (iNat behind the same protocol, added second), the wide-union merge, the model-coverage relationship, and the default-subset behaviour — all behind unit tests with a stubbed source (no network in CI). Fold a name-mismatch audit of the uncovered Q&V labels into the Phase 1 work, since it gates how much to invest in fuzzy matching.

docs/claude/analysis/phase0_regional_coverage_spike.py

@@ -0,0 +1,198 @@
+"""Phase 0 spike (#1364) — verify GBIF/iNat regional endpoints and measure real
+classifier label coverage. Throwaway measurement script, run inside the django
+container so it has DB + network + the repo's create_session().
+
+    docker compose exec -T django python docs/claude/analysis/phase0_regional_coverage_spike.py
+
+Region: Vermont (matches the Quebec & Vermont classifier). Taxon scope: Lepidoptera.
+Nothing is written to the DB. Output is a coverage report on stdout.
+"""
+
+from __future__ import annotations
+
+import concurrent.futures
+import sys
+import time
+
+from ami.ml.models.algorithm import Algorithm
+from ami.utils.requests import create_session
+
+GBIF = "https://api.gbif.org/v1"
+INAT = "https://api.inaturalist.org/v1"
+LEP_GBIF = 797  # Lepidoptera, GBIF backbone taxonKey
+LEP_INAT = 47157  # Lepidoptera, iNaturalist taxon_id
+VT_POINT = (44.26, -72.58)  # Montpelier, VT — for reverse-geocode (also exercises A3)
+QV_ALGORITHM_PK = 10  # "Quebec & Vermont Species Classifier - Apr 2024" (2497 labels)
+
+
+def norm(name: str) -> str:
+    return (name or "").strip()
+
+
+def load_classifier_labels() -> set[str]:
+    algo = Algorithm.objects.get(pk=QV_ALGORITHM_PK)
+    labels = algo.category_map.labels or []
+    print(f"[classifier] {algo.name!r} — {len(labels)} labels")
+    return {norm(x) for x in labels if norm(x)}
+
+
+def gbif_reverse_gadm1(session, lat, lng) -> str | None:
+    """A3 check: reverse-geocode a point to its GADM level-1 gid."""
+    r = session.get(f"{GBIF}/geocode/reverse", params={"lat": lat, "lng": lng}, timeout=30)
+    r.raise_for_status()
+    gadm1 = None
+    for item in r.json():
+        gid = item.get("id", "")
+        # GADM level-1 ids look like "USA.46_1"; level-0 "USA", level-2 "USA.46.14_1"
+        if item.get("source", "").lower().find("gadm") >= 0 or gid.startswith("USA"):
+            parts = gid.split(".")
+            if len(parts) == 2 and gid.endswith("_1"):
+                gadm1 = gid
+                print(f"[gbif] reverse-geocode {lat},{lng} -> GADM1 {gid} ({item.get('title')})")
+    return gadm1
+
+
+def gbif_vt_species(session, gadm_gid: str, cap: int = 3000) -> set[str]:
+    """Distinct Lepidoptera species in the region via speciesKey facet, resolved to names."""
+    keys: list[int] = []
+    offset = 0
+    page = 1000
+    while True:
+        r = session.get(
+            f"{GBIF}/occurrence/search",
+            params={
+                "taxonKey": LEP_GBIF,
+                "gadmGid": gadm_gid,
+                "hasCoordinate": "true",
+                "facet": "speciesKey",
+                "facetLimit": page,
+                "facetOffset": offset,
+                "limit": 0,
+            },
+            timeout=60,
+        )
+        r.raise_for_status()
+        facets = r.json().get("facets", [])
+        counts = facets[0]["counts"] if facets else []
+        if not counts:
+            break
+        keys.extend(int(c["name"]) for c in counts)
+        offset += page
+        if len(counts) < page or len(keys) >= cap:
+            break
+    keys = keys[:cap]
+    print(f"[gbif] {len(keys)} distinct Lepidoptera speciesKeys in {gadm_gid} (cap={cap}); resolving names...")
+
+    def resolve(key: int) -> str | None:
+        s = create_session()
+        try:
+            rr = s.get(f"{GBIF}/species/{key}", timeout=30)
+            if rr.status_code != 200:
+                return None
+            d = rr.json()
+            return norm(d.get("canonicalName") or d.get("species") or d.get("scientificName"))
+        except Exception:
+            return None
+
+    names: set[str] = set()
+    t0 = time.time()
+    with concurrent.futures.ThreadPoolExecutor(max_workers=16) as ex:
+        for nm in ex.map(resolve, keys):
+            if nm:
+                names.add(nm)
+    print(f"[gbif] resolved {len(names)} species names in {time.time() - t0:.0f}s")
+    return names
+
+
+def inat_place_id(session, q: str) -> int | None:
+    r = session.get(f"{INAT}/places/autocomplete", params={"q": q, "per_page": 10}, timeout=30)
+    r.raise_for_status()
+    results = r.json().get("results", [])
+    # admin_level 10 == state/province in iNat
+    for p in results:
+        if p.get("admin_level") == 10 and norm(p.get("name")) == q:
+            print(f"[inat] place {q!r} -> place_id {p['id']} (admin_level 10)")
+            return p["id"]
+    if results:
+        p = results[0]
+        print(f"[inat] place {q!r} -> place_id {p['id']} ({p.get('name')}, admin_level {p.get('admin_level')}) [fallback]")
+        return p["id"]
+    return None
+
+
+def inat_vt_species(session, place_id: int) -> set[str]:
+    """Species-level Lepidoptera in the place; names come back inline (fast)."""
+    names: set[str] = set()
+    page = 1
+    while True:
+        r = session.get(
+            f"{INAT}/observations/species_counts",
+            params={
+                "place_id": place_id,
+                "taxon_id": LEP_INAT,
+                "quality_grade": "research",
+                "hrank": "species",
+                "per_page": 500,
+                "page": page,
+            },
+            timeout=60,
+        )
+        r.raise_for_status()
+        data = r.json()
+        results = data.get("results", [])
+        for row in results:
+            t = row.get("taxon", {})
+            if t.get("rank") == "species" and t.get("name"):
+                names.add(norm(t["name"]))
+        total = data.get("total_results", 0)
+        if page * 500 >= total or not results:
+            break
+        page += 1
+    print(f"[inat] {len(names)} species-rank Lepidoptera in place {place_id}")
+    return names
+
+
+def pct(a: int, b: int) -> str:
+    return f"{(100.0 * a / b):.1f}%" if b else "n/a"
+
+
+def main() -> int:
+    session = create_session()
+    labels = load_classifier_labels()
+
+    gadm1 = gbif_reverse_gadm1(session, *VT_POINT)
+    gbif_names: set[str] = set()
+    if gadm1:
+        gbif_names = gbif_vt_species(session, gadm1)
+    else:
+        print("[gbif] could not resolve GADM1 gid; skipping GBIF")
+
+    place = inat_place_id(session, "Vermont")
+    inat_names = inat_vt_species(session, place) if place else set()
+
+    union = gbif_names | inat_names
+    print("\n==================== COVERAGE REPORT ====================")
+    print(f"Q&V classifier labels        : {len(labels)}")
+    print(f"GBIF VT Lepidoptera species  : {len(gbif_names)}")
+    print(f"iNat VT Lepidoptera species  : {len(inat_names)}")
+    print(f"Region union (GBIF ∪ iNat)   : {len(union)}")
+    print("--------------------------------------------------------")
+    cov_g = labels & gbif_names
+    cov_i = labels & inat_names
+    cov_u = labels & union
+    print(f"Q&V ∩ GBIF   : {len(cov_g):5d}  ({pct(len(cov_g), len(labels))} of labels)")
+    print(f"Q&V ∩ iNat   : {len(cov_i):5d}  ({pct(len(cov_i), len(labels))} of labels)")
+    print(f"Q&V ∩ union  : {len(cov_u):5d}  ({pct(len(cov_u), len(labels))} of labels)  <-- default masking-list size")
+    print("--------------------------------------------------------")
+    only_region = union - labels
+    print(f"Region species NOT in Q&V labels (no-model-coverage bucket): {len(only_region)}")
+    print(f"  ({pct(len(only_region), len(union))} of the regional union — these are the include_uncovered opt-in rows)")
+    print("========================================================")
+    # A few examples of the intersection and the uncovered set for eyeballing name-match quality.
+    print("sample Q&V∩union :", sorted(list(cov_u))[:8])
+    print("sample region-only:", sorted(list(only_region))[:8])
+    return 0
+
+
+if __name__ == "__main__":
+    sys.exit(main())