From bbac84cf29a4ede4d32b5c3b024f11a646079b73 Mon Sep 17 00:00:00 2001 From: koki3070 Date: Wed, 17 Jun 2026 11:02:03 +0900 Subject: [PATCH 1/2] =?UTF-8?q?fix:=20issue=20#59=20=E2=80=94=20ad-hoc=20?= =?UTF-8?q?=E5=BE=AE=E5=B0=8F=E6=AD=A3=E5=AE=9A=E6=95=B0=E3=81=AE=E6=95=B4?= =?UTF-8?q?=E7=90=86=E3=81=A8=E8=AB=96=E6=96=87=E5=A4=96=E3=83=A2=E3=83=87?= =?UTF-8?q?=E3=83=AA=E3=83=B3=E3=82=B0=E5=BA=8A=E3=81=AE=E5=89=8A=E9=99=A4?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit 数値安定化定数を numerical_eps.py に集約し、encoder clamp(0.2) と legacy +1e-4 を削除。 sqrt(0) 勾配爆発を sqrt(dist_sq + NUMERICAL_EPS) で修正し、TopologicalLoss の全バッチ失敗を RuntimeError にする。 Closes #59 Co-authored-by: Cursor --- REPRODUCIBILITY.md | 55 +++++++++++++++++++++++++ configs/issue59_verify_mahalanobis.yaml | 48 +++++++++++++++++++++ tda_ml/losses.py | 7 +++- tda_ml/models.py | 16 ++++--- tda_ml/numerical_eps.py | 13 ++++++ tda_ml/topology.py | 19 +++++---- tests/test_models.py | 10 +++++ tests/test_topology.py | 10 +++-- tests/test_topology_gradients.py | 32 ++++++++++++++ 9 files changed, 192 insertions(+), 18 deletions(-) create mode 100644 configs/issue59_verify_mahalanobis.yaml create mode 100644 tda_ml/numerical_eps.py create mode 100644 tests/test_topology_gradients.py diff --git a/REPRODUCIBILITY.md b/REPRODUCIBILITY.md index aa5fb6d..2b243b8 100644 --- a/REPRODUCIBILITY.md +++ b/REPRODUCIBILITY.md @@ -58,6 +58,61 @@ uv run python experiments/run_backend_multiseed.py \ したがって、`run_backend_multiseed.py` で同じ YAML を回しても、**位相損失が見ている距離空間はバックエンド間で同一ではありません**。ここでは「同一のデータ・スケジュール・設定表面での再現パイプライン比較」を意図しており、**両バックエンドが数学的に完全に同型の重み付き距離目的関数を共有する**という読み方はしません。`ellphi` 側に Mahalanobis の確率重みに相当する項を無理に足す予定はなく、比較の解釈は本節および `README.md` の英語節（*Backend comparison: outlier-probability weighting*）に従ってください。 +## 教師あり学習の目的関数（論文 Methods 用） + +本線 `tda_ml/` の学習は **点ラベル BCE** と **clean 点群の $H_1$ 持久図との Wasserstein 教師** を併用します（`configs/reproduce.yaml` 系）。 + +\[ +\mathcal{L} += w_{\mathrm{class}}\mathcal{L}_{\mathrm{BCE}} ++ w_{\mathrm{topo}}\, W_2^2\!\bigl(\mathrm{PD}_{H_1}(D^{\theta,p}),\,\mathrm{PD}_{H_1}(X_{\mathrm{clean}})\bigr) ++ w_{\mathrm{size}}\mathcal{L}_{\mathrm{size}} ++ w_{\mathrm{aniso}}\mathcal{L}_{\mathrm{aniso}}. +\] + +**楕円パラメータ（`axis_param=legacy`、主表の既定）** — 局所 PCA の $a_{i,\mathrm{base}}, b_{i,\mathrm{base}}, \theta_{i,\mathrm{base}}$ に対し: + +\[ +a_i = a_{i,\mathrm{base}}\, e^{\Delta a_i},\quad +b_i = b_{i,\mathrm{base}}\, e^{\Delta b_i},\quad +\theta_i = \theta_{i,\mathrm{base}} + \tanh(\Delta\theta_i)\frac{\pi}{2}. +\] + +**正則化（`tda_ml/losses.py`）** + +\[ +\mathcal{L}_{\mathrm{size}} = \frac{1}{N}\sum_i (M_i^2 + m_i^2),\quad +M_i=\max(a_i,b_i),\; m_i=\min(a_i,b_i). +\] + +主表 config では `aniso_mode: linear` により +$\mathcal{L}_{\mathrm{aniso}} = \frac{1}{N}\sum_i R_i$（$R_i=M_i/m_i$）。 +図用 ablation（`ablation_localscale_try2` 等）では +`aniso_mode: barrier` として +$\mathcal{L}_{\mathrm{aniso}} = \frac{10}{N}\sum_i \mathrm{ReLU}(R_i-\tau)^2$。 + +**Mahalanobis 距離**（`tda_ml/topology.py`）は outlier 確率 $p_i$ により二乗距離を +$1/\bigl((1-p_i)(1-p_j)\bigr)$ で重み付け（`INLIER_PROB_MIN` で下限クリップ）。 + +**数値安定化のみの定数**（モデリング床ではない）は `tda_ml/numerical_eps.py` に集約し、付録で列挙します。encoder の `clamp(0.2)` や legacy の `+10^{-4}` といった**論文に無い床は削除済み**（issue #59）。 + +### issue #59 検証（早期打ち切り + 診断） + +床削除後の教師あり本線を確認するには: + +```bash +uv run python experiments/issue59_verify_mahalanobis.py +``` + +- 設定: `configs/issue59_verify_mahalanobis.yaml`（`reproduce_1week_tuned` 相当・**mahalanobis**） +- **3 epoch 後**も `val_mcc` / `train_mcc` が閾値（既定 0.05 / 0.02）を超えない、または全 outlier 予測（`val_recall≈1`）なら **学習を打ち切り** +- 打ち切り時の成果物（`/logs/`）: + - `issue59_abort_report.json` / `.md` — 楕円軸・encoder PCA・距離行列・分類の統計と**原因仮説リスト** + - `abort_checkpoint.pth` — 打ち切り時点の重み + - `run_manifest.json` — コミット SHA・early_abort 設定 + +CLI: `--epochs 12 --seed 42 --out-base outputs/issue59_verify` + ## 図・定性出力 `docs/paper/` 以下の LaTeX の **すべての図を一括で出す単一スクリプトはありません**。原稿専用のアセットもあります。下の表は **コードに近い** 図の流れを示すものです。論文の図を増やしたら表も追記してください。 diff --git a/configs/issue59_verify_mahalanobis.yaml b/configs/issue59_verify_mahalanobis.yaml new file mode 100644 index 0000000..ac1264c --- /dev/null +++ b/configs/issue59_verify_mahalanobis.yaml @@ -0,0 +1,48 @@ +# Issue #59 verification: reproduce_1week_tuned + mahalanobis + early MCC abort + diagnostics. +# Mahalanobis only; aborts if MCC stays ~0 after probe_epochs (default 8 — size reg needs epochs). + +meta: + config_id: "issue59_verify_mahalanobis" + +model: + point_dim: 2 + feature_dim: 128 + threshold: 0.5 + topology_loss: + distance_backend: "mahalanobis" + ellphi_differentiable: true + +loss: + w_class: 1.0 + w_topo: 0.2 + w_aniso: 0.01099204345474479 + w_size: 0.0005578582308620256 + pos_weight: 1.0 + aniso_mode: "linear" + +training: + lr: 0.0004897466143769238 + epochs: 12 + grad_clip_value: 1.0 + visualize_every: 1 + warmup_epochs: 0 + use_amp: true + early_abort: + enabled: false + +data: + max_points: 100 + num_outliers: 20 + seed: 42 + test_size: 1000 + num_workers: 4 + batch_size: 16 + pin_memory: true + +performance: + cudnn_benchmark: true + enable_tf32: true + matmul_precision: high + +outputs: + base_dir: "outputs" diff --git a/tda_ml/losses.py b/tda_ml/losses.py index 8a07423..0ae9bfc 100644 --- a/tda_ml/losses.py +++ b/tda_ml/losses.py @@ -6,6 +6,7 @@ from torch_topological.nn import VietorisRipsComplex, WassersteinDistance from tda_ml.distance_backend import compute_distance_matrix_batch +from tda_ml.numerical_eps import NUMERICAL_EPS from tda_ml.topology import compute_anisotropic_distance_matrix # Distance-mode aliases for topology loss (kept for config/test compatibility). @@ -109,7 +110,7 @@ def forward(self, params): major_axis = axes.max(dim=-1)[0] minor_axis = axes.min(dim=-1)[0] - aspect_ratios = major_axis / (minor_axis + 1e-6) + aspect_ratios = major_axis / (minor_axis + NUMERICAL_EPS) if self.mode == 'barrier': barrier_term = F.relu(aspect_ratios - self.barrier_threshold).pow(2).mean() @@ -193,6 +194,8 @@ def forward(self, points, params, logits, clean_pd_info): ) if valid_samples == 0: - return torch.tensor(0.0, device=points.device, requires_grad=True) + raise RuntimeError( + "TopologicalLoss: all batch items failed; no valid Wasserstein terms" + ) return self.weight * (total_loss / valid_samples) diff --git a/tda_ml/models.py b/tda_ml/models.py index 92190df..e26e1c0 100644 --- a/tda_ml/models.py +++ b/tda_ml/models.py @@ -1,6 +1,8 @@ import torch import torch.nn as nn +from tda_ml.numerical_eps import EIGENVALUE_FLOOR, PCA_RIDGE_EPS + class DecoupledGeometricEncoder(nn.Module): r""" Decoupled encoder: local neighborhood PCA / eigendecomposition plus separate MLP heads. @@ -66,17 +68,21 @@ def forward(self, x): centered = relative_coords - mean_neighbor cov = torch.matmul(centered.transpose(-1, -2), centered) / (self.k - 1) - - e, v = torch.linalg.eigh(cov + torch.eye(2, device=x.device) * 1e-6) + + # eigh is not implemented for float16; local PCA in float32 is standard. + cov32 = cov.float() + eye2 = torch.eye(2, device=x.device, dtype=torch.float32) + e, v = torch.linalg.eigh(cov32 + eye2 * PCA_RIDGE_EPS) v1 = v[:, :, :, 1] base_angle = torch.atan2(v1[:, :, 1], v1[:, :, 0]).unsqueeze(-1) - base_axes = torch.sqrt(torch.clamp(e, min=1e-6)) + base_axes = torch.sqrt(torch.clamp(e, min=EIGENVALUE_FLOOR)) base_axes = torch.flip(base_axes, dims=[-1]) - base_axes = base_axes / (base_axes.max(dim=-1, keepdim=True)[0] + 1e-6) - base_axes = torch.clamp(base_axes, min=0.2) + base_axes = base_axes / ( + base_axes.max(dim=-1, keepdim=True)[0] + EIGENVALUE_FLOOR + ) local_coords_canonical = torch.matmul(centered, v) diff --git a/tda_ml/numerical_eps.py b/tda_ml/numerical_eps.py new file mode 100644 index 0000000..fc6f49b --- /dev/null +++ b/tda_ml/numerical_eps.py @@ -0,0 +1,13 @@ +"""Named numerical floors (document in paper Methods / supplement).""" + +# Denominators in metric tensor, probability weighting, aspect ratio. +NUMERICAL_EPS = 1e-8 + +# Ridge on local 2x2 covariance before ``eigh``. +PCA_RIDGE_EPS = 1e-6 + +# Floor before ``sqrt`` on PCA eigenvalues (non-negative spectrum). +EIGENVALUE_FLOOR = 1e-6 + +# Lower clamp on inlier probability in Mahalanobis distance weighting. +INLIER_PROB_MIN = 1e-4 diff --git a/tda_ml/topology.py b/tda_ml/topology.py index 99efa7b..f2f529c 100644 --- a/tda_ml/topology.py +++ b/tda_ml/topology.py @@ -2,6 +2,8 @@ import torch +from tda_ml.numerical_eps import INLIER_PROB_MIN, NUMERICAL_EPS + def compute_anisotropic_metric( axes: torch.Tensor, @@ -20,15 +22,15 @@ def compute_anisotropic_metric( Returns: tuple: (m00, m11, m01) components of the metric tensors. """ - a = axes[..., 0] + 1e-6 - b = axes[..., 1] + 1e-6 + a = axes[..., 0] + b = axes[..., 1] theta = angles.squeeze(-1) cos_t = torch.cos(theta) sin_t = torch.sin(theta) - inv_a2 = 1.0 / (a**2 + 1e-8) - inv_b2 = 1.0 / (b**2 + 1e-8) + inv_a2 = 1.0 / (a**2 + NUMERICAL_EPS) + inv_b2 = 1.0 / (b**2 + NUMERICAL_EPS) # Components of M = R diag(a^-2, b^-2) R^T m00 = cos_t**2 * inv_a2 + sin_t**2 * inv_b2 @@ -106,9 +108,9 @@ def compute_anisotropic_distance_matrix( # Weight by inlier probabilities if provided: divide squared distance by p_in[i]*p_in[j]. if probs is not None: inlier_probs = 1.0 - probs - inlier_probs = torch.clamp(inlier_probs, min=1e-4) + inlier_probs = torch.clamp(inlier_probs, min=INLIER_PROB_MIN) prob_matrix = inlier_probs.unsqueeze(2) * inlier_probs.unsqueeze(1) # (B, N, N) - dist_sq_ij = dist_sq_ij / (prob_matrix + 1e-8) + dist_sq_ij = dist_sq_ij / (prob_matrix + NUMERICAL_EPS) # Symmetrization dist_sq_ji = dist_sq_ij.transpose(-1, -2) @@ -122,5 +124,6 @@ def compute_anisotropic_distance_matrix( f"symmetrize must be 'max' or 'min', got {symmetrize!r}" ) - dist_sq = torch.clamp(dist_sq, min=0.0) - return torch.sqrt(dist_sq) + # sqrt'(0) is undefined; diagonal entries are exactly zero, so add NUMERICAL_EPS + # inside sqrt for a bounded derivative (declared L2 floor, not a modeling prior). + return torch.sqrt(dist_sq + NUMERICAL_EPS) diff --git a/tests/test_models.py b/tests/test_models.py index e05a496..cda866d 100644 --- a/tests/test_models.py +++ b/tests/test_models.py @@ -24,5 +24,15 @@ def test_forward_pass(self): self.assertIsInstance(outlier_logits, torch.Tensor) self.assertIsInstance(ellipse_params, torch.Tensor) + def test_legacy_axes_positive_without_modeling_floor(self): + """Legacy forward keeps a,b > 0 without encoder clamp or +1e-4 offset.""" + torch.manual_seed(3) + model = AnisotropicOutlierClassifier() + x = torch.rand(1, 24, 2) * 0.8 - 0.4 + with torch.no_grad(): + _, params = model(x) + axes = params[0, :, 0:2] + self.assertTrue((axes > 0).all()) + if __name__ == '__main__': unittest.main() diff --git a/tests/test_topology.py b/tests/test_topology.py index da5afcf..8031254 100644 --- a/tests/test_topology.py +++ b/tests/test_topology.py @@ -1,7 +1,9 @@ +import math import unittest import torch +from tda_ml.numerical_eps import NUMERICAL_EPS from tda_ml.topology import compute_anisotropic_distance_matrix @@ -17,12 +19,13 @@ def _make_batch(self): ) return points, params - def test_zero_diagonal_same_point(self): + def test_diagonal_uses_numerical_sqrt_floor(self): points, params = self._make_batch() dist = compute_anisotropic_distance_matrix(points, params, symmetrize="max") + expected = math.sqrt(NUMERICAL_EPS) n = points.shape[1] for i in range(n): - self.assertAlmostEqual(dist[0, i, i].item(), 0.0, places=8) + self.assertAlmostEqual(dist[0, i, i].item(), expected, places=8) def test_symmetry_max_and_min(self): points, params = self._make_batch() @@ -50,7 +53,8 @@ def test_probs_changes_off_diagonal(self): atol=1e-8, ) ) - self.assertAlmostEqual(dist_probs[0, 0, 0].item(), 0.0, places=8) + expected_diag = math.sqrt(NUMERICAL_EPS) + self.assertAlmostEqual(dist_probs[0, 0, 0].item(), expected_diag, places=8) def test_invalid_symmetrize_raises(self): points, params = self._make_batch() diff --git a/tests/test_topology_gradients.py b/tests/test_topology_gradients.py new file mode 100644 index 0000000..a8715fa --- /dev/null +++ b/tests/test_topology_gradients.py @@ -0,0 +1,32 @@ +import torch + +from tda_ml.topology import compute_anisotropic_distance_matrix + + +def test_anisotropic_distance_backward_finite_at_zero_diagonal(): + """Diagonal distances are zero; backward must not emit inf/nan gradients.""" + batch_size, num_points = 1, 5 + points = torch.randn(batch_size, num_points, 2, dtype=torch.float32) + params = torch.tensor( + [ + [ + [1.0, 0.5, 0.0], + [1.0, 0.5, 0.0], + [1.0, 0.5, 0.0], + [1.0, 0.5, 0.0], + [1.0, 0.5, 0.0], + ] + ], + dtype=torch.float32, + requires_grad=True, + ) + probs = torch.full((batch_size, num_points), 0.2, dtype=torch.float32) + + dist = compute_anisotropic_distance_matrix( + points, params, probs=probs, symmetrize="max" + ) + assert torch.all(dist.diagonal(dim1=-2, dim2=-1) > 0) + + dist.sum().backward() + assert params.grad is not None + assert torch.isfinite(params.grad).all() From 82f9528accaa4bd00f2c2b02c71cee2cb99952b8 Mon Sep 17 00:00:00 2001 From: koki3070 Date: Thu, 18 Jun 2026 10:18:49 +0900 Subject: [PATCH 2/2] =?UTF-8?q?fix:=20PR=20#62=20=E3=83=AC=E3=83=93?= =?UTF-8?q?=E3=83=A5=E3=83=BC=E5=AF=BE=E5=BF=9C=20=E2=80=94=20legacy=20+1e?= =?UTF-8?q?-4=20=E5=89=8A=E9=99=A4=E3=81=A8=20issue59=20=E6=A4=9C=E8=A8=BC?= =?UTF-8?q?=E3=82=B9=E3=82=AF=E3=83=AA=E3=83=97=E3=83=88=E8=BF=BD=E5=8A=A0?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Codex 指摘のドキュメント/実装不一致を解消。early_abort・診断を main に統合し REPRODUCIBILITY.md が案内する experiments/issue59_verify_mahalanobis.py を追加。 Co-authored-by: Cursor --- REPRODUCIBILITY.md | 4 +- configs/issue59_verify_mahalanobis.yaml | 6 +- experiments/issue59_verify_mahalanobis.py | 101 ++++++ tda_ml/main.py | 119 ++++++- tda_ml/models.py | 2 +- tda_ml/supervised_diagnostics.py | 396 ++++++++++++++++++++++ tests/test_supervised_diagnostics.py | 87 +++++ 7 files changed, 704 insertions(+), 11 deletions(-) create mode 100644 experiments/issue59_verify_mahalanobis.py create mode 100644 tda_ml/supervised_diagnostics.py create mode 100644 tests/test_supervised_diagnostics.py diff --git a/REPRODUCIBILITY.md b/REPRODUCIBILITY.md index 2b243b8..37ac3a9 100644 --- a/REPRODUCIBILITY.md +++ b/REPRODUCIBILITY.md @@ -4,7 +4,7 @@ ## リポジトリに含まれる範囲（目安） -- **含む:** `tda_ml/`、**`configs/` 直下の正本 YAML**（`base.yaml` と `reproduce` / `dev` / `prod` / `test_fast` の各ファイル）、`tests/`、追跡されている `scripts/`、`experiments/run_backend_multiseed.py`、および `README.md` / `REPRODUCIBILITY.md` / `pyproject.toml` / `uv.lock` / `LICENSE` / `CITATION.cff` などのメタデータ。 +- **含む:** `tda_ml/`、**`configs/` 直下の正本 YAML**（`base.yaml` と `reproduce` / `dev` / `prod` / `test_fast` の各ファイル）、`tests/`、追跡されている `scripts/`、`experiments/run_backend_multiseed.py`、`experiments/issue59_verify_mahalanobis.py`、および `README.md` / `REPRODUCIBILITY.md` / `pyproject.toml` / `uv.lock` / `LICENSE` / `CITATION.cff` などのメタデータ。 - **含めない:** `docs/` 以下、`configs/archive/`（履歴用 YAML を置く場合は **ローカルのみ**）、`outputs/`、`data/` など。`load_config("archive/...")` は、手元に `configs/archive/*.yaml` を置いた場合にのみ使えます。 ## 環境 @@ -105,7 +105,7 @@ uv run python experiments/issue59_verify_mahalanobis.py ``` - 設定: `configs/issue59_verify_mahalanobis.yaml`（`reproduce_1week_tuned` 相当・**mahalanobis**） -- **3 epoch 後**も `val_mcc` / `train_mcc` が閾値（既定 0.05 / 0.02）を超えない、または全 outlier 予測（`val_recall≈1`）なら **学習を打ち切り** +- **`probe_epochs` 後**（`configs/issue59_verify_mahalanobis.yaml` では既定 8）も `val_mcc` / `train_mcc` が閾値（0.05 / 0.02）を超えない、または全 outlier 予測（`val_recall≈1`）なら **学習を打ち切り** - 打ち切り時の成果物（`/logs/`）: - `issue59_abort_report.json` / `.md` — 楕円軸・encoder PCA・距離行列・分類の統計と**原因仮説リスト** - `abort_checkpoint.pth` — 打ち切り時点の重み diff --git a/configs/issue59_verify_mahalanobis.yaml b/configs/issue59_verify_mahalanobis.yaml index ac1264c..a430163 100644 --- a/configs/issue59_verify_mahalanobis.yaml +++ b/configs/issue59_verify_mahalanobis.yaml @@ -28,7 +28,11 @@ training: warmup_epochs: 0 use_amp: true early_abort: - enabled: false + enabled: true + probe_epochs: 8 + min_val_mcc: 0.05 + min_train_mcc: 0.02 + max_val_recall_for_abort: 0.99 data: max_points: 100 diff --git a/experiments/issue59_verify_mahalanobis.py b/experiments/issue59_verify_mahalanobis.py new file mode 100644 index 0000000..9f98386 --- /dev/null +++ b/experiments/issue59_verify_mahalanobis.py @@ -0,0 +1,101 @@ +#!/usr/bin/env python3 +""" +Issue #59 supervised verification (mahalanobis, reproduce_1week_tuned hyperparameters). + +Runs training with early abort when val/train MCC stay near zero after probe_epochs. +On abort, writes ``logs/issue59_abort_report.{json,md}`` with ellipse / encoder / +distance / classification diagnostics and ranked failure hypotheses. + +Usage:: + + uv run python experiments/issue59_verify_mahalanobis.py + uv run python experiments/issue59_verify_mahalanobis.py --epochs 12 --seed 42 +""" + +from __future__ import annotations + +import argparse +import json +from pathlib import Path + +from tda_ml.config import deep_update, load_config +from tda_ml.main import main as train_main +from tda_ml.supervised_diagnostics import git_revision + +REPO = Path(__file__).resolve().parents[1] + + +def preflight(config_name: str) -> dict: + cfg = load_config(config_name, project_root=REPO) + data = cfg.get("data", {}) + training = cfg.get("training", {}) + early = training.get("early_abort", {}) + if not early.get("enabled", False): + raise ValueError( + f"{config_name}: training.early_abort.enabled must be true for verification" + ) + backend = ( + cfg.get("model", {}).get("topology_loss", {}).get("distance_backend", "") + ) + if backend != "mahalanobis": + raise ValueError( + f"issue59 verification expects mahalanobis backend; got {backend!r}" + ) + data_root = REPO / "data" + if not data_root.exists(): + raise FileNotFoundError( + f"MNIST data root missing: {data_root}. Run training once to download." + ) + return { + "config_id": cfg["meta"].get("config_id"), + "source_revision": git_revision(REPO), + "epochs": training.get("epochs"), + "seed": data.get("seed"), + "early_abort": early, + "distance_backend": backend, + } + + +def main() -> int: + parser = argparse.ArgumentParser(description=__doc__) + parser.add_argument( + "--config", + default="issue59_verify_mahalanobis", + help="Config name under configs/ (default: issue59_verify_mahalanobis)", + ) + parser.add_argument("--epochs", type=int, default=None) + parser.add_argument("--seed", type=int, default=None) + parser.add_argument( + "--out-base", + type=Path, + default=REPO / "outputs" / "issue59_verify", + ) + args = parser.parse_args() + + manifest_preview = preflight(args.config) + print("=== Preflight OK ===") + print(json.dumps(manifest_preview, indent=2, ensure_ascii=True)) + + cfg = load_config(args.config, project_root=REPO) + overrides: dict = {"outputs": {"base_dir": str(args.out_base)}} + if args.epochs is not None: + overrides.setdefault("training", {})["epochs"] = args.epochs + if args.seed is not None: + overrides.setdefault("data", {})["seed"] = args.seed + cfg = deep_update(cfg, overrides) + + result = train_main(config=cfg) + status = result.get("status", "completed") + print("\n=== Verification result ===") + print(f"status: {status}") + print(f"run_dir: {result.get('run_dir')}") + print(f"best_val_mcc: {result.get('best_val_mcc', result.get('val_mcc')):.4f}") + if result.get("abort_report"): + print(f"abort_report: {result['abort_report']}") + print("\nReview hypotheses in issue59_abort_report.md before long runs.") + return 2 + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) diff --git a/tda_ml/main.py b/tda_ml/main.py index f79de25..e3ce9e9 100644 --- a/tda_ml/main.py +++ b/tda_ml/main.py @@ -12,6 +12,12 @@ from tda_ml.models import AnisotropicOutlierClassifier from tda_ml.seed_utils import set_global_seed from tda_ml.runtime_profile import build_runtime_profile +from tda_ml.supervised_diagnostics import ( + git_revision, + run_abort_diagnostics, + should_early_abort, + write_abort_report, +) from tda_ml.trainer import Trainer logger = logging.getLogger(__name__) @@ -198,6 +204,26 @@ def main(config_name=None, config=None, trial=None, config_overrides=None): with open(runtime_profile_path, "w", encoding="utf-8") as f: json.dump(runtime_profile, f, ensure_ascii=True, indent=2) logger.info("Runtime profile saved: %s", runtime_profile_path) + + manifest = { + "timestamp_utc": datetime.datetime.now(datetime.timezone.utc).isoformat(), + "source_revision": git_revision(), + "config_id": config_id, + "command_entry": "tda_ml.main", + "seed": seed, + "epochs_planned": config["training"]["epochs"], + "distance_backend": config.get("model", {}) + .get("topology_loss", {}) + .get("distance_backend", "mahalanobis"), + "early_abort": config.get("training", {}).get("early_abort"), + "run_dir": run_dir, + "fallback_status": "not applicable", + } + manifest_path = os.path.join(log_dir, "run_manifest.json") + with open(manifest_path, "w", encoding="utf-8") as f: + json.dump(manifest, f, ensure_ascii=True, indent=2) + logger.info("Run manifest saved: %s", manifest_path) + with open(metrics_path, 'w', newline='') as f: writer = csv.writer(f) writer.writerow(['epoch', 'train_loss', 'train_class_loss', 'train_topo_loss', 'train_aniso_loss', 'train_size_loss', 'val_loss', 'val_recall', 'val_mcc', 'val_aniso', 'val_size']) @@ -205,16 +231,36 @@ def main(config_name=None, config=None, trial=None, config_overrides=None): # --- Training Loop --- epochs = config['training']['epochs'] save_every = config['outputs'].get('save_every', 10) - best_val_mcc = -1.0 - + best_val_mcc = -1.0 + early_abort_cfg = config.get("training", {}).get("early_abort", {}) + metrics_history: list[dict] = [] + final_status = "completed" + abort_report_path = None + for epoch in range(1, epochs + 1): # res returns (avg_loss, class_loss, topo_loss, aniso_loss, size_loss, ...) - res = trainer.train_epoch(data_loader, epoch) + res = trainer.train_epoch(data_loader, epoch) val_res = trainer.validate(val_loader) - + val_mcc = val_res[4] # MCC is at index 4 - print(f"Epoch {epoch}: Val MCC={val_mcc:.4f}, Aniso={val_res[5]:.4f}") - + train_mcc = res[10] + val_recall = val_res[1] + print(f"Epoch {epoch}: Val MCC={val_mcc:.4f}, Aniso={val_res[5]:.4f}") + + metrics_history.append( + { + "epoch": epoch, + "val_mcc": float(val_mcc), + "train_mcc": float(train_mcc), + "val_recall": float(val_recall), + "val_specificity": float(val_res[2]), + "val_loss": float(val_res[0]), + "train_loss": float(res[0]), + "val_size": float(val_res[6]), + "val_aniso": float(val_res[5]), + } + ) + # Save best model logic if val_mcc > best_val_mcc: best_val_mcc = val_mcc @@ -238,7 +284,63 @@ def main(config_name=None, config=None, trial=None, config_overrides=None): with open(metrics_path, 'a', newline='') as f: writer = csv.writer(f) writer.writerow([epoch, res[0], res[1], res[2], res[3], res[4], val_res[0], val_res[1], val_res[4], val_res[5], val_res[6]]) - + + do_abort, abort_reason = should_early_abort( + epoch=epoch, + best_val_mcc=best_val_mcc, + val_recall=val_recall, + val_mcc=val_mcc, + train_mcc=train_mcc, + early_abort_cfg=early_abort_cfg, + ) + if do_abort: + abort_ckpt = os.path.join(run_dir, "abort_checkpoint.pth") + torch.save( + { + "epoch": epoch, + "model_state_dict": model.state_dict(), + "val_mcc": val_mcc, + "best_val_mcc": best_val_mcc, + "abort_reason": abort_reason, + }, + abort_ckpt, + ) + logger.warning("Early abort at epoch %s: %s", epoch, abort_reason) + report = run_abort_diagnostics( + trainer=trainer, + model=model, + val_loader=val_loader, + device=device, + epoch=epoch, + metrics_history=metrics_history, + abort_reason=abort_reason, + ) + abort_report_path = write_abort_report(log_dir, report) + logger.warning("Abort diagnostics: %s", abort_report_path) + manifest["final_status"] = "early-aborted" + manifest["abort_epoch"] = epoch + manifest["abort_reason"] = abort_reason + manifest["abort_report"] = str(abort_report_path) + with open(manifest_path, "w", encoding="utf-8") as f: + json.dump(manifest, f, ensure_ascii=True, indent=2) + final_status = "early-aborted" + break + + if final_status == "completed": + manifest["final_status"] = "completed" + with open(manifest_path, "w", encoding="utf-8") as f: + json.dump(manifest, f, ensure_ascii=True, indent=2) + + if final_status == "early-aborted": + return { + "val_gmean": val_res[3], + "val_mcc": val_res[4], + "run_dir": run_dir, + "status": final_status, + "best_val_mcc": best_val_mcc, + "abort_report": str(abort_report_path) if abort_report_path else None, + } + final_model_path = os.path.join(run_dir, 'final_model.pth') torch.save(model.state_dict(), final_model_path) print(f"Saved final model to {final_model_path}") @@ -279,6 +381,9 @@ def main(config_name=None, config=None, trial=None, config_overrides=None): "val_gmean": val_res[3], "val_mcc": val_res[4], "run_dir": run_dir, + "status": final_status, + "best_val_mcc": best_val_mcc, + "abort_report": None, } if __name__ == "__main__": diff --git a/tda_ml/models.py b/tda_ml/models.py index e26e1c0..d3c4f14 100644 --- a/tda_ml/models.py +++ b/tda_ml/models.py @@ -146,7 +146,7 @@ def forward(self, x): raw = self.topology_head(topo_feats) axes_scale = torch.exp(raw[:, :, 0:2]) - axes = axes_scale * base_axes + 1e-4 + axes = axes_scale * base_axes angle_delta = torch.tanh(raw[:, :, 2:3]) * (torch.pi / 2) angle = base_angle + angle_delta ellipse_params = torch.cat([axes, angle], dim=2) diff --git a/tda_ml/supervised_diagnostics.py b/tda_ml/supervised_diagnostics.py new file mode 100644 index 0000000..42dc6d7 --- /dev/null +++ b/tda_ml/supervised_diagnostics.py @@ -0,0 +1,396 @@ +"""Diagnostics when supervised training MCC stays near zero (issue #59 verification).""" + +from __future__ import annotations + +import json +import subprocess +from datetime import datetime, timezone +from pathlib import Path +from typing import Any + +import torch + +from tda_ml.distance_backend import compute_distance_matrix_batch +from tda_ml.metrics import compute_recall_specificity_gmean_mcc +from tda_ml.numerical_eps import NUMERICAL_EPS +from tda_ml.topology import compute_anisotropic_metric + + +def git_revision(repo_root: Path | None = None) -> str: + root = repo_root or Path(__file__).resolve().parents[1] + try: + return ( + subprocess.check_output( + ["git", "rev-parse", "HEAD"], + cwd=root, + stderr=subprocess.DEVNULL, + text=True, + ) + .strip() + ) + except (subprocess.CalledProcessError, FileNotFoundError): + return "unknown" + + +def _finite_stats(t: torch.Tensor) -> dict[str, float | int]: + flat = t.detach().float().reshape(-1) + finite = torch.isfinite(flat) + n = int(flat.numel()) + n_finite = int(finite.sum().item()) + if n_finite == 0: + return { + "count": n, + "finite_count": 0, + "min": float("nan"), + "median": float("nan"), + "max": float("nan"), + "nan_count": n - n_finite, + } + vals = flat[finite] + return { + "count": n, + "finite_count": n_finite, + "min": float(vals.min().item()), + "median": float(vals.median().item()), + "max": float(vals.max().item()), + "nan_count": n - n_finite, + } + + +def collect_ellipse_stats(params: torch.Tensor) -> dict[str, Any]: + if params.shape[-1] not in (3, 5): + raise ValueError(f"params last dim must be 3 or 5; got {params.shape}") + axes = params[..., 2:4] if params.shape[-1] == 5 else params[..., 0:2] + a = axes[..., 0] + b = axes[..., 1] + major = torch.maximum(a, b) + minor = torch.minimum(a, b) + ratio = major / (minor + NUMERICAL_EPS) + return { + "a": _finite_stats(a), + "b": _finite_stats(b), + "major": _finite_stats(major), + "minor": _finite_stats(minor), + "aspect_ratio": _finite_stats(ratio), + "frac_below_1e-4": float((axes < 1e-4).float().mean().item()), + "frac_below_1e-3": float((axes < 1e-3).float().mean().item()), + } + + +def collect_classification_health( + logits: torch.Tensor, + labels: torch.Tensor, + *, + threshold: float, +) -> dict[str, Any]: + probs = torch.sigmoid(logits.squeeze(-1)).detach().float() + labels_flat = labels.detach().float().reshape(-1) + preds = (probs > threshold).long() + pred_np = preds.reshape(-1).cpu().numpy() + label_np = labels_flat.cpu().numpy() + recall, specificity, gmean, mcc = compute_recall_specificity_gmean_mcc( + label_np, pred_np + ) + return { + "logits": _finite_stats(logits), + "probs": _finite_stats(probs), + "pred_outlier_fraction": float(pred_np.mean()), + "label_outlier_fraction": float(label_np.mean()), + "recall": float(recall), + "specificity": float(specificity), + "gmean": float(gmean), + "mcc": float(mcc), + "threshold": threshold, + } + + +def collect_encoder_health(model, data: torch.Tensor) -> dict[str, Any]: + with torch.no_grad(): + _, _, base_angle, base_axes = model.encoder(data) + ratio = base_axes[..., 0] / (base_axes[..., 1] + NUMERICAL_EPS) + return { + "base_axes_major": _finite_stats(base_axes[..., 0]), + "base_axes_minor": _finite_stats(base_axes[..., 1]), + "base_pca_aspect_ratio": _finite_stats(ratio), + "base_angle": _finite_stats(base_angle), + "frac_base_minor_below_0.2": float( + (base_axes[..., 1] < 0.2).float().mean().item() + ), + } + + +def collect_metric_tensor_health(params: torch.Tensor) -> dict[str, Any]: + metric_params = params[..., 2:5] if params.shape[-1] == 5 else params + m00, m11, m01 = compute_anisotropic_metric( + metric_params[..., 0:2], metric_params[..., 2:3] + ) + return { + "m00": _finite_stats(m00), + "m11": _finite_stats(m11), + "m01": _finite_stats(m01), + } + + +def collect_distance_health( + points: torch.Tensor, + params: torch.Tensor, + logits: torch.Tensor, + *, + distance_backend: str, + ellphi_differentiable: bool, + threshold: float, +) -> dict[str, Any]: + probs = torch.sigmoid(logits.squeeze(-1)) + d = compute_distance_matrix_batch( + points, + params, + probs=probs, + symmetrize="max", + backend=distance_backend, + ellphi_differentiable=ellphi_differentiable, + ) + d0 = d[0] + sym_err = (d0 - d0.T).abs().max().item() + off_diag = d0[~torch.eye(d0.shape[0], dtype=torch.bool, device=d0.device)] + return { + "distance_matrix": _finite_stats(d0), + "off_diagonal": _finite_stats(off_diag), + "max_symmetry_error": float(sym_err), + "backend": distance_backend, + } + + +def _first_val_batch(val_loader, device: torch.device): + for data, labels, clean_pc in val_loader: + return ( + data.to(device), + labels.to(device), + clean_pc.to(device), + ) + raise RuntimeError("validation loader is empty") + + +def infer_failure_hypotheses(report: dict[str, Any]) -> list[dict[str, str]]: + hypotheses: list[dict[str, str]] = [] + cls = report.get("classification", {}) + ell = report.get("ellipse", {}) + enc = report.get("encoder", {}) + dist = report.get("distance", {}) + hist = report.get("metrics_history", []) + + if cls.get("pred_outlier_fraction", 0.0) > 0.95: + hypotheses.append( + { + "id": "all_outlier_predictions", + "detail": ( + "Almost all points classified as outlier " + f"(fraction={cls['pred_outlier_fraction']:.3f}). " + "Check BCE signal, threshold, logit scale." + ), + } + ) + if cls.get("pred_outlier_fraction", 0.0) < 0.05: + hypotheses.append( + { + "id": "all_inlier_predictions", + "detail": ( + "Almost no outlier predictions; specificity may be high but " + "recall near zero." + ), + } + ) + + frac_tiny = ell.get("frac_below_1e-4", 0.0) + if frac_tiny > 0.1: + hypotheses.append( + { + "id": "axis_collapse", + "detail": ( + f"{frac_tiny:.1%} of axis values below 1e-4 " + "(fixpullback-style collapse). Check w_size / removed encoder floors." + ), + } + ) + + a_max = ell.get("a", {}).get("max", 0.0) + if a_max > 100.0: + hypotheses.append( + { + "id": "axis_explosion", + "detail": f"Large axis lengths detected (a_max={a_max:.3g}). Check w_aniso / w_size.", + } + ) + + if dist.get("distance_matrix", {}).get("nan_count", 0) > 0: + hypotheses.append( + { + "id": "non_finite_distances", + "detail": "Distance matrix contains NaN/Inf; inspect metric tensor and axes.", + } + ) + + logit_std = cls.get("logits", {}).get("max", 0.0) - cls.get("logits", {}).get("min", 0.0) + if logit_std < 0.05: + hypotheses.append( + { + "id": "logit_collapse", + "detail": ( + f"Logit dynamic range ~{logit_std:.4f}; classifier head may not be learning." + ), + } + ) + + if hist: + val_mccs = [float(row.get("val_mcc", 0.0)) for row in hist] + if max(val_mccs) < 0.05: + hypotheses.append( + { + "id": "mcc_never_rose", + "detail": ( + f"val_mcc stayed below 0.05 for all {len(val_mccs)} logged epochs " + f"(max={max(val_mccs):.4f})." + ), + } + ) + + minor_med = enc.get("base_axes_minor", {}).get("median", float("nan")) + if minor_med == minor_med and minor_med < 0.05: + hypotheses.append( + { + "id": "degenerate_pca_minor_axis", + "detail": ( + f"Encoder PCA minor axis median={minor_med:.4g}; " + "local neighborhoods may be nearly collinear." + ), + } + ) + + if not hypotheses: + hypotheses.append( + { + "id": "unclear", + "detail": ( + "MCC remained low but no single dominant failure mode was flagged. " + "Inspect metrics_history and saved checkpoint." + ), + } + ) + return hypotheses + + +def run_abort_diagnostics( + *, + trainer, + model, + val_loader, + device: torch.device, + epoch: int, + metrics_history: list[dict[str, Any]], + abort_reason: str, +) -> dict[str, Any]: + data, labels, _clean_pc = _first_val_batch(val_loader, device) + with torch.no_grad(): + logits, params = model(data) + + report: dict[str, Any] = { + "status": "early-aborted", + "abort_reason": abort_reason, + "abort_epoch": epoch, + "timestamp_utc": datetime.now(timezone.utc).isoformat(), + "config_id": trainer.config["meta"].get("config_id"), + "distance_backend": trainer.distance_backend, + "metrics_history": metrics_history, + "classification": collect_classification_health( + logits, labels, threshold=trainer.threshold + ), + "ellipse": collect_ellipse_stats(params), + "encoder": collect_encoder_health(model, data), + "metric_tensor": collect_metric_tensor_health(params), + "distance": collect_distance_health( + data, + params, + logits, + distance_backend=trainer.distance_backend, + ellphi_differentiable=trainer.ellphi_differentiable, + threshold=trainer.threshold, + ), + } + report["hypotheses"] = infer_failure_hypotheses(report) + return report + + +def write_abort_report(log_dir: str | Path, report: dict[str, Any]) -> Path: + log_dir = Path(log_dir) + log_dir.mkdir(parents=True, exist_ok=True) + json_path = log_dir / "issue59_abort_report.json" + with json_path.open("w", encoding="utf-8") as f: + json.dump(report, f, ensure_ascii=True, indent=2) + + md_path = log_dir / "issue59_abort_report.md" + lines = [ + "# Supervised training early abort report", + "", + f"- **status**: {report.get('status')}", + f"- **reason**: {report.get('abort_reason')}", + f"- **epoch**: {report.get('abort_epoch')}", + f"- **backend**: {report.get('distance_backend')}", + "", + "## Hypotheses", + "", + ] + for h in report.get("hypotheses", []): + lines.append(f"- **{h['id']}**: {h['detail']}") + lines.extend(["", "## Classification", ""]) + cls = report.get("classification", {}) + lines.append( + f"- pred_outlier_fraction={cls.get('pred_outlier_fraction')}, " + f"mcc={cls.get('mcc')}, recall={cls.get('recall')}, " + f"specificity={cls.get('specificity')}" + ) + lines.extend(["", "## Ellipse axes", ""]) + ell = report.get("ellipse", {}) + lines.append( + f"- a_median={ell.get('a', {}).get('median')}, " + f"b_median={ell.get('b', {}).get('median')}, " + f"frac_below_1e-4={ell.get('frac_below_1e-4')}" + ) + lines.extend(["", "## Metrics history", ""]) + for row in report.get("metrics_history", []): + lines.append( + f"- epoch {row.get('epoch')}: val_mcc={row.get('val_mcc')}, " + f"train_mcc={row.get('train_mcc')}, val_recall={row.get('val_recall')}" + ) + md_path.write_text("\n".join(lines) + "\n", encoding="utf-8") + return json_path + + +def should_early_abort( + *, + epoch: int, + best_val_mcc: float, + val_recall: float, + val_mcc: float, + train_mcc: float, + early_abort_cfg: dict[str, Any], +) -> tuple[bool, str]: + if not early_abort_cfg.get("enabled", False): + return False, "" + probe_epochs = int(early_abort_cfg.get("probe_epochs", 3)) + if epoch < probe_epochs: + return False, "" + + min_val_mcc = float(early_abort_cfg.get("min_val_mcc", 0.05)) + min_train_mcc = float(early_abort_cfg.get("min_train_mcc", 0.02)) + max_val_recall = float(early_abort_cfg.get("max_val_recall_for_abort", 0.99)) + + if best_val_mcc < min_val_mcc and train_mcc < min_train_mcc: + return True, ( + f"after epoch {epoch}: best_val_mcc={best_val_mcc:.4f} < {min_val_mcc} " + f"and train_mcc={train_mcc:.4f} < {min_train_mcc}" + ) + if val_mcc < min_val_mcc and val_recall >= max_val_recall: + return True, ( + f"after epoch {epoch}: val_mcc={val_mcc:.4f} with val_recall={val_recall:.4f} " + f"(likely all-outlier predictions)" + ) + return False, "" diff --git a/tests/test_supervised_diagnostics.py b/tests/test_supervised_diagnostics.py new file mode 100644 index 0000000..3a4b8a8 --- /dev/null +++ b/tests/test_supervised_diagnostics.py @@ -0,0 +1,87 @@ +"""Tests for supervised training early-abort diagnostics.""" + +import unittest + +import torch + +from tda_ml.models import AnisotropicOutlierClassifier +from tda_ml.supervised_diagnostics import ( + collect_classification_health, + collect_ellipse_stats, + infer_failure_hypotheses, + should_early_abort, +) + + +class TestSupervisedDiagnostics(unittest.TestCase): + def test_should_early_abort_after_probe(self): + cfg = { + "enabled": True, + "probe_epochs": 3, + "min_val_mcc": 0.05, + "min_train_mcc": 0.02, + } + ok, _ = should_early_abort( + epoch=2, + best_val_mcc=0.0, + val_recall=1.0, + val_mcc=0.0, + train_mcc=0.0, + early_abort_cfg=cfg, + ) + self.assertFalse(ok) + + abort, reason = should_early_abort( + epoch=3, + best_val_mcc=0.0, + val_recall=1.0, + val_mcc=0.0, + train_mcc=0.0, + early_abort_cfg=cfg, + ) + self.assertTrue(abort) + self.assertIn("best_val_mcc", reason) + + def test_should_not_abort_when_mcc_rises(self): + cfg = {"enabled": True, "probe_epochs": 3, "min_val_mcc": 0.05} + abort, _ = should_early_abort( + epoch=3, + best_val_mcc=0.12, + val_recall=0.5, + val_mcc=0.12, + train_mcc=0.15, + early_abort_cfg=cfg, + ) + self.assertFalse(abort) + + def test_all_outlier_hypothesis(self): + report = { + "classification": {"pred_outlier_fraction": 0.99, "logits": {"min": 0.0, "max": 0.01}}, + "ellipse": {"frac_below_1e-4": 0.0, "a": {"max": 1.0}}, + "encoder": {"base_axes_minor": {"median": 0.5}}, + "distance": {"distance_matrix": {"nan_count": 0}}, + "metrics_history": [{"val_mcc": 0.0}], + } + hyps = infer_failure_hypotheses(report) + ids = {h["id"] for h in hyps} + self.assertIn("all_outlier_predictions", ids) + self.assertIn("mcc_never_rose", ids) + + def test_collect_ellipse_stats_on_model_output(self): + model = AnisotropicOutlierClassifier() + x = torch.rand(2, 16, 2) + with torch.no_grad(): + _, params = model(x) + stats = collect_ellipse_stats(params) + self.assertGreater(stats["a"]["median"], 0.0) + self.assertIn("aspect_ratio", stats) + + def test_collect_classification_health(self): + logits = torch.tensor([[[5.0], [-5.0], [0.0], [2.0]]]) + labels = torch.tensor([[1.0, 0.0, 0.0, 1.0]]) + health = collect_classification_health(logits, labels, threshold=0.5) + self.assertEqual(health["pred_outlier_fraction"], 0.5) + + +if __name__ == "__main__": + unittest.main()