make QA scoring reproducible and threshold-independent

The greedy IoU matcher gave a reference to whichever candidate was offered
first. Row order decided that, and every detection in a run shares one
transaction timestamp, so ordering by created_at left the assignment
undefined: the same QA run over the same data produced different mean IoU,
and the geometry shown to a reviewer as a false positive could be the better
of two detections. Candidates are now ranked by confidence with feature
identity as tiebreaker, which is also the COCO/PASCAL rule.

A single precision/recall/F1 triple describes one operating point, so two
models cannot be compared from it: a conservatively calibrated model looks
worse at a low confidence cut and better at a high one without detecting
anything differently. DetectionMetricsService adds the full curve, average
precision and the threshold where F1 actually peaks.

Also:
- report the population the metrics were computed over, so
  matches + false_positives equals candidate_feature_count even under an
  area filter; raw dataset totals move to the _raw fields;
- state whether candidates are axis-aligned boxes or footprint polygons.
  A box can never reach IoU 1 against a rotated building, so the strict
  score has a ceiling that has nothing to do with detection quality.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

backend/app/services/detection_metrics_service.py

@@ -0,0 +1,175 @@
+"""Threshold-independent detection metrics.
+
+A single precision/recall/F1 triple describes one operating point. Which point
+that is depends on the confidence threshold the operator typed, so two models
+cannot be compared from it: a conservatively calibrated model looks worse at a
+low cut and better at a high one without detecting anything differently.
+
+This service produces the standard alternative — the full precision/recall
+curve over every confidence value that occurs in the run, the average
+precision derived from it, and the threshold where F1 actually peaks — using
+the same greedy IoU matching rule as the rest of QA so the numbers stay
+comparable with the persisted quality checks.
+"""
+
+from __future__ import annotations
+
+from typing import Any
+
+from shapely.geometry.base import BaseGeometry
+from shapely.strtree import STRtree
+
+from app.services.qa_service import QaService
+
+
+class DetectionMetricsService:
+    SUPPORTED_GEOMETRY_TYPES = QaService.SUPPORTED_GEOMETRY_TYPES
+
+    @staticmethod
+    def _rank_candidates(
+        candidates: list[tuple[dict[str, Any], BaseGeometry]],
+    ) -> list[tuple[float, str, dict[str, Any], BaseGeometry]]:
+        """Order candidates by confidence, highest first, deterministically."""
+
+        ranked: list[tuple[float, str, dict[str, Any], BaseGeometry]] = []
+        for index, (feature, geometry) in enumerate(candidates):
+            if geometry.is_empty or geometry.geom_type not in DetectionMetricsService.SUPPORTED_GEOMETRY_TYPES:
+                continue
+            if geometry.area <= 0:
+                continue
+            confidence = QaService._feature_confidence(feature)
+            identifier = QaService._feature_identifier(feature, "candidate", index)
+            ranked.append((float(confidence if confidence is not None else 0.0), identifier, feature, geometry))
+        ranked.sort(key=lambda item: (-item[0], item[1]))
+        return ranked
+
+    @staticmethod
+    def _greedy_hits(
+        ranked: list[tuple[float, str, dict[str, Any], BaseGeometry]],
+        references: list[tuple[dict[str, Any], BaseGeometry]],
+        iou_threshold: float,
+    ) -> tuple[list[bool], int]:
+        """Mark each ranked candidate as a hit or a miss, best score first.
+
+        Walking the ranking once and consuming references as they are claimed
+        is exactly the COCO/PASCAL rule, and it is what makes the result
+        independent of the order rows came out of the database.
+        """
+
+        supported = [
+            (index, geometry)
+            for index, (_, geometry) in enumerate(references)
+            if geometry.geom_type in DetectionMetricsService.SUPPORTED_GEOMETRY_TYPES
+            and not geometry.is_empty
+            and geometry.area > 0
+        ]
+        reference_count = len(supported)
+        if not supported:
+            return [False] * len(ranked), 0
+
+        tree = STRtree([geometry for _, geometry in supported])
+        claimed: set[int] = set()
+        hits: list[bool] = []
+
+        for _, _, _, geometry in ranked:
+            best_iou = 0.0
+            best_index: int | None = None
+            for position in sorted(int(value) for value in tree.query(geometry)):
+                if position in claimed:
+                    continue
+                _, reference_geometry = supported[position]
+                intersection_area = geometry.intersection(reference_geometry).area
+                if intersection_area <= 0:
+                    continue
+                union_area = geometry.area + reference_geometry.area - intersection_area
+                if union_area <= 0:
+                    continue
+                iou = intersection_area / union_area
+                if iou > best_iou:
+                    best_iou = iou
+                    best_index = position
+            if best_index is not None and best_iou >= iou_threshold:
+                claimed.add(best_index)
+                hits.append(True)
+            else:
+                hits.append(False)
+
+        return hits, reference_count
+
+    @staticmethod
+    def _average_precision(points: list[dict[str, Any]]) -> float:
+        """Area under the precision/recall curve, with precision made monotone.
+
+        Interpolating precision to its running maximum from the right is the
+        VOC/COCO convention; without it the sawtooth from individual false
+        positives shows up as noise in the score.
+        """
+
+        if not points:
+            return 0.0
+
+        recalls = [0.0] + [point["recall"] for point in points]
+        precisions = [1.0] + [point["precision"] for point in points]
+        for index in range(len(precisions) - 2, -1, -1):
+            precisions[index] = max(precisions[index], precisions[index + 1])
+
+        area = 0.0
+        for index in range(1, len(recalls)):
+            area += (recalls[index] - recalls[index - 1]) * precisions[index]
+        return area
+
+    @staticmethod
+    def precision_recall_curve(
+        candidates: list[tuple[dict[str, Any], BaseGeometry]],
+        references: list[tuple[dict[str, Any], BaseGeometry]],
+        *,
+        iou_threshold: float,
+    ) -> dict[str, Any]:
+        """Sweep every confidence value present and describe the whole curve."""
+
+        ranked = DetectionMetricsService._rank_candidates(candidates)
+        hits, reference_count = DetectionMetricsService._greedy_hits(ranked, references, iou_threshold)
+
+        points: list[dict[str, Any]] = []
+        true_positives = 0
+        for position, hit in enumerate(hits):
+            if hit:
+                true_positives += 1
+            false_positives = position + 1 - true_positives
+            precision = true_positives / (position + 1)
+            recall = true_positives / reference_count if reference_count else 0.0
+            f1 = (2 * precision * recall / (precision + recall)) if precision + recall > 0 else 0.0
+            points.append(
+                {
+                    "confidence_threshold": ranked[position][0],
+                    "candidate_count": position + 1,
+                    "true_positives": true_positives,
+                    "false_positives": false_positives,
+                    "false_negatives": max(0, reference_count - true_positives),
+                    "precision": precision,
+                    "recall": recall,
+                    "f1_score": f1,
+                }
+            )
+
+        # Keep one point per distinct threshold: the last one, which is the
+        # complete tally for everything at or above that confidence.
+        deduplicated: list[dict[str, Any]] = []
+        for point in points:
+            if deduplicated and deduplicated[-1]["confidence_threshold"] == point["confidence_threshold"]:
+                deduplicated[-1] = point
+            else:
+                deduplicated.append(point)
+
+        best = max(deduplicated, key=lambda point: (point["f1_score"], point["confidence_threshold"]), default=None)
+        return {
+            "iou_threshold": iou_threshold,
+            "reference_count": reference_count,
+            "candidate_count": len(ranked),
+            "average_precision": DetectionMetricsService._average_precision(points),
+            "best_f1": best["f1_score"] if best else 0.0,
+            "best_f1_threshold": best["confidence_threshold"] if best else None,
+            "best_f1_precision": best["precision"] if best else None,
+            "best_f1_recall": best["recall"] if best else None,
+            "points": deduplicated,
+        }

