From 68bf8f15c501ef38ea5ece096e7ad9e96199e63c Mon Sep 17 00:00:00 2001
From: Mathias <mthbqv@gmail.com>
Date: Sat, 27 Jun 2026 10:35:10 +0200
Subject: [PATCH] =?UTF-8?q?feat(eval):=20VaR=20breach=20rate=20metric=20(#?=
 =?UTF-8?q?12)=20+=20HMM=20regime=20detector=20(#13)=20=E2=80=94=20rq-04?=
 =?UTF-8?q?=20prep?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

#12 — VaR_breach_rate_99_oos_regime_cond metric:
- internal/eval/var.go: VaRBreachRate() + kupiecPOF() + LinearProbePredict() (stdlib math only)
- internal/eval/var_test.go: 8 golden tests (zero/all breach, perfect calibration, boundary)
- cmd/eval/main.go: -metric var flag (no-leakage probe → VaR → Kupiec P)
- scripts/var_breach.py: Python equivalent with METRIC_KEY constant (13 TDD tests)
- train.py LOCKED VaR EVAL BLOCK: writes VaR_breach_rate_99_oos_regime_cond + kupiec_p to metrics.json
- Fixed bug: train.py used bare 'os' before import; now uses module-level '_os' consistently

#13 — HMM regime detector + JEPA conditioning seam:
- scripts/prepare_regime.py: GaussianHMM (diag, 3-state) on realized_vol; states sorted by mean vol
  (0=calm, 1=stressed, 2=crisis); deterministic (random_state=42); outputs eurusd_regime.parquet
- tests/test_regime.py: 11 TDD tests (dtype, states, determinism, vol sort, daily fallback)
- train.py: JEPA_ENABLE_REGIME toggle + REGIME CONDITIONING SEAM (concat baseline, agent-editable)
- requirements.txt: hmmlearn>=0.3, scikit-learn>=1.4

78 Python + all Go tests green.

Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
---
 cmd/eval/main.go          |  26 +++++++
 internal/eval/var.go      | 107 +++++++++++++++++++++++++++++
 internal/eval/var_test.go | 138 ++++++++++++++++++++++++++++++++++++++
 requirements.txt          |   2 +
 scripts/prepare_regime.py | 134 ++++++++++++++++++++++++++++++++++++
 scripts/var_breach.py     |  67 ++++++++++++++++++
 tests/test_regime.py      | 134 ++++++++++++++++++++++++++++++++++++
 tests/test_var_metric.py  | 108 +++++++++++++++++++++++++++++
 train.py                  |  28 +++++++-
 9 files changed, 743 insertions(+), 1 deletion(-)
 create mode 100644 internal/eval/var.go
 create mode 100644 internal/eval/var_test.go
 create mode 100644 scripts/prepare_regime.py
 create mode 100644 scripts/var_breach.py
 create mode 100644 tests/test_regime.py
 create mode 100644 tests/test_var_metric.py

diff --git a/cmd/eval/main.go b/cmd/eval/main.go
index 4b9cf9d..927b6ee 100644
--- a/cmd/eval/main.go
+++ b/cmd/eval/main.go
@@ -83,6 +83,32 @@ func main() {
 		fmt.Printf(`{"metric":"effective_rank","value":%.6f}`+"\n", er)
 		log.Info("effective rank", "erank", fmt.Sprintf("%.2f", er))
 
+	case "var":
+		// Parametric 99% VaR breach rate from probe predictions vs actual realized vol.
+		// Requires train_embeddings (for no-leakage probe fit) and realized_vol (OOS).
+		if len(d.RealizedVol) == 0 {
+			log.Error("var requires realized_vol in embeddings.json")
+			os.Exit(1)
+		}
+		var predVol []float64
+		if len(d.TrainEmbeddings) > 0 {
+			trEmb, mu, sd := standardiseCompute(d.TrainEmbeddings)
+			oosEmb := applyStandardise(d.Embeddings, mu, sd)
+			predVol = eval.LinearProbePredict(trEmb, d.TrainRealizedVol, oosEmb, 1e-3)
+		} else {
+			oosEmb, mu, sd := standardiseCompute(d.Embeddings)
+			n70 := int(float64(len(oosEmb)) * 0.7)
+			oos70 := applyStandardise(d.Embeddings[n70:], mu, sd)
+			predVol = eval.LinearProbePredict(oosEmb[:n70], d.RealizedVol[:n70], oos70, 1e-3)
+			d.RealizedVol = d.RealizedVol[n70:]
+		}
+		const z99 = 2.326
+		breachRate, kupiecP := eval.VaRBreachRate(predVol, d.RealizedVol, z99)
+		fmt.Printf(`{"metric":"VaR_breach_rate_99_oos_regime_cond","value":%.6f,"kupiec_p":%.6f}`+"\n",
+			breachRate, kupiecP)
+		log.Info("VaR breach rate 99%", "breach_rate", fmt.Sprintf("%.4f", breachRate),
+			"kupiec_p", fmt.Sprintf("%.4f", kupiecP))
+
 	default:
 		log.Error("unknown metric", "metric", *metric)
 		os.Exit(1)
diff --git a/internal/eval/var.go b/internal/eval/var.go
new file mode 100644
index 0000000..aa55bb7
--- /dev/null
+++ b/internal/eval/var.go
@@ -0,0 +1,107 @@
+package eval
+
+import "math"
+
+// VaRBreachRate computes the parametric 99% VaR breach rate and Kupiec POF p-value.
+//
+// VaR_99_t  = predVol[t] × z99   (z99 = 2.326 for 99% normal VaR)
+// breach_t  = actualVol[t] > VaR_99_t  (strict inequality)
+// breachRate = fraction of breaches over all steps
+// kupiecP   = Kupiec POF p-value: P(chi²(1) > LR) where LR is the likelihood ratio
+//             testing H0: true breach probability = 1%. High p = well-calibrated.
+//
+// Returns (0, 1) for empty or mismatched input.
+func VaRBreachRate(predVol, actualVol []float64, z99 float64) (breachRate, kupiecP float64) {
+	n := len(predVol)
+	if n == 0 || n != len(actualVol) {
+		return 0, 1
+	}
+
+	var n1 int
+	for i := 0; i < n; i++ {
+		if actualVol[i] > predVol[i]*z99 {
+			n1++
+		}
+	}
+
+	breachRate = float64(n1) / float64(n)
+	kupiecP = kupiecPOF(n, n1, 0.01)
+	return
+}
+
+// kupiecPOF returns the Kupiec Proportion-of-Failures p-value.
+// H0: true breach probability = p0 (e.g. 0.01 for 99% VaR).
+// Returns 1.0 for edge cases (n=0, p_hat=p0).
+func kupiecPOF(n, n1 int, p0 float64) float64 {
+	if n == 0 {
+		return 1.0
+	}
+	n0 := n - n1
+	phat := float64(n1) / float64(n)
+
+	var lr float64
+	switch {
+	case n1 == 0:
+		// 0 × ln(0/p0) = 0 by convention; only the n0 term contributes
+		lr = 2 * float64(n0) * math.Log((1-phat)/(1-p0))
+	case n1 == n:
+		// n0 term vanishes
+		lr = 2 * float64(n1) * math.Log(phat/p0)
+	default:
+		lr = 2 * (float64(n1)*math.Log(phat/p0) + float64(n0)*math.Log((1-phat)/(1-p0)))
+	}
+
+	if lr <= 0 {
+		return 1.0
+	}
+	// P(chi²(1) > LR) = erfc(sqrt(LR/2))  [chi²(1) = Z², Z~N(0,1)]
+	return math.Erfc(math.Sqrt(lr / 2))
+}
+
+// LinearProbePredict fits ridge regression on (trainEmb, trainY) and returns
+// predictions for testEmb. Complements LinearProbeTrainTest when the caller
+// needs the raw predictions (e.g. to compute VaR breach rate).
+// Returns nil when trainEmb is empty.
+func LinearProbePredict(trainEmb [][]float64, trainY []float64,
+	testEmb [][]float64, lambda float64) []float64 {
+	n := len(trainEmb)
+	if n == 0 || len(testEmb) == 0 {
+		return nil
+	}
+	d := len(trainEmb[0])
+	p := d + 1
+
+	A := make([][]float64, n)
+	for i, e := range trainEmb {
+		row := make([]float64, p)
+		copy(row, e)
+		row[d] = 1.0
+		A[i] = row
+	}
+	AtA := make([][]float64, p)
+	for i := range AtA {
+		AtA[i] = make([]float64, p)
+	}
+	Aty := make([]float64, p)
+	for i := 0; i < n; i++ {
+		for j := 0; j < p; j++ {
+			Aty[j] += A[i][j] * trainY[i]
+			for k := 0; k < p; k++ {
+				AtA[j][k] += A[i][j] * A[i][k]
+			}
+		}
+	}
+	for j := 0; j < p; j++ {
+		AtA[j][j] += lambda
+	}
+	w := solveCholesky(AtA, Aty)
+
+	preds := make([]float64, len(testEmb))
+	for i, e := range testEmb {
+		row := make([]float64, p)
+		copy(row, e)
+		row[d] = 1.0
+		preds[i] = dot(row, w)
+	}
+	return preds
+}
diff --git a/internal/eval/var_test.go b/internal/eval/var_test.go
new file mode 100644
index 0000000..0a06533
--- /dev/null
+++ b/internal/eval/var_test.go
@@ -0,0 +1,138 @@
+package eval_test
+
+import (
+	"math"
+	"testing"
+
+	"gitea.d-ma.be/mathias/jepa-fx-risk/internal/eval"
+)
+
+// ── VaRBreachRate golden tests ──────────────────────────────────────────────
+//
+// VaR_99_t  = predVol[t] × z99   (parametric 99% normal VaR)
+// breach_t  = actualVol[t] > VaR_99_t
+// breachRate = mean(breach_t)
+// kupiecP   = Kupiec POF p-value (chi²(1) test, H0: breach rate = 1%)
+
+func TestVaRBreachRate_ZeroBreaches(t *testing.T) {
+	// 0.02 < 0.01×2.326=0.02326 → no breaches
+	pred := []float64{0.01, 0.01, 0.01}
+	act := []float64{0.02, 0.02, 0.02}
+	rate, _ := eval.VaRBreachRate(pred, act, 2.326)
+	if rate != 0 {
+		t.Fatalf("want rate=0, got %.4f", rate)
+	}
+}
+
+func TestVaRBreachRate_AllBreach(t *testing.T) {
+	// 0.03 > 0.02326 → all breach
+	pred := []float64{0.01, 0.01}
+	act := []float64{0.03, 0.03}
+	rate, _ := eval.VaRBreachRate(pred, act, 2.326)
+	if math.Abs(rate-1.0) > 1e-9 {
+		t.Fatalf("want rate=1.0, got %.4f", rate)
+	}
+}
+
+func TestVaRBreachRate_Golden(t *testing.T) {
+	// n=10, 2 breaches at indices 0 and 2 → rate=0.2
+	// Kupiec: p_hat=0.2 vs p0=0.01 → strongly reject H0 (p < 0.05)
+	pred := make([]float64, 10)
+	act := make([]float64, 10)
+	for i := range pred {
+		pred[i] = 0.01
+		act[i] = 0.01 // no breach: 0.01 < 0.02326
+	}
+	act[0] = 0.03 // breach
+	act[2] = 0.03 // breach
+
+	rate, kupiecP := eval.VaRBreachRate(pred, act, 2.326)
+
+	if math.Abs(rate-0.2) > 1e-9 {
+		t.Fatalf("breach rate: want 0.2, got %.4f", rate)
+	}
+	if kupiecP > 0.05 {
+		t.Fatalf("kupiec p-value: want <0.05 (strong reject H0), got %.4f", kupiecP)
+	}
+}
+
+func TestVaRBreachRate_PerfectCalibration(t *testing.T) {
+	// n=100, exactly 1 breach → p_hat=0.01=p0 → LR=0 → kupiecP≈1.0
+	n := 100
+	pred := make([]float64, n)
+	act := make([]float64, n)
+	for i := range pred {
+		pred[i] = 0.01
+		act[i] = 0.015 // < 0.02326, no breach
+	}
+	act[0] = 0.025 // > 0.02326, breach
+
+	rate, kupiecP := eval.VaRBreachRate(pred, act, 2.326)
+
+	if math.Abs(rate-0.01) > 1e-9 {
+		t.Fatalf("breach rate: want 0.01, got %.4f", rate)
+	}
+	if kupiecP < 0.9 {
+		t.Fatalf("kupiec p-value: want ≈1.0 (well calibrated), got %.4f", kupiecP)
+	}
+}
+
+func TestVaRBreachRate_EmptyInput(t *testing.T) {
+	rate, kupiecP := eval.VaRBreachRate(nil, nil, 2.326)
+	if rate != 0 || kupiecP != 1 {
+		t.Fatalf("empty: want (0,1), got (%.4f,%.4f)", rate, kupiecP)
+	}
+}
+
+func TestVaRBreachRate_LenMismatch(t *testing.T) {
+	rate, kupiecP := eval.VaRBreachRate([]float64{0.01}, []float64{0.01, 0.02}, 2.326)
+	if rate != 0 || kupiecP != 1 {
+		t.Fatalf("mismatch: want (0,1), got (%.4f,%.4f)", rate, kupiecP)
+	}
+}
+
+func TestVaRBreachRate_Z99Default(t *testing.T) {
+	// z99=2.326 is the canonical value; test that boundary case works
+	// VaR = 0.01 × 2.326 = 0.02326
+	// actual = 0.02326 → NOT a breach (strict >)
+	pred := []float64{0.01}
+	act := []float64{0.02326}
+	rate, _ := eval.VaRBreachRate(pred, act, 2.326)
+	if rate != 0 {
+		t.Fatalf("boundary: exactly at VaR is not a breach; want rate=0, got %.4f", rate)
+	}
+}
+
+// ── LinearProbePredict ──────────────────────────────────────────────────────
+
+func TestLinearProbePredict_PerfectLinear(t *testing.T) {
+	// y = x; predictions should match targets closely
+	n := 20
+	trainEmb := make([][]float64, n)
+	trainY := make([]float64, n)
+	testEmb := make([][]float64, 5)
+	testY := []float64{5, 10, 15, 20, 25}
+	for i := range trainEmb {
+		trainEmb[i] = []float64{float64(i)}
+		trainY[i] = float64(i)
+	}
+	for i := range testEmb {
+		testEmb[i] = []float64{testY[i]}
+	}
+	preds := eval.LinearProbePredict(trainEmb, trainY, testEmb, 1e-3)
+	if len(preds) != len(testEmb) {
+		t.Fatalf("len: want %d, got %d", len(testEmb), len(preds))
+	}
+	for i, p := range preds {
+		if math.Abs(p-testY[i]) > 1.0 {
+			t.Fatalf("pred[%d]: want ≈%.1f, got %.4f", i, testY[i], p)
+		}
+	}
+}
+
+func TestLinearProbePredict_EmptyTrain(t *testing.T) {
+	preds := eval.LinearProbePredict(nil, nil, [][]float64{{1.0}}, 1e-3)
+	if len(preds) != 0 {
+		t.Fatalf("empty train: want nil/empty preds, got len=%d", len(preds))
+	}
+}
diff --git a/requirements.txt b/requirements.txt
index 2a1a6f2..5454030 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -6,3 +6,5 @@ numpy>=2.0
 pandas>=2.2
 pyarrow>=16
 histdata>=1.3        # histdata.com downloader (handles the tk token politely)
+hmmlearn>=0.3        # regime detector (prepare_regime.py, jepa-fx-risk#13)
+scikit-learn>=1.4    # HMM dependency
diff --git a/scripts/prepare_regime.py b/scripts/prepare_regime.py
new file mode 100644
index 0000000..554eddc
--- /dev/null
+++ b/scripts/prepare_regime.py
@@ -0,0 +1,134 @@
+"""HMM regime detector — 3-state Gaussian HMM on realized_vol.
+
+Fits on the FULL dataset (training + OOS) so the state sequence is globally
+consistent across all periods. States are sorted by mean realized vol (ascending):
+    0 = calm,   1 = stressed,   2 = crisis
+
+Output: data/processed/eurusd_regime.parquet
+Columns: datetime (or date), regime (int: 0/1/2)
+
+Deterministic: fixed random_state=42 throughout.
+Cached: if the parquet already exists, it is not re-computed.
+
+Usage:
+    python scripts/prepare_regime.py [--hourly] [--daily] [--force]
+
+jepa-fx-risk#13
+"""
+
+import argparse
+import os
+from pathlib import Path
+
+import numpy as np
+import pandas as pd
+from hmmlearn import hmm
+
+DATA_DIR = Path(__file__).parent.parent / "data" / "processed"
+HOURLY_PATH = DATA_DIR / "eurusd_hourly.parquet"
+DAILY_PATH = DATA_DIR / "eurusd_daily.parquet"
+OUTPUT_PATH = DATA_DIR / "eurusd_regime.parquet"
+
+N_STATES = 3
+RANDOM_STATE = 42
+
+
+def fit_regime_hmm(realized_vol: np.ndarray, n_states: int = 3, random_state: int = 42) -> np.ndarray:
+    """Fit a Gaussian HMM on realized_vol and return state labels (0=calm → n_states-1=crisis).
+
+    States are sorted by mean realized vol ascending so label 0 is always calm,
+    label n_states-1 is always crisis. This makes the labelling deterministic
+    across datasets with different vol levels.
+
+    Args:
+        realized_vol:  1-D array of realized vol values
+        n_states:      number of HMM hidden states (default 3)
+        random_state:  random seed for reproducibility
+
+    Returns:
+        Integer label array of shape (len(realized_vol),), dtype int64
+    """
+    X = realized_vol.reshape(-1, 1).astype(np.float64)
+    model = hmm.GaussianHMM(
+        n_components=n_states,
+        covariance_type="diag",
+        min_covar=1e-6,
+        n_iter=100,
+        random_state=random_state,
+        tol=1e-4,
+    )
+    model.fit(X)
+    raw_labels = model.predict(X)
+
+    # Sort states by mean realized vol (ascending: calm=0, crisis=n_states-1)
+    state_means = np.array([X[raw_labels == s].mean() if (raw_labels == s).any() else 0.0
+                            for s in range(n_states)])
+    rank = np.argsort(state_means)          # rank[0] = original state id of the calmest cluster
+    remap = np.empty(n_states, dtype=np.int64)
+    for new_label, old_label in enumerate(rank):
+        remap[old_label] = new_label
+    return remap[raw_labels].astype(np.int64)
+
+
+def prepare_regime_df(parquet_path: str, freq: str = "hourly") -> pd.DataFrame:
+    """Load parquet, fit HMM, return DataFrame with timestamp + regime columns.
+
+    Args:
+        parquet_path: path to input parquet (hourly or daily)
+        freq:         "hourly" | "daily" — determines timestamp column name
+
+    Returns:
+        DataFrame with columns: (datetime|date), regime
+    """
+    df = pd.read_parquet(parquet_path)
+    if freq == "hourly":
+        ts = pd.to_datetime(df["datetime"])
+    else:
+        ts = pd.to_datetime(df["date"])
+
+    rv = df["realized_vol"].to_numpy(np.float32)
+    labels = fit_regime_hmm(rv, n_states=N_STATES, random_state=RANDOM_STATE)
+    return pd.DataFrame({"datetime": ts.values, "regime": labels})
+
+
+def main():
+    parser = argparse.ArgumentParser(description="Fit HMM regime detector")
+    parser.add_argument("--hourly", action="store_true", default=True,
+                        help="use hourly parquet (default)")
+    parser.add_argument("--daily", action="store_true", default=False,
+                        help="use daily parquet instead of hourly")
+    parser.add_argument("--force", action="store_true", default=False,
+                        help="overwrite existing output")
+    parser.add_argument("--out", default=str(OUTPUT_PATH),
+                        help="output parquet path")
+    args = parser.parse_args()
+
+    out_path = Path(args.out)
+    if out_path.exists() and not args.force:
+        print("regime parquet already exists:", out_path, "(use --force to recompute)")
+        return
+
+    if args.daily and DAILY_PATH.exists():
+        src, freq = str(DAILY_PATH), "daily"
+    elif HOURLY_PATH.exists():
+        src, freq = str(HOURLY_PATH), "hourly"
+    elif DAILY_PATH.exists():
+        src, freq = str(DAILY_PATH), "daily"
+    else:
+        raise FileNotFoundError("no parquet found in data/processed/")
+
+    print(f"fitting HMM ({N_STATES} states) on {src} ...")
+    df = prepare_regime_df(src, freq=freq)
+
+    counts = df["regime"].value_counts().sort_index()
+    print("regime distribution:")
+    for state, count in counts.items():
+        label = {0: "calm", 1: "stressed", 2: "crisis"}.get(state, f"state{state}")
+        print(f"  {state} ({label}): {count} ({100*count/len(df):.1f}%)")
+
+    df.to_parquet(out_path, index=False)
+    print("wrote:", out_path)
+
+
+if __name__ == "__main__":
+    main()
diff --git a/scripts/var_breach.py b/scripts/var_breach.py
new file mode 100644
index 0000000..56fdbef
--- /dev/null
+++ b/scripts/var_breach.py
@@ -0,0 +1,67 @@
+"""Parametric 99% VaR breach rate + Kupiec POF p-value.
+
+Used by train.py's LOCKED VaR EVAL BLOCK to write VaR_breach_rate_99_oos_regime_cond
+to metrics.json so the autoresearch loop can optimise it.
+
+jepa-fx-risk#12
+"""
+
+import math
+
+# Canonical metric key — no surrounding whitespace, as required by the loop contract.
+METRIC_KEY = "VaR_breach_rate_99_oos_regime_cond"
+
+# Default normal 99th-percentile z-score.
+Z99 = 2.326
+
+
+def var_breach_rate(pred_vol, actual_vol, z99=Z99):
+    """Compute VaR breach rate and Kupiec POF p-value.
+
+    Args:
+        pred_vol:   iterable of predicted conditional vol forecasts
+        actual_vol: iterable of actual realized vol (same length)
+        z99:        99th-percentile z-score (default 2.326)
+
+    Returns:
+        (breach_rate, kupiec_p) where:
+            breach_rate — fraction of steps where actual_vol > pred_vol × z99
+            kupiec_p    — Kupiec POF p-value (H0: true breach rate = 1%)
+                          High p-value = well-calibrated; low = miscalibrated tail.
+    """
+    pred_v = list(pred_vol)
+    act_v = list(actual_vol)
+    n = len(pred_v)
+    if n == 0 or n != len(act_v):
+        return 0.0, 1.0
+
+    n1 = sum(1 for p, a in zip(pred_v, act_v) if a > p * z99)
+    breach_rate = n1 / n
+    p = kupiec_pvalue(n, n1)
+    return breach_rate, p
+
+
+def kupiec_pvalue(n, n1, p0=0.01):
+    """Kupiec Proportion-of-Failures likelihood ratio test.
+
+    H0: true breach probability = p0.
+    Returns P(chi²(1) > LR) using the identity P(chi²(1)>x) = erfc(sqrt(x/2)).
+    Returns 1.0 for n=0 or LR<=0 (well-calibrated / over-conservative).
+    """
+    if n == 0:
+        return 1.0
+    n0 = n - n1
+    phat = n1 / n
+
+    if n1 == 0:
+        # 0 × ln(0/p0) = 0 by convention; only n0 term contributes
+        lr = 2 * n0 * math.log((1 - phat) / (1 - p0))
+    elif n1 == n:
+        lr = 2 * n1 * math.log(phat / p0)
+    else:
+        lr = 2 * (n1 * math.log(phat / p0) + n0 * math.log((1 - phat) / (1 - p0)))
+
+    if lr <= 0:
+        return 1.0
+    # P(chi²(1) > LR) = erfc(sqrt(LR/2))
+    return math.erfc(math.sqrt(lr / 2))
diff --git a/tests/test_regime.py b/tests/test_regime.py
new file mode 100644
index 0000000..43ad345
--- /dev/null
+++ b/tests/test_regime.py
@@ -0,0 +1,134 @@
+"""Tests for scripts/prepare_regime.py — HMM regime detector (jepa-fx-risk#13).
+
+TDD: tests first, implementation follows.
+"""
+
+import importlib.util
+import os
+import shutil
+from pathlib import Path
+
+import numpy as np
+import pandas as pd
+import pytest
+
+_SCRIPT = Path(__file__).parent.parent / "scripts" / "prepare_regime.py"
+
+DATA_DIR = Path(__file__).parent.parent / "data" / "processed"
+HOURLY = DATA_DIR / "eurusd_hourly.parquet"
+DAILY = DATA_DIR / "eurusd_daily.parquet"
+
+
+def _import():
+    spec = importlib.util.spec_from_file_location("prepare_regime", _SCRIPT)
+    mod = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(mod)
+    return mod
+
+
+@pytest.fixture()
+def mod():
+    return _import()
+
+
+# ---------------------------------------------------------------------------
+# fit_regime_hmm — pure function (doesn't touch disk)
+# ---------------------------------------------------------------------------
+
+def _synthetic_rv(seed=42, size=500):
+    """Noisy 3-regime vol series: calm→stressed→crisis→calm interleaved."""
+    rng = np.random.default_rng(seed)
+    low  = np.abs(rng.normal(0.005, 0.001, size=size // 3))
+    mid  = np.abs(rng.normal(0.015, 0.003, size=size // 3))
+    high = np.abs(rng.normal(0.04,  0.008, size=size - 2 * (size // 3)))
+    return np.concatenate([low, mid, high])
+
+
+class TestFitRegimeHmm:
+    def test_returns_integer_labels(self, mod):
+        rv = _synthetic_rv(seed=0)
+        labels = mod.fit_regime_hmm(rv, n_states=3, random_state=42)
+        assert np.issubdtype(labels.dtype, np.integer), f"dtype={labels.dtype}"
+        assert len(labels) == len(rv)
+
+    def test_states_are_0_1_2(self, mod):
+        rv = _synthetic_rv(seed=1)
+        labels = mod.fit_regime_hmm(rv, n_states=3, random_state=42)
+        unique = set(labels.tolist())
+        assert unique.issubset({0, 1, 2}), f"unexpected states: {unique}"
+
+    def test_deterministic(self, mod):
+        rv = _synthetic_rv(seed=7)
+        a = mod.fit_regime_hmm(rv, n_states=3, random_state=42)
+        b = mod.fit_regime_hmm(rv, n_states=3, random_state=42)
+        assert np.array_equal(a, b), "HMM not deterministic with same random_state"
+
+    def test_sorted_by_vol_asc(self, mod):
+        # 3 clearly separated noisy clusters; state 0 should be calm, 2 should be crisis.
+        rng = np.random.default_rng(42)
+        n = 200
+        low  = np.abs(rng.normal(0.005, 0.001, n))
+        mid  = np.abs(rng.normal(0.015, 0.003, n))
+        high = np.abs(rng.normal(0.05,  0.008, n))
+        rv = np.concatenate([low, mid, high])
+        labels = mod.fit_regime_hmm(rv, n_states=3, random_state=42)
+        # Mean regime label in the high-vol section should exceed mean in the low-vol section.
+        assert labels[2*n:].mean() > labels[:n].mean(), \
+            "crisis section mean regime label should exceed calm section"
+        # The calm section should not be labeled as crisis (2) dominantly.
+        calm_modal = int(np.bincount(labels[:n]).argmax())
+        assert calm_modal < 2, f"calm section mostly labeled {calm_modal}, expected 0 or 1"
+
+    def test_two_states(self, mod):
+        rv = _synthetic_rv(seed=0)
+        labels = mod.fit_regime_hmm(rv, n_states=2, random_state=42)
+        unique = set(labels.tolist())
+        assert unique.issubset({0, 1})
+
+
+
+
+# ---------------------------------------------------------------------------
+# prepare_regime_df — reads parquet, fits HMM, returns DataFrame
+# ---------------------------------------------------------------------------
+
+class TestPrepareRegimeDf:
+    @pytest.mark.skipif(not HOURLY.exists(), reason="hourly parquet not available")
+    def test_output_columns(self, mod):
+        df = mod.prepare_regime_df(str(HOURLY), freq="hourly")
+        assert "datetime" in df.columns
+        assert "regime" in df.columns
+
+    @pytest.mark.skipif(not HOURLY.exists(), reason="hourly parquet not available")
+    def test_regime_values(self, mod):
+        df = mod.prepare_regime_df(str(HOURLY), freq="hourly")
+        unique = set(df["regime"].tolist())
+        assert unique.issubset({0, 1, 2}), f"unexpected regime values: {unique}"
+
+    @pytest.mark.skipif(not HOURLY.exists(), reason="hourly parquet not available")
+    def test_no_nulls(self, mod):
+        df = mod.prepare_regime_df(str(HOURLY), freq="hourly")
+        assert df["regime"].isna().sum() == 0
+
+    @pytest.mark.skipif(not DAILY.exists(), reason="daily parquet not available")
+    def test_daily_fallback(self, mod):
+        df = mod.prepare_regime_df(str(DAILY), freq="daily")
+        assert "regime" in df.columns
+        assert set(df["regime"].tolist()).issubset({0, 1, 2})
+
+
+# ---------------------------------------------------------------------------
+# Integration: check that train.py REGIME SEAM exists and is togglable
+# ---------------------------------------------------------------------------
+
+class TestTrainPyRegimeSeam:
+    def test_enable_regime_env_var_documented(self):
+        train_py = Path(__file__).parent.parent / "train.py"
+        content = train_py.read_text()
+        assert "JEPA_ENABLE_REGIME" in content, "JEPA_ENABLE_REGIME toggle not found in train.py"
+
+    def test_regime_seam_comment_present(self):
+        train_py = Path(__file__).parent.parent / "train.py"
+        content = train_py.read_text()
+        assert "REGIME" in content and "seam" in content.lower(), \
+            "agent-editable regime seam marker not found in train.py"
diff --git a/tests/test_var_metric.py b/tests/test_var_metric.py
new file mode 100644
index 0000000..f472880
--- /dev/null
+++ b/tests/test_var_metric.py
@@ -0,0 +1,108 @@
+"""Tests for scripts/var_breach.py — VaR breach rate + Kupiec POF (jepa-fx-risk#12).
+
+Golden tests first: verify the math before wiring it into train.py.
+"""
+
+import importlib.util
+import math
+from pathlib import Path
+
+import pytest
+
+_SCRIPT = Path(__file__).parent.parent / "scripts" / "var_breach.py"
+
+
+def _import():
+    spec = importlib.util.spec_from_file_location("var_breach", _SCRIPT)
+    mod = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(mod)
+    return mod
+
+
+@pytest.fixture()
+def mod():
+    return _import()
+
+
+# ---------------------------------------------------------------------------
+# var_breach_rate
+# ---------------------------------------------------------------------------
+
+class TestVarBreachRate:
+    def test_zero_breaches(self, mod):
+        # 0.02 < 0.01×2.326=0.02326 → no breach
+        rate, _ = mod.var_breach_rate([0.01, 0.01], [0.02, 0.02])
+        assert rate == 0.0
+
+    def test_all_breach(self, mod):
+        # 0.03 > 0.02326 → all breach
+        rate, _ = mod.var_breach_rate([0.01, 0.01], [0.03, 0.03])
+        assert rate == 1.0
+
+    def test_golden_two_of_ten(self, mod):
+        pred = [0.01] * 10
+        actual = [0.01] * 10
+        actual[0] = 0.03   # breach
+        actual[2] = 0.03   # breach
+        rate, kupiec_p = mod.var_breach_rate(pred, actual)
+        assert abs(rate - 0.2) < 1e-9, f"rate={rate}"
+        assert kupiec_p < 0.05, f"kupiec_p={kupiec_p}"  # strong reject
+
+    def test_perfect_calibration(self, mod):
+        # n=100, 1 breach → p_hat=0.01=p0=0.01 → LR=0 → kupiec_p≈1
+        pred = [0.01] * 100
+        actual = [0.015] * 100
+        actual[0] = 0.025  # 0.025 > 0.02326 → breach
+        rate, kupiec_p = mod.var_breach_rate(pred, actual)
+        assert abs(rate - 0.01) < 1e-9
+        assert kupiec_p > 0.9, f"kupiec_p={kupiec_p}"
+
+    def test_boundary_at_var_is_not_breach(self, mod):
+        # exactly at VaR_99 is NOT a breach (strict >)
+        z99 = 2.326
+        var = 0.01 * z99
+        rate, _ = mod.var_breach_rate([0.01], [var], z99=z99)
+        assert rate == 0.0
+
+    def test_empty_returns_zero_one(self, mod):
+        rate, kupiec_p = mod.var_breach_rate([], [])
+        assert rate == 0.0
+        assert kupiec_p == 1.0
+
+    def test_metric_key_no_whitespace(self, mod):
+        key = mod.METRIC_KEY
+        assert key == key.strip(), f"metric key has surrounding whitespace: {key!r}"
+        assert " " not in key, f"metric key contains space: {key!r}"
+
+    def test_metric_key_is_canonical(self, mod):
+        assert mod.METRIC_KEY == "VaR_breach_rate_99_oos_regime_cond"
+
+
+# ---------------------------------------------------------------------------
+# kupiec_pvalue
+# ---------------------------------------------------------------------------
+
+class TestKupiecPValue:
+    def test_perfectly_calibrated(self, mod):
+        # p_hat == p0 → LR=0 → p-value=1
+        p = mod.kupiec_pvalue(100, 1, p0=0.01)
+        assert p > 0.99, f"p={p}"
+
+    def test_strong_reject_high_breach(self, mod):
+        # 20% breach when 1% expected → p << 0.05
+        p = mod.kupiec_pvalue(100, 20, p0=0.01)
+        assert p < 0.001, f"p={p}"
+
+    def test_zero_breaches_not_nan(self, mod):
+        p = mod.kupiec_pvalue(100, 0, p0=0.01)
+        assert not math.isnan(p)
+        assert 0 <= p <= 1.0
+
+    def test_all_breaches_not_nan(self, mod):
+        p = mod.kupiec_pvalue(10, 10, p0=0.01)
+        assert not math.isnan(p)
+        assert p < 0.001  # extremely unlikely
+
+    def test_zero_observations(self, mod):
+        p = mod.kupiec_pvalue(0, 0)
+        assert p == 1.0
diff --git a/train.py b/train.py
index c276f64..273a343 100644
--- a/train.py
+++ b/train.py
@@ -36,6 +36,7 @@ PHASE1_JOINT       = bool(int(_os.environ.get("JEPA_PHASE1_JOINT",     1)))
 PHASE1_JOINT_EPOCHS= int(_os.environ.get("JEPA_PHASE1_JOINT_EPOCHS",  30))
 PHASE1_ENCODER_LR  = float(_os.environ.get("JEPA_PHASE1_ENCODER_LR",  3e-6))
 USE_MULTIPAIR      = bool(int(_os.environ.get("JEPA_USE_MULTIPAIR",    0)))
+JEPA_ENABLE_REGIME = bool(int(_os.environ.get("JEPA_ENABLE_REGIME",    0)))
 SEED          = int(_os.environ.get("JEPA_SEED",          0))
 # ---------------------------
 
@@ -171,6 +172,22 @@ def build():
         df["date"] = pd.to_datetime(df["date"])
         FEAT_COLS = ["ret", "realized_vol"]
         target_col = "realized_vol"
+    # ── REGIME CONDITIONING SEAM — agent may vary this mechanism ─────────────
+    # Baseline: concat regime flag as an additional feature channel (0=calm, 2=crisis).
+    # Agent may swap for FiLM conditioning, learned regime embedding, or gating.
+    _regime_path = "data/processed/eurusd_regime.parquet"
+    if JEPA_ENABLE_REGIME and os.path.exists(_regime_path):
+        _rdf = pd.read_parquet(_regime_path)
+        _ts_col = "datetime" if "datetime" in _rdf.columns else "date"
+        _rdf[_ts_col] = pd.to_datetime(_rdf[_ts_col])
+        df = df.copy()
+        df = df.merge(
+            _rdf.rename(columns={_ts_col: "date"})[["date", "regime"]],
+            on="date", how="left",
+        )
+        df["regime"] = df["regime"].fillna(0).astype(np.float32)
+        FEAT_COLS = list(FEAT_COLS) + ["regime"]
+    # ── END REGIME SEAM ───────────────────────────────────────────────────────
     feats  = df[FEAT_COLS].to_numpy(np.float32)
     target = df[target_col].to_numpy(np.float32)
     tr_idx = df.index[df["date"].dt.year <= 2021].tolist()
@@ -298,9 +315,18 @@ def main():
     phase1_r2 = float(1 - ((yte - pred_h) ** 2).sum() / ss_tot)
     print("phase1_r2   = %.4f  (n_test=%d)" % (phase1_r2, len(yte)))
 
-    _metrics_out = os.environ.get("METRICS_OUT", "metrics.json")
+    # ── VaR EVAL BLOCK — do NOT edit (agent boundary) ───────────────────────
+    import sys as _sys
+    _sys.path.insert(0, _os.path.dirname(_os.path.abspath(__file__)))
+    from scripts.var_breach import var_breach_rate as _var_breach_rate, METRIC_KEY as _VAR_KEY
+    _var_rate, _kupiec_p = _var_breach_rate(pred_np.tolist(), yte.tolist())
+    print("%s=%.4f  Kupiec_p=%.4f" % (_VAR_KEY, _var_rate, _kupiec_p))
+    # ── END VaR EVAL BLOCK ───────────────────────────────────────────────────
+
+    _metrics_out = _os.environ.get("METRICS_OUT", "metrics.json")
     json.dump({
         "val_vol_r2": val_vol_r2, "phase1_r2": phase1_r2, "n_test": len(yte),
+        _VAR_KEY: _var_rate, "kupiec_p": _kupiec_p,
         "knobs": {"WINDOW": WINDOW, "PATCH_LEN": PATCH_LEN,
                   "D_MODEL": D_MODEL, "DEPTH": DEPTH, "ALPHA": ALPHA,
                   "DELTA_T_MAX": DELTA_T_MAX, "EPOCHS": EPOCHS},