backend/app/services/detection_qa_service.py

@@ -168,24 +168,61 @@ class DetectionQaService:
             clipped_boundary_count=clipped_boundary_count,
         )
 
+    # A polygon whose area is within this fraction of its own bounding box is
+    # an axis-aligned rectangle for practical purposes.
+    RECTANGULAR_AREA_RATIO = 0.99
+
+    @staticmethod
+    def candidate_geometry_mode(geometries: list[tuple[dict[str, Any], BaseGeometry]]) -> str:
+        """Say whether the candidates are detector boxes or true footprints.
+
+        It matters for reading the score. An axis-aligned box can never reach
+        IoU 1 against a rotated or L-shaped building footprint, so a strict
+        footprint IoU understates a box detector by a fixed amount that has
+        nothing to do with whether it found the building.
+        """
+
+        polygonal = [
+            geometry
+            for _, geometry in geometries
+            if geometry.geom_type in {"Polygon", "MultiPolygon"} and geometry.area > 0
+        ]
+        if not polygonal:
+            return "unknown"
+        rectangular = sum(
+            1
+            for geometry in polygonal
+            if geometry.area / geometry.envelope.area >= DetectionQaService.RECTANGULAR_AREA_RATIO
+        )
+        return "axis_aligned_boxes" if rectangular == len(polygonal) else "footprint_polygons"
+
     @staticmethod
     def box_to_footprint_diagnostics(
         strict_evidence: QaMatchEvidence,
         envelope_evidence: QaMatchEvidence,
         *,
         iou_threshold: float,
+        candidate_geometry_mode: str = "unknown",
     ) -> dict[str, Any]:
         envelope_metrics = DetectionQaService._metrics(envelope_evidence)
-        return {
+        diagnostics = {
             "diagnostic_only": True,
             "canonical_method": "candidate_polygon_vs_reference_footprint_iou",
             "diagnostic_method": "candidate_polygon_vs_reference_envelope_iou",
             "iou_threshold": iou_threshold,
+            "candidate_geometry_mode": candidate_geometry_mode,
             "strict_matches": strict_evidence.matches,
             "envelope_matches": envelope_evidence.matches,
             "possible_box_to_footprint_mismatch_count": max(0, envelope_evidence.matches - strict_evidence.matches),
             **envelope_metrics,
         }
+        if candidate_geometry_mode == "axis_aligned_boxes":
+            diagnostics["interpretation"] = (
+                "Candidates are axis-aligned detector boxes. The strict footprint IoU therefore has a "
+                "ceiling below 1 for rotated or non-rectangular buildings; the envelope figures isolate "
+                "detection quality from that shape mismatch."
+            )
+        return diagnostics
 
     @staticmethod
     def _metrics(evidence: QaMatchEvidence) -> dict[str, Any]:

backend/app/services/qa_service.py

@@ -123,6 +123,51 @@ class QaService:
                 return str(value)
         return f"{fallback_prefix}-{index + 1}"
 
+    @staticmethod
+    def _feature_confidence(feature: dict[str, Any]) -> float | None:
+        """Read a detector confidence from a QA feature, if the source has one.
+
+        Vector-vs-vector comparisons have no confidence at all; those fall back
+        to identity ordering so the result stays reproducible either way.
+        """
+
+        candidates: list[Any] = [feature.get("confidence")]
+        properties = feature.get("properties")
+        if isinstance(properties, dict):
+            candidates.append(properties.get("confidence"))
+        for value in candidates:
+            if value is None or isinstance(value, bool):
+                continue
+            try:
+                confidence = float(value)
+            except (TypeError, ValueError):
+                continue
+            if confidence == confidence:  # reject NaN
+                return confidence
+        return None
+
+    @staticmethod
+    def _candidate_match_order(
+        source_supported: list[tuple[int, dict[str, Any], BaseGeometry]],
+    ) -> list[tuple[int, dict[str, Any], BaseGeometry]]:
+        """Order candidates the way detection benchmarks do: best score first.
+
+        Greedy IoU matching gives the reference to whichever candidate is
+        offered first, so the input order decides both the score and which
+        geometry an operator sees as false-positive evidence. Database row
+        order is not a defensible answer to that question — every detection in
+        a run shares one transaction timestamp — so candidates are ranked by
+        confidence, with feature identity as a stable tiebreaker.
+        """
+
+        def order_key(entry: tuple[int, dict[str, Any], BaseGeometry]) -> tuple[float, str, int]:
+            index, feature, _ = entry
+            confidence = QaService._feature_confidence(feature)
+            identifier = QaService._feature_identifier(feature, "candidate", index)
+            return (-(confidence if confidence is not None else 0.0), identifier, index)
+
+        return sorted(source_supported, key=order_key)
+
     @staticmethod
     def _match_io_u_evidence(
         source_geometries: list[tuple[dict[str, Any], BaseGeometry]],
@@ -153,7 +198,7 @@ class QaService:
                 unsupported=True,
                 false_positive_evidence=[
                     {"candidate_feature_id": QaService._feature_identifier(feature, "candidate", source_index)}
-                    for source_index, feature, _ in source_supported
+                    for source_index, feature, _ in QaService._candidate_match_order(source_supported)
                 ],
                 false_negative_evidence=[
                     {"reference_feature_id": QaService._feature_identifier(feature, "reference", reference_index)}
@@ -167,7 +212,7 @@ class QaService:
         }
         evidence = QaMatchEvidence(warnings=[f"Unsupported geometry types: {unsupported}"] if unsupported else [], unsupported=bool(unsupported))
 
-        for source_index, source_feature, source_geom in source_supported:
+        for source_index, source_feature, source_geom in QaService._candidate_match_order(source_supported):
             source_feature_id = QaService._feature_identifier(source_feature, "candidate", source_index)
             if source_geom.area <= 0:
                 evidence.false_positives += 1
@@ -278,6 +323,9 @@ class QaService:
             dataset_ids=(candidate_dataset_id, reference_dataset_id),
         )
 
+        candidate_feature_count_raw = len(candidate_payload.get("features", [])) if isinstance(candidate_payload, dict) else 0
+        reference_feature_count_raw = len(reference_payload.get("features", [])) if isinstance(reference_payload, dict) else 0
+
         if area_geometry is not None:
             candidate_geometries = QaService._apply_area_filter(candidate_geometries, area_geometry, dataset_id=candidate_dataset.id)
             reference_geometries = QaService._apply_area_filter(reference_geometries, area_geometry, dataset_id=reference_dataset.id)
@@ -288,8 +336,11 @@ class QaService:
             iou_threshold,
         )
 
-        candidate_feature_count = len(candidate_payload.get("features", [])) if isinstance(candidate_payload, dict) else 0
-        reference_feature_count = len(reference_payload.get("features", [])) if isinstance(reference_payload, dict) else 0
+        # Report the population the metrics were computed over, not the raw
+        # dataset totals: with an area filter the two differ, and a count that
+        # disagrees with matches + false positives is unreadable as evidence.
+        candidate_feature_count = len(candidate_geometries)
+        reference_feature_count = len(reference_geometries)
         mean_iou = None if not evidence.match_iou_values else sum(evidence.match_iou_values) / len(evidence.match_iou_values)
 
         precision = None
@@ -310,6 +361,8 @@ class QaService:
             warnings=_extract_crs_warnings(candidate_dataset, reference_dataset) + evidence.warnings,
             candidate_feature_count=candidate_feature_count,
             reference_feature_count=reference_feature_count,
+            candidate_feature_count_raw=candidate_feature_count_raw,
+            reference_feature_count_raw=reference_feature_count_raw,
             matches=evidence.matches,
             false_positives=evidence.false_positives,
             false_negatives=evidence.false_negatives,