fix(eval): export block uses next-period RV target (t+1) to match Python probe

Go harness reported 0.42 vs Python 0.36 because export used realized_vol[t]
(current) while Python probe used realized_vol[t+1] (next-period). Fix adds
t+1 < len(df2) guard and uses iloc[t+1] as target. Go now matches Python: 0.3585.

Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>

.gitignore

@@ -28,3 +28,9 @@ go.work.sum
 # Project-specific
 bin/
 *.templ.go
+
+# python venv (autoresearch loop)
+.venv/
+
+# downloaded + processed market data (track via DVC/MinIO, #10 — not git)
+data/

STATUS.md

@@ -0,0 +1,18 @@
+# Autoresearch STATUS
+
+| iter | val_vol_r2 | delta | action | secs | gpu | change |
+|------|-----------|-------|--------|------|-----|--------|
+| 1 | 0.3749 | +0.0928 | KEEP | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter1 |
+| 1 | 0.3011 | +0.0776 | KEEP | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter1 |
+| 2 | 0.3032 | +0.0021 | KEEP | 2s | gpu=0% vram=10054/12227MiB temp=35°C | iter2 |
+| 1 | 0.2759 | -0.0273 | revert | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter1 |
+| 2 | 0.3442 | +0.0410 | KEEP | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter2 |
+| 3 | 0.3371 | -0.0071 | revert | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter3 |
+| 4 | 0.3143 | -0.0299 | revert | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter4 |
+| 5 | 0.3355 | -0.0087 | revert | 2s | gpu=0% vram=10054/12227MiB temp=34°C | iter5 |
+| 6 | 0.2377 | -0.1065 | revert | 2s | gpu=0% vram=10054/12227MiB temp=35°C | iter6 |
+| 1 | -0.1247 | +0.0296 | KEEP | 4s | gpu=0% vram=10054/12227MiB temp=35°C | iter1 |
+| 2 | -0.1203 | +0.0044 | KEEP | 4s | gpu=0% vram=10054/12227MiB temp=35°C | iter2 |
+| 3 | -0.0716 | +0.0487 | KEEP | 4s | gpu=0% vram=10054/12227MiB temp=36°C | iter3 |
+| 4 | 0.0590 | +0.1306 | KEEP | 5s | gpu=0% vram=10054/12227MiB temp=36°C | iter4 |
+| 5 | 0.0599 | +0.0009 | KEEP | 5s | gpu=0% vram=10054/12227MiB temp=37°C | iter5 |

Taskfile.yml

@@ -5,16 +5,48 @@ tasks:
     desc: Run templ generate
     cmds: [templ generate]
   build:
-    desc: Build the binary
+    desc: Build all binaries
     deps: [generate]
-    cmds: [go build -o bin/hostexecutor ./cmd/hostexecutor]
+    cmds:
+      - go build -o bin/jepa-fx-risk ./cmd/jepa-fx-risk
+      - go build -o bin/eval ./cmd/eval
   run:
     deps: [build]
-    cmds: [./bin/hostexecutor]
+    cmds: [./bin/jepa-fx-risk]
   test:
     desc: Run all tests
     deps: [generate]
     cmds: [go test ./... -race]
+
+  data:fetch:
+    desc: "Download EUR/USD M1 from histdata (set YEARS env var)"
+    cmds: [.venv/bin/python scripts/fetch_data.py]
+  data:fetch:historical:
+    desc: "Download EUR/USD M1 2008-2018 from histdata"
+    cmds:
+      - YEARS=2008,2009,2010,2011,2012,2013,2014,2015,2016,2017,2018 .venv/bin/python scripts/fetch_data.py
+  data:prepare:daily:
+    desc: "Rebuild eurusd_daily.parquet from all M1 zips"
+    cmds: [.venv/bin/python scripts/prepare_data.py]
+  data:prepare:hourly:
+    desc: "Build eurusd_hourly.parquet from all M1 zips"
+    cmds: [.venv/bin/python scripts/prepare_hourly.py]
+  data:prepare:all:
+    desc: "Build both daily and hourly parquets"
+    deps: [data:prepare:daily, data:prepare:hourly]
+  data:test:
+    desc: "Run Python data pipeline tests"
+    cmds: [.venv/bin/python -m pytest tests/test_prepare_hourly.py tests/test_hepa.py -v]
+
+  eval:probe:
+    desc: "Run linear-probe (val_vol_r2) on embeddings from metrics.json"
+    cmds: [./bin/eval -metric probe]
+  eval:silhouette:
+    desc: "Run silhouette on embeddings vs binary HV labels"
+    cmds: [./bin/eval -metric silhouette]
+  eval:collapse:
+    desc: "Run effective-rank collapse diagnostic"
+    cmds: [./bin/eval -metric erank]
   lint:
     cmds: [golangci-lint run ./...]
   check:

cmd/eval/main.go

@@ -0,0 +1,168 @@
+// cmd/eval — CLI driver for the jepa-fx-risk evaluation harness.
+//
+//	./bin/eval -metric probe|silhouette|erank [-emb embeddings.json]
+//
+// embeddings.json format (from train.py EXPORT_EMBEDDINGS=1):
+//
+//	{
+//	  "embeddings":         [[...], ...],  // OOS frozen embeddings
+//	  "realized_vol":       [...],         // OOS target (next-day RV)
+//	  "hv_label":           [...],         // binary HV label (top-33%)
+//	  "train_embeddings":   [[...], ...],  // train-set frozen embeddings
+//	  "train_realized_vol": [...]          // train-set RV targets
+//	}
+//
+// eval:probe standardises both sets using train statistics (no leakage).
+// Falls back to internal 70/30 split of OOS if train_embeddings absent.
+package main
+
+import (
+	"encoding/json"
+	"flag"
+	"fmt"
+	"log/slog"
+	"math"
+	"os"
+
+	"gitea.d-ma.be/mathias/jepa-fx-risk/internal/eval"
+)
+
+func main() {
+	metric := flag.String("metric", "probe", "probe | silhouette | erank")
+	embFile := flag.String("emb", "embeddings.json", "path to embeddings JSON")
+	flag.Parse()
+
+	log := slog.New(slog.NewJSONHandler(os.Stdout, nil))
+
+	d, err := readJSON(*embFile)
+	if err != nil {
+		log.Error("load embeddings", "err", err)
+		os.Exit(1)
+	}
+	log.Info("loaded", "oos", len(d.Embeddings), "dim", len(d.Embeddings[0]),
+		"train", len(d.TrainEmbeddings), "metric", *metric)
+
+	switch *metric {
+	case "probe":
+		var r2 float64
+		if len(d.TrainEmbeddings) > 0 {
+			// standardise both sets using train statistics to prevent leakage
+			trEmb, mu, sd := standardiseCompute(d.TrainEmbeddings)
+			oosEmb := applyStandardise(d.Embeddings, mu, sd)
+			r2 = eval.LinearProbeTrainTest(trEmb, d.TrainRealizedVol, oosEmb, d.RealizedVol, 1e-3)
+			log.Info("probe mode", "fit_on", "train_embeddings", "eval_on", "oos")
+		} else {
+			// fallback: internal 70/30 split of OOS embeddings
+			oosEmb, mu, sd := standardiseCompute(d.Embeddings)
+			n70 := int(float64(len(oosEmb)) * 0.7)
+			oos70 := applyStandardise(d.Embeddings[n70:], mu, sd)
+			r2 = eval.LinearProbeTrainTest(oosEmb[:n70], d.RealizedVol[:n70],
+				oos70, d.RealizedVol[n70:], 1e-3)
+			log.Info("probe mode", "fit_on", "oos[0:70%]", "eval_on", "oos[70%:]")
+		}
+		fmt.Printf(`{"metric":"val_vol_r2","value":%.6f}`+"\n", r2)
+		log.Info("linear probe", "val_vol_r2", fmt.Sprintf("%.4f", r2))
+
+	case "silhouette":
+		if len(d.HVLabel) == 0 {
+			log.Error("silhouette requires hv_label in embeddings.json")
+			os.Exit(1)
+		}
+		oosEmb := standardise(d.Embeddings)
+		sil, err := eval.Silhouette(oosEmb, d.HVLabel)
+		if err != nil {
+			log.Error("silhouette", "err", err)
+			os.Exit(1)
+		}
+		fmt.Printf(`{"metric":"silhouette","value":%.6f}`+"\n", sil)
+		log.Info("silhouette", "score", fmt.Sprintf("%.4f", sil))
+
+	case "erank":
+		oosEmb := standardise(d.Embeddings)
+		er := eval.EffectiveRank(oosEmb)
+		fmt.Printf(`{"metric":"effective_rank","value":%.6f}`+"\n", er)
+		log.Info("effective rank", "erank", fmt.Sprintf("%.2f", er))
+
+	default:
+		log.Error("unknown metric", "metric", *metric)
+		os.Exit(1)
+	}
+}
+
+type embJSON struct {
+	Embeddings       [][]float64 `json:"embeddings"`
+	Dates            []string    `json:"dates"`
+	RealizedVol      []float64   `json:"realized_vol"`
+	HVLabel          []int       `json:"hv_label"`
+	TrainEmbeddings  [][]float64 `json:"train_embeddings"`
+	TrainRealizedVol []float64   `json:"train_realized_vol"`
+}
+
+func readJSON(path string) (*embJSON, error) {
+	f, err := os.Open(path)
+	if err != nil {
+		return nil, fmt.Errorf("open %s: %w", path, err)
+	}
+	defer func() { _ = f.Close() }()
+	var d embJSON
+	if err := json.NewDecoder(f).Decode(&d); err != nil {
+		return nil, fmt.Errorf("decode: %w", err)
+	}
+	if len(d.Embeddings) == 0 {
+		return nil, fmt.Errorf("empty embeddings in %s", path)
+	}
+	return &d, nil
+}
+
+// standardise centres + scales to zero mean / unit std; returns normalised rows.
+func standardise(rows [][]float64) [][]float64 {
+	out, _, _ := standardiseCompute(rows)
+	return out
+}
+
+// standardiseCompute centres + scales and returns (normalised, mu, sd) for reuse.
+func standardiseCompute(rows [][]float64) ([][]float64, []float64, []float64) {
+	if len(rows) == 0 {
+		return rows, nil, nil
+	}
+	n, dim := len(rows), len(rows[0])
+	mu := make([]float64, dim)
+	for _, r := range rows {
+		for j, v := range r {
+			mu[j] += v
+		}
+	}
+	for j := range mu {
+		mu[j] /= float64(n)
+	}
+	sd := make([]float64, dim)
+	for _, r := range rows {
+		for j, v := range r {
+			diff := v - mu[j]
+			sd[j] += diff * diff
+		}
+	}
+	for j := range sd {
+		sd[j] = math.Sqrt(sd[j]/float64(n)) + 1e-8
+	}
+	out := make([][]float64, n)
+	for i, r := range rows {
+		out[i] = make([]float64, dim)
+		for j, v := range r {
+			out[i][j] = (v - mu[j]) / sd[j]
+		}
+	}
+	return out, mu, sd
+}
+
+// applyStandardise normalises rows using pre-computed mu and sd.
+func applyStandardise(rows [][]float64, mu, sd []float64) [][]float64 {
+	out := make([][]float64, len(rows))
+	for i, r := range rows {
+		out[i] = make([]float64, len(r))
+		for j, v := range r {
+			out[i][j] = (v - mu[j]) / sd[j]
+		}
+	}
+	return out
+}

cmd/jepa-fx-risk/main.go

@@ -5,7 +5,7 @@ import (
 	"net/http"
 	"os"
 
-	"gitea.d-ma.be/mathias/hostexecutor/internal/web"
+	"gitea.d-ma.be/mathias/jepa-fx-risk/internal/web"
 )
 
 func main() {

go.mod

@@ -1,7 +1,5 @@
-module gitea.d-ma.be/mathias/hostexecutor
+module gitea.d-ma.be/mathias/jepa-fx-risk
 
 go 1.26
 
-require (
+require github.com/a-h/templ v0.3.1020
-	github.com/a-h/templ v0.2.778
-)

go.sum

@@ -0,0 +1,4 @@
+github.com/a-h/templ v0.3.1020 h1:ypAT/L5ySWEnZ6Zft/5yfoWXYYkhFNvEFOeeqecg4tw=
+github.com/a-h/templ v0.3.1020/go.mod h1:A2DlK61v+K+NRoGnhmYbNYVmtYHcFO5/AisMvBdDxTM=
+github.com/google/go-cmp v0.6.0 h1:ofyhxvXcZhMsU5ulbFiLKl/XBFqE1GSq7atu8tAmTRI=
+github.com/google/go-cmp v0.6.0/go.mod h1:17dUlkBOakJ0+DkrSSNjCkIjxS6bF9zb3elmeNGIjoY=

internal/eval/eval.go

@@ -0,0 +1,383 @@
+// Package eval implements the Go evaluation harness for jepa-fx-risk (#4).
+// Three diagnostics on frozen embeddings exported from train.py:
+//   - LinearProbe  — val_vol_r2: OOS R² of a ridge probe predicting next-day realized vol
+//   - Silhouette   — mean silhouette score of embeddings vs a binary label (HV regime)
+//   - EffectiveRank — Roy's effective rank: exp(H(σ²)) where H is entropy of normalised singular values
+package eval
+
+import (
+	"errors"
+	"math"
+)
+
+// LinearProbeTrainTest fits ridge regression on (trainEmb, trainY) and evaluates
+// on (testEmb, testY). Returns OOS R². Use this for proper held-out evaluation.
+func LinearProbeTrainTest(trainEmb [][]float64, trainY []float64,
+	testEmb [][]float64, testY []float64, lambda float64) float64 {
+	n := len(trainEmb)
+	if n == 0 || len(testEmb) == 0 {
+		return 0
+	}
+	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)
+
+	yMean := mean(testY)
+	var ssRes, ssTot float64
+	for i, e := range testEmb {
+		row := make([]float64, p)
+		copy(row, e)
+		row[d] = 1.0
+		pred := dot(row, w)
+		ssRes += (testY[i] - pred) * (testY[i] - pred)
+		ssTot += (testY[i] - yMean) * (testY[i] - yMean)
+	}
+	if ssTot == 0 {
+		return 0
+	}
+	return 1 - ssRes/ssTot
+}
+
+// LinearProbe fits a ridge regression (closed-form) on (emb, y) with regularisation λ
+// and returns R² on the same data. Call with train embeddings; probe on held-out by
+// splitting before calling.
+//
+// emb[i] is the embedding vector for sample i; y[i] is the scalar target.
+func LinearProbe(emb [][]float64, y []float64, lambda float64) float64 {
+	n := len(emb)
+	if n == 0 {
+		return 0
+	}
+	d := len(emb[0])
+
+	// Build augmented design matrix A = [emb | 1] (n × d+1)
+	A := make([][]float64, n)
+	for i, e := range emb {
+		row := make([]float64, d+1)
+		copy(row, e)
+		row[d] = 1.0
+		A[i] = row
+	}
+
+	// Normal equations: (AᵀA + λI) w = Aᵀy  (ridge)
+	p := d + 1
+	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] * y[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)
+
+	// R² = 1 - SS_res / SS_tot
+	yMean := mean(y)
+	var ssRes, ssTot float64
+	for i := 0; i < n; i++ {
+		pred := dot(A[i], w)
+		ssRes += (y[i] - pred) * (y[i] - pred)
+		ssTot += (y[i] - yMean) * (y[i] - yMean)
+	}
+	if ssTot == 0 {
+		return 0
+	}
+	return 1 - ssRes/ssTot
+}
+
+// Silhouette returns the mean silhouette coefficient of the embeddings with respect
+// to the given integer labels. Distances are Euclidean. Returns an error if fewer
+// than 2 distinct labels are present.
+func Silhouette(emb [][]float64, labels []int) (float64, error) {
+	n := len(emb)
+	if n == 0 {
+		return 0, errors.New("eval: empty embeddings")
+	}
+
+	// count distinct labels
+	labelSet := map[int]struct{}{}
+	for _, l := range labels {
+		labelSet[l] = struct{}{}
+	}
+	if len(labelSet) < 2 {
+		return 0, errors.New("eval: silhouette requires at least 2 distinct labels")
+	}
+
+	// group indices by label
+	groups := map[int][]int{}
+	for i, l := range labels {
+		groups[l] = append(groups[l], i)
+	}
+
+	var total float64
+	for i := 0; i < n; i++ {
+		li := labels[i]
+
+		// a(i) = mean intra-cluster distance
+		var aSum float64
+		inGroup := groups[li]
+		for _, j := range inGroup {
+			if j != i {
+				aSum += euclidean(emb[i], emb[j])
+			}
+		}
+		var a float64
+		if len(inGroup) > 1 {
+			a = aSum / float64(len(inGroup)-1)
+		}
+
+		// b(i) = min mean inter-cluster distance
+		b := math.MaxFloat64
+		for l, idxs := range groups {
+			if l == li {
+				continue
+			}
+			var dSum float64
+			for _, j := range idxs {
+				dSum += euclidean(emb[i], emb[j])
+			}
+			avg := dSum / float64(len(idxs))
+			if avg < b {
+				b = avg
+			}
+		}
+
+		s := (b - a) / math.Max(a, b)
+		total += s
+	}
+	return total / float64(n), nil
+}
+
+// EffectiveRank computes Roy's effective rank of the embedding matrix:
+// exp(H) where H = -∑ pᵢ log(pᵢ) is the Shannon entropy of the normalised
+// squared singular values. Returns 1 for a rank-1 matrix and ≈ dim for
+// a full-rank isotropic matrix.
+func EffectiveRank(emb [][]float64) float64 {
+	n := len(emb)
+	if n == 0 {
+		return 0
+	}
+	d := len(emb[0])
+
+	// Compute covariance-like matrix CᵀC where C is mean-centered embedding.
+	mu := make([]float64, d)
+	for _, e := range emb {
+		for j, v := range e {
+			mu[j] += v
+		}
+	}
+	for j := range mu {
+		mu[j] /= float64(n)
+	}
+
+	// C = emb - mu  (n × d); compute CᵀC (d × d)
+	CtC := make([][]float64, d)
+	for i := range CtC {
+		CtC[i] = make([]float64, d)
+	}
+	for _, e := range emb {
+		for j := 0; j < d; j++ {
+			cj := e[j] - mu[j]
+			for k := 0; k < d; k++ {
+				CtC[j][k] += cj * (e[k] - mu[k])
+			}
+		}
+	}
+
+	// Eigenvalues of CᵀC via power iteration approximation isn't great;
+	// use the Frobenius / trace approach: σᵢ² ∝ eigenvalues of CᵀC.
+	// For a pure-Go impl without LAPACK: use the fact that the normalised
+	// squared singular values equal normalised eigenvalues of CᵀC.
+	// Compute them via Jacobi iteration for small d, or use the analytical
+	// formula for 2×2, or use iterative QR for general d.
+	eigs := jacobiEigenvalues(CtC)
+
+	// normalise to sum-1 distribution
+	var sumEig float64
+	for _, v := range eigs {
+		if v > 0 {
+			sumEig += v
+		}
+	}
+	if sumEig == 0 {
+		return 1
+	}
+	var H float64
+	for _, v := range eigs {
+		if v > 0 {
+			p := v / sumEig
+			H -= p * math.Log(p)
+		}
+	}
+	return math.Exp(H)
+}
+
+// ── internal helpers ──────────────────────────────────────────────────────────
+
+func euclidean(a, b []float64) float64 {
+	var s float64
+	for i := range a {
+		d := a[i] - b[i]
+		s += d * d
+	}
+	return math.Sqrt(s)
+}
+
+func dot(a, b []float64) float64 {
+	var s float64
+	for i := range a {
+		s += a[i] * b[i]
+	}
+	return s
+}
+
+func mean(y []float64) float64 {
+	var s float64
+	for _, v := range y {
+		s += v
+	}
+	return s / float64(len(y))
+}
+
+// solveCholesky solves Ax = b for symmetric positive-definite A via
+// Cholesky decomposition. Falls back to pseudo-inverse on failure.
+func solveCholesky(A [][]float64, b []float64) []float64 {
+	n := len(A)
+	// Cholesky decomposition: A = LLᵀ
+	L := make([][]float64, n)
+	for i := range L {
+		L[i] = make([]float64, n)
+	}
+	for i := 0; i < n; i++ {
+		for j := 0; j <= i; j++ {
+			s := A[i][j]
+			for k := 0; k < j; k++ {
+				s -= L[i][k] * L[j][k]
+			}
+			if i == j {
+				if s <= 0 {
+					s = 1e-12
+				}
+				L[i][j] = math.Sqrt(s)
+			} else {
+				L[i][j] = s / L[j][j]
+			}
+		}
+	}
+	// Forward substitution Ly = b
+	y := make([]float64, n)
+	for i := 0; i < n; i++ {
+		s := b[i]
+		for k := 0; k < i; k++ {
+			s -= L[i][k] * y[k]
+		}
+		y[i] = s / L[i][i]
+	}
+	// Back substitution Lᵀx = y
+	x := make([]float64, n)
+	for i := n - 1; i >= 0; i-- {
+		s := y[i]
+		for k := i + 1; k < n; k++ {
+			s -= L[k][i] * x[k]
+		}
+		x[i] = s / L[i][i]
+	}
+	return x
+}
+
+// jacobiEigenvalues returns eigenvalues of a symmetric matrix via Jacobi iteration.
+func jacobiEigenvalues(A [][]float64) []float64 {
+	n := len(A)
+	// copy
+	a := make([][]float64, n)
+	for i := range a {
+		a[i] = make([]float64, n)
+		copy(a[i], A[i])
+	}
+
+	const maxIter = 100
+	const tol = 1e-10
+	for iter := 0; iter < maxIter; iter++ {
+		// find largest off-diagonal element
+		p, q, amax := 0, 1, 0.0
+		for i := 0; i < n; i++ {
+			for j := i + 1; j < n; j++ {
+				if v := math.Abs(a[i][j]); v > amax {
+					amax = v
+					p, q = i, j
+				}
+			}
+		}
+		if amax < tol {
+			break
+		}
+		// Jacobi rotation
+		theta := 0.5 * math.Atan2(2*a[p][q], a[q][q]-a[p][p])
+		c, s := math.Cos(theta), math.Sin(theta)
+		// apply rotation
+		newA := make([][]float64, n)
+		for i := range newA {
+			newA[i] = make([]float64, n)
+			copy(newA[i], a[i])
+		}
+		app := c*c*a[p][p] + 2*c*s*a[p][q] + s*s*a[q][q]
+		aqq := s*s*a[p][p] - 2*c*s*a[p][q] + c*c*a[q][q]
+		apq := 0.0
+		newA[p][p] = app
+		newA[q][q] = aqq
+		newA[p][q] = apq
+		newA[q][p] = apq
+		for r := 0; r < n; r++ {
+			if r == p || r == q {
+				continue
+			}
+			arp := c*a[r][p] + s*a[r][q]
+			arq := -s*a[r][p] + c*a[r][q]
+			newA[r][p] = arp
+			newA[p][r] = arp
+			newA[r][q] = arq
+			newA[q][r] = arq
+		}
+		a = newA
+	}
+
+	eigs := make([]float64, n)
+	for i := range eigs {
+		eigs[i] = a[i][i]
+	}
+	return eigs
+}

internal/eval/eval_test.go

@@ -0,0 +1,138 @@
+package eval_test
+
+import (
+	"math"
+	"math/rand"
+	"testing"
+
+	"gitea.d-ma.be/mathias/jepa-fx-risk/internal/eval"
+)
+
+func seededRNG(seed int64) *rand.Rand {
+	return rand.New(rand.NewSource(seed))
+}
+
+// ── LinearProbe (val_vol_r2) ──────────────────────────────────────────────────
+
+func TestLinearProbe_Perfect(t *testing.T) {
+	n := 50
+	emb := make([][]float64, n)
+	y := make([]float64, n)
+	for i := range emb {
+		emb[i] = []float64{float64(i)}
+		y[i] = float64(i)
+	}
+	r2 := eval.LinearProbe(emb, y, 1e-3)
+	if r2 < 0.99 {
+		t.Fatalf("perfect predictor: want R²≥0.99, got %.4f", r2)
+	}
+}
+
+func TestLinearProbe_ConstantTarget(t *testing.T) {
+	n := 40
+	emb := make([][]float64, n)
+	y := make([]float64, n)
+	for i := range emb {
+		emb[i] = []float64{float64(i), float64(i * i)}
+		y[i] = 3.0
+	}
+	r2 := eval.LinearProbe(emb, y, 1e-3)
+	if r2 > 0.01 {
+		t.Fatalf("constant target: want R²≤0.01, got %.4f", r2)
+	}
+}
+
+func TestLinearProbe_NoiseEmbedding(t *testing.T) {
+	rng := seededRNG(42)
+	n := 80
+	emb := make([][]float64, n)
+	y := make([]float64, n)
+	for i := range emb {
+		emb[i] = []float64{rng.NormFloat64(), rng.NormFloat64()}
+		y[i] = float64(i)
+	}
+	r2 := eval.LinearProbe(emb, y, 1e-3)
+	if r2 > 0.10 {
+		t.Fatalf("noise embedding: want R²<0.10, got %.4f", r2)
+	}
+}
+
+// ── Silhouette ────────────────────────────────────────────────────────────────
+
+func TestSilhouette_PerfectClusters(t *testing.T) {
+	emb := make([][]float64, 40)
+	labels := make([]int, 40)
+	for i := range emb {
+		if i < 20 {
+			emb[i] = []float64{0.0, 0.0}
+			labels[i] = 0
+		} else {
+			emb[i] = []float64{1000.0, 1000.0}
+			labels[i] = 1
+		}
+	}
+	sil, err := eval.Silhouette(emb, labels)
+	if err != nil {
+		t.Fatal(err)
+	}
+	if sil < 0.95 {
+		t.Fatalf("perfect clusters: want sil≥0.95, got %.4f", sil)
+	}
+}
+
+func TestSilhouette_SingleLabel(t *testing.T) {
+	emb := [][]float64{{1, 2}, {3, 4}, {5, 6}}
+	labels := []int{0, 0, 0}
+	_, err := eval.Silhouette(emb, labels)
+	if err == nil {
+		t.Fatal("expected error for single-label input")
+	}
+}
+
+func TestSilhouette_RandomClusters(t *testing.T) {
+	rng := seededRNG(7)
+	n := 60
+	emb := make([][]float64, n)
+	labels := make([]int, n)
+	for i := range emb {
+		emb[i] = []float64{rng.NormFloat64(), rng.NormFloat64()}
+		labels[i] = i % 2
+	}
+	sil, err := eval.Silhouette(emb, labels)
+	if err != nil {
+		t.Fatal(err)
+	}
+	if math.Abs(sil) > 0.30 {
+		t.Fatalf("random clusters: want |sil|≤0.30, got %.4f", sil)
+	}
+}
+
+// ── EffectiveRank ─────────────────────────────────────────────────────────────
+
+func TestEffectiveRank_Rank1(t *testing.T) {
+	emb := make([][]float64, 30)
+	for i := range emb {
+		emb[i] = []float64{1.0, 2.0, 3.0, 4.0}
+	}
+	er := eval.EffectiveRank(emb)
+	if er > 1.5 {
+		t.Fatalf("rank-1 matrix: want erank≤1.5, got %.4f", er)
+	}
+}
+
+func TestEffectiveRank_FullRank(t *testing.T) {
+	rng := seededRNG(99)
+	dim := 8
+	emb := make([][]float64, 200)
+	for i := range emb {
+		row := make([]float64, dim)
+		for j := range row {
+			row[j] = rng.NormFloat64()
+		}
+		emb[i] = row
+	}
+	er := eval.EffectiveRank(emb)
+	if er < float64(dim)*0.7 {
+		t.Fatalf("full-rank: want erank≥%.1f, got %.4f", float64(dim)*0.7, er)
+	}
+}

internal/web/index_templ.go

@@ -0,0 +1,100 @@
+// Code generated by templ - DO NOT EDIT.
+
+// templ: version: v0.3.1020
+package web
+
+//lint:file-ignore SA4006 This context is only used if a nested component is present.
+
+import "github.com/a-h/templ"
+import templruntime "github.com/a-h/templ/runtime"
+
+func Index() templ.Component {
+	return templruntime.GeneratedTemplate(func(templ_7745c5c3_Input templruntime.GeneratedComponentInput) (templ_7745c5c3_Err error) {
+		templ_7745c5c3_W, ctx := templ_7745c5c3_Input.Writer, templ_7745c5c3_Input.Context
+		if templ_7745c5c3_CtxErr := ctx.Err(); templ_7745c5c3_CtxErr != nil {
+			return templ_7745c5c3_CtxErr
+		}
+		templ_7745c5c3_Buffer, templ_7745c5c3_IsBuffer := templruntime.GetBuffer(templ_7745c5c3_W)
+		if !templ_7745c5c3_IsBuffer {
+			defer func() {
+				templ_7745c5c3_BufErr := templruntime.ReleaseBuffer(templ_7745c5c3_Buffer)
+				if templ_7745c5c3_Err == nil {
+					templ_7745c5c3_Err = templ_7745c5c3_BufErr
+				}
+			}()
+		}
+		ctx = templ.InitializeContext(ctx)
+		templ_7745c5c3_Var1 := templ.GetChildren(ctx)
+		if templ_7745c5c3_Var1 == nil {
+			templ_7745c5c3_Var1 = templ.NopComponent
+		}
+		ctx = templ.ClearChildren(ctx)
+		templ_7745c5c3_Var2 := templruntime.GeneratedTemplate(func(templ_7745c5c3_Input templruntime.GeneratedComponentInput) (templ_7745c5c3_Err error) {
+			templ_7745c5c3_W, ctx := templ_7745c5c3_Input.Writer, templ_7745c5c3_Input.Context
+			templ_7745c5c3_Buffer, templ_7745c5c3_IsBuffer := templruntime.GetBuffer(templ_7745c5c3_W)
+			if !templ_7745c5c3_IsBuffer {
+				defer func() {
+					templ_7745c5c3_BufErr := templruntime.ReleaseBuffer(templ_7745c5c3_Buffer)
+					if templ_7745c5c3_Err == nil {
+						templ_7745c5c3_Err = templ_7745c5c3_BufErr
+					}
+				}()
+			}
+			ctx = templ.InitializeContext(ctx)
+			templ_7745c5c3_Err = templruntime.WriteString(templ_7745c5c3_Buffer, 1, "<h1 class=\"text-3xl font-semibold mb-6\">hostexecutor</h1><button hx-get=\"/api/hello\" hx-target=\"#out\" class=\"px-4 py-2 bg-slate-900 text-white rounded-md hover:bg-slate-700\">Say hello</button><div id=\"out\" class=\"mt-6 text-slate-700\"></div>")
+			if templ_7745c5c3_Err != nil {
+				return templ_7745c5c3_Err
+			}
+			return nil
+		})
+		templ_7745c5c3_Err = Layout("hostexecutor").Render(templ.WithChildren(ctx, templ_7745c5c3_Var2), templ_7745c5c3_Buffer)
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		return nil
+	})
+}
+
+func Hello(name string) templ.Component {
+	return templruntime.GeneratedTemplate(func(templ_7745c5c3_Input templruntime.GeneratedComponentInput) (templ_7745c5c3_Err error) {
+		templ_7745c5c3_W, ctx := templ_7745c5c3_Input.Writer, templ_7745c5c3_Input.Context
+		if templ_7745c5c3_CtxErr := ctx.Err(); templ_7745c5c3_CtxErr != nil {
+			return templ_7745c5c3_CtxErr
+		}
+		templ_7745c5c3_Buffer, templ_7745c5c3_IsBuffer := templruntime.GetBuffer(templ_7745c5c3_W)
+		if !templ_7745c5c3_IsBuffer {
+			defer func() {
+				templ_7745c5c3_BufErr := templruntime.ReleaseBuffer(templ_7745c5c3_Buffer)
+				if templ_7745c5c3_Err == nil {
+					templ_7745c5c3_Err = templ_7745c5c3_BufErr
+				}
+			}()
+		}
+		ctx = templ.InitializeContext(ctx)
+		templ_7745c5c3_Var3 := templ.GetChildren(ctx)
+		if templ_7745c5c3_Var3 == nil {
+			templ_7745c5c3_Var3 = templ.NopComponent
+		}
+		ctx = templ.ClearChildren(ctx)
+		templ_7745c5c3_Err = templruntime.WriteString(templ_7745c5c3_Buffer, 2, "<p>Hello, ")
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		var templ_7745c5c3_Var4 string
+		templ_7745c5c3_Var4, templ_7745c5c3_Err = templ.JoinStringErrs(name)
+		if templ_7745c5c3_Err != nil {
+			return templ.Error{Err: templ_7745c5c3_Err, FileName: `internal/web/index.templ`, Line: 15, Col: 17}
+		}
+		_, templ_7745c5c3_Err = templ_7745c5c3_Buffer.WriteString(templ.EscapeString(templ_7745c5c3_Var4))
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		templ_7745c5c3_Err = templruntime.WriteString(templ_7745c5c3_Buffer, 3, "!</p>")
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		return nil
+	})
+}
+
+var _ = templruntime.GeneratedTemplate

internal/web/layout_templ.go

@@ -0,0 +1,61 @@
+// Code generated by templ - DO NOT EDIT.
+
+// templ: version: v0.3.1020
+package web
+
+//lint:file-ignore SA4006 This context is only used if a nested component is present.
+
+import "github.com/a-h/templ"
+import templruntime "github.com/a-h/templ/runtime"
+
+func Layout(title string) templ.Component {
+	return templruntime.GeneratedTemplate(func(templ_7745c5c3_Input templruntime.GeneratedComponentInput) (templ_7745c5c3_Err error) {
+		templ_7745c5c3_W, ctx := templ_7745c5c3_Input.Writer, templ_7745c5c3_Input.Context
+		if templ_7745c5c3_CtxErr := ctx.Err(); templ_7745c5c3_CtxErr != nil {
+			return templ_7745c5c3_CtxErr
+		}
+		templ_7745c5c3_Buffer, templ_7745c5c3_IsBuffer := templruntime.GetBuffer(templ_7745c5c3_W)
+		if !templ_7745c5c3_IsBuffer {
+			defer func() {
+				templ_7745c5c3_BufErr := templruntime.ReleaseBuffer(templ_7745c5c3_Buffer)
+				if templ_7745c5c3_Err == nil {
+					templ_7745c5c3_Err = templ_7745c5c3_BufErr
+				}
+			}()
+		}
+		ctx = templ.InitializeContext(ctx)
+		templ_7745c5c3_Var1 := templ.GetChildren(ctx)
+		if templ_7745c5c3_Var1 == nil {
+			templ_7745c5c3_Var1 = templ.NopComponent
+		}
+		ctx = templ.ClearChildren(ctx)
+		templ_7745c5c3_Err = templruntime.WriteString(templ_7745c5c3_Buffer, 1, "<!doctype html><html lang=\"en\"><head><meta charset=\"utf-8\"><meta name=\"viewport\" content=\"width=device-width,initial-scale=1\"><title>")
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		var templ_7745c5c3_Var2 string
+		templ_7745c5c3_Var2, templ_7745c5c3_Err = templ.JoinStringErrs(title)
+		if templ_7745c5c3_Err != nil {
+			return templ.Error{Err: templ_7745c5c3_Err, FileName: `internal/web/layout.templ`, Line: 9, Col: 17}
+		}
+		_, templ_7745c5c3_Err = templ_7745c5c3_Buffer.WriteString(templ.EscapeString(templ_7745c5c3_Var2))
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		templ_7745c5c3_Err = templruntime.WriteString(templ_7745c5c3_Buffer, 2, "</title><script src=\"https://unpkg.com/htmx.org@2.0.0\"></script><script src=\"https://cdn.tailwindcss.com\"></script></head><body class=\"min-h-screen bg-slate-50 text-slate-900 antialiased\"><main class=\"max-w-3xl mx-auto px-6 py-12\">")
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		templ_7745c5c3_Err = templ_7745c5c3_Var1.Render(ctx, templ_7745c5c3_Buffer)
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		templ_7745c5c3_Err = templruntime.WriteString(templ_7745c5c3_Buffer, 3, "</main></body></html>")
+		if templ_7745c5c3_Err != nil {
+			return templ_7745c5c3_Err
+		}
+		return nil
+	})
+}
+
+var _ = templruntime.GeneratedTemplate

loop.py

@@ -0,0 +1,205 @@
+"""loop.py — Karpathy-style autoresearch loop for jepa-fx-risk.
+
+Agent (on iguana/berget — NOT koala, whose GPU is reserved for train.py) reads
+program.md + train.py + STATUS.md, proposes ONE change to train.py, we run it,
+keep if val_vol_r2 improved else git-revert. Appends per-iter record to STATUS.md.
+
+    LITELLM_KEY=xxx python loop.py [--iters N] [--model MODEL]
+
+Env:
+  LITELLM_KEY      — LiteLLM master key (required)
+  LITELLM_BASE     — default http://localhost:30401/v1
+  LOOP_MODEL       — default berget/gemma4-31b (non-thinking; iguana/berget only)
+  LOOP_ITERS       — default 3
+  TRAIN_TIMEOUT    — seconds per train.py run, default 120
+"""
+import argparse
+import json
+import os
+import subprocess
+import sys
+import time
+import textwrap
+from pathlib import Path
+
+import urllib.request
+
+LITELLM_BASE = os.environ.get("LITELLM_BASE", "http://localhost:30401/v1")
+LITELLM_KEY  = os.environ.get("LITELLM_KEY", "")
+LOOP_MODEL   = os.environ.get("LOOP_MODEL", "berget/gemma4-31b")
+LOOP_ITERS   = int(os.environ.get("LOOP_ITERS", "3"))
+TRAIN_TIMEOUT = int(os.environ.get("TRAIN_TIMEOUT", "120"))
+STATUS_MD    = Path("STATUS.md")
+METRICS_JSON = Path("metrics.json")
+TRAIN_PY     = Path("train.py")
+
+AGENT_SYSTEM = textwrap.dedent("""\
+    You are the autoresearch agent for jepa-fx-risk. Your job: propose ONE small,
+    targeted change to train.py to improve val_vol_r2 (OOS R² predicting 1-day
+    realized vol from frozen embeddings). Higher is better.
+
+    Rules:
+    - Return ONLY the full new content of train.py — nothing else, no explanation,
+      no markdown fence. Raw Python only.
+    - Change ONE thing at a time (one knob, one structural idea).
+    - Do NOT touch prepare_data.py, loop.py, or the data pipeline — only train.py.
+    - Do NOT add new data sources or new files.
+    - The metric is computed externally from your frozen embeddings; trust it.
+""")
+
+
+def read_file(p: Path) -> str:
+    return p.read_text() if p.exists() else ""
+
+
+def gpu_snapshot() -> str:
+    try:
+        out = subprocess.check_output(
+            ["nvidia-smi", "--query-gpu=utilization.gpu,memory.used,memory.total,temperature.gpu",
+             "--format=csv,noheader,nounits"], timeout=5, text=True
+        ).strip()
+        util, mem_used, mem_total, temp = [x.strip() for x in out.split(",")]
+        return "gpu=%s%% vram=%s/%sMiB temp=%s°C" % (util, mem_used, mem_total, temp)
+    except Exception:
+        return "gpu=N/A"
+
+
+def read_metric() -> float | None:
+    if not METRICS_JSON.exists():
+        return None
+    try:
+        return float(json.loads(METRICS_JSON.read_text())["val_vol_r2"])
+    except Exception:
+        return None
+
+
+def run_train() -> tuple[float | None, float, str]:
+    """Run train.py. Returns (val_vol_r2 or None, wall_secs, stderr_tail)."""
+    t0 = time.time()
+    gpu_before = gpu_snapshot()
+    try:
+        r = subprocess.run(
+            [sys.executable, "train.py"],
+            capture_output=True, text=True, timeout=TRAIN_TIMEOUT,
+        )
+        elapsed = time.time() - t0
+        if r.returncode != 0:
+            return None, elapsed, (r.stderr or r.stdout)[-300:]
+        metric = read_metric()
+        return metric, elapsed, ""
+    except subprocess.TimeoutExpired:
+        return None, TRAIN_TIMEOUT, "TIMEOUT"
+
+
+def call_agent(iteration: int, best_so_far: float | None) -> str:
+    """Ask the LLM agent to edit train.py. Returns new train.py content."""
+    context = "\n\n".join([
+        "# program.md\n" + read_file(Path("program.md")),
+        "# train.py (current)\n" + read_file(TRAIN_PY),
+        "# STATUS.md (history)\n" + read_file(STATUS_MD)[-2000:],
+        "# metrics.json (last run)\n" + read_file(METRICS_JSON),
+        "Iteration %d. Best val_vol_r2 so far: %s. Improve it." % (
+            iteration, "%.4f" % best_so_far if best_so_far is not None else "none yet"
+        ),
+    ])
+    payload = json.dumps({
+        "model": LOOP_MODEL,
+        "messages": [
+            {"role": "system", "content": AGENT_SYSTEM},
+            {"role": "user",   "content": context},
+        ],
+        "temperature": 0.7,
+        "max_tokens": 4096,
+    }).encode()
+    req = urllib.request.Request(
+        LITELLM_BASE + "/chat/completions",
+        data=payload,
+        headers={"Authorization": "Bearer " + LITELLM_KEY,
+                 "Content-Type": "application/json"},
+        method="POST",
+    )
+    resp = urllib.request.urlopen(req, timeout=60)
+    data = json.load(resp)
+    return data["choices"][0]["message"]["content"]
+
+
+def revert_train(original_content: str):
+    TRAIN_PY.write_text(original_content)
+
+
+def append_status(line: str):
+    with open(STATUS_MD, "a") as f:
+        f.write(line + "\n")
+
+
+def main():
+    if not LITELLM_KEY:
+        print("ERROR: set LITELLM_KEY"); sys.exit(1)
+
+    if not STATUS_MD.exists():
+        STATUS_MD.write_text("# Autoresearch STATUS\n\n| iter | val_vol_r2 | delta | action | secs | gpu | change |\n|------|-----------|-------|--------|------|-----|--------|\n")
+
+    # establish baseline
+    baseline = read_metric()
+    if baseline is None:
+        print("No metrics.json — running train.py for baseline...")
+        m, secs, err = run_train()
+        if m is None:
+            print("Baseline run failed:", err); sys.exit(1)
+        baseline = m
+        print("Baseline: val_vol_r2 = %.4f (%.1fs)" % (baseline, secs))
+
+    best = baseline
+    print("Starting loop | model=%s | iters=%d | baseline=%.4f" % (LOOP_MODEL, LOOP_ITERS, best))
+
+    for i in range(1, LOOP_ITERS + 1):
+        print("\n--- iter %d/%d ---" % (i, LOOP_ITERS))
+        original = TRAIN_PY.read_text()
+
+        print("  calling agent (%s)..." % LOOP_MODEL)
+        t_agent = time.time()
+        try:
+            new_code = call_agent(i, best)
+        except Exception as e:
+            print("  agent call failed:", e)
+            append_status("| %d | ERR | — | agent-fail | — | — | %s |" % (i, str(e)[:60]))
+            continue
+        agent_secs = time.time() - t_agent
+        print("  agent replied in %.1fs" % agent_secs)
+
+        # strip accidental markdown fences
+        if new_code.strip().startswith("```"):
+            lines = new_code.strip().splitlines()
+            new_code = "\n".join(lines[1:-1] if lines[-1].strip() == "```" else lines[1:])
+
+        TRAIN_PY.write_text(new_code)
+
+        gpu = gpu_snapshot()
+        print("  running train.py [%s]..." % gpu)
+        metric, secs, err = run_train()
+
+        if metric is None:
+            print("  train.py FAILED — reverting. err:", err[:100])
+            revert_train(original)
+            append_status("| %d | FAIL | — | revert | %.0fs | %s | run error |" % (i, secs, gpu))
+            continue
+
+        delta = metric - best
+        if metric > best:
+            best = metric
+            action = "KEEP"
+        else:
+            revert_train(original)
+            action = "revert"
+
+        summary = "| %d | %.4f | %+.4f | %s | %.0fs | %s | iter%d |" % (
+            i, metric, delta, action, secs, gpu, i)
+        append_status(summary)
+        print("  val_vol_r2=%.4f  delta=%+.4f  action=%s  [%.0fs]" % (metric, delta, action, secs))
+
+    print("\nDone. Best val_vol_r2 = %.4f (baseline was %.4f, delta %+.4f)" % (best, baseline, best - baseline))
+    print("STATUS.md updated.")
+
+
+if __name__ == "__main__":
+    main()

metrics.json

@@ -0,0 +1,14 @@
+{
+  "val_vol_r2": 0.05988483092470609,
+  "n_test": 263,
+  "knobs": {
+    "WINDOW": 60,
+    "PATCH_LEN": 5,
+    "STRIDE": 5,
+    "D_MODEL": 64,
+    "DEPTH": 2,
+    "MASK_FRAC": 0.5,
+    "SIGREG_LAM": 0.01,
+    "EPOCHS": 300
+  }
+}

requirements.txt

@@ -0,0 +1,8 @@
+# Python deps for the autoresearch loop (train.py + scripts). Install torch from
+# the cu130 index FIRST (koala Blackwell sm_120, torch 2.12.1+cu130 verified):
+#   pip install torch --index-url https://download.pytorch.org/whl/cu130
+#   pip install -r requirements.txt
+numpy>=2.0
+pandas>=2.2
+pyarrow>=16
+histdata>=1.3        # histdata.com downloader (handles the tk token politely)

results/summaries/phase-0-metrics.json

@@ -0,0 +1,13 @@
+{
+  "label": "null",
+  "mean_sil": 0.018206419112781685,
+  "pca_sil": 0.13589094579219818,
+  "spread": 0.8673340065023978,
+  "pc1_hv_corr": 0.525803392278542,
+  "per_seed": [
+    0.026790648698806763,
+    0.010999602265655994,
+    0.016829006373882294
+  ],
+  "passed": false
+}

results/summaries/phase-0-null.md

@@ -0,0 +1,25 @@
+# Phase-0 SSL feasibility gate — null
+
+**Date:** 2026-06-24
+**Path B deviation:** Daily 2019-2023 (not hourly 2008-2022); Python harness
+(not Go #4); gate metric adapted from silhouette-on-embedding to match
+available data. Go harness (#4) remains open for production experiments.
+
+## Data
+- Train: EUR/USD daily 2019-2021 (907 windows)
+- OOS: EUR/USD daily 2022-2023 (593 windows)
+- HV label: top-33% realized-vol days = high-volatility (196 days)
+
+## Results
+| | Value | Gate |
+|---|---|---|
+| TS-JEPA mean silhouette (3 seeds) | 0.0182 | > 0.20 → **False** |
+| Beats PCA baseline (0.1359) | 0.0182 | > PCA → **False** |
+| Seed stability (spread) | 86.73% | < 10% → **False** |
+| PC1/HV correlation | 0.5258 | < 0.95 → **True** |
+
+Per-seed: ['0.0268', '0.0110', '0.0168']
+
+## Verdict: **NULL**
+
+One or more gate criteria not met. See null result protocol in #5.

scripts/check_gpu.py

@@ -0,0 +1,21 @@
+"""Phase-0 compute gate (brain wiki/jepa-fx/facts/autoresearch-integration-phase1):
+PyTorch cu130 must see the koala Blackwell GPU and compute before any experiment.
+
+    python scripts/check_gpu.py   # exits 0 if the GPU is usable, 1 otherwise
+
+Note: koala shares this 12GB card with the llama-swap LLM stack. The autoresearch
+agent should run on iguana/berget models so koala's GPU stays free for train.py.
+"""
+import sys
+import torch
+
+print("torch", torch.__version__)
+if not torch.cuda.is_available():
+    print("CUDA NOT AVAILABLE — gate BLOCKED")
+    sys.exit(1)
+print("device:", torch.cuda.get_device_name(0))
+print("capability: sm_%d%d" % torch.cuda.get_device_capability(0))
+x = torch.randn(2000, 2000, device="cuda")
+(x @ x).sum().item()
+torch.cuda.synchronize()
+print("GPU matmul OK — Phase-0 compute gate GREEN")

scripts/fetch_data.py

@@ -0,0 +1,31 @@
+"""Fetch EUR/USD M1 bars from histdata.com (free, research use).
+
+Polite: one request per year, spaced; past years query month=None. Uses the
+maintained `histdata` package which handles histdata's anti-hotlink tk token.
+Output: data/raw/DAT_ASCII_EURUSD_M1_<year>.zip
+
+    YEARS=2019,2020,2021 python scripts/fetch_data.py
+"""
+import os
+import time
+
+from histdata import download_hist_data
+from histdata.api import Platform as P, TimeFrame as T
+
+YEARS = [y.strip() for y in os.environ.get("YEARS", "2019,2020,2021").split(",")]
+
+
+def main():
+    os.makedirs("data/raw", exist_ok=True)
+    for yr in YEARS:
+        f = download_hist_data(
+            year=yr, month=None, pair="eurusd",
+            platform=P.GENERIC_ASCII, time_frame=T.ONE_MINUTE,
+            output_directory="data/raw",
+        )
+        print("fetched", yr, "->", f)
+        time.sleep(2)  # be a good citizen
+
+
+if __name__ == "__main__":
+    main()

scripts/phase0_gate.py

@@ -0,0 +1,237 @@
+"""Phase-0 SSL feasibility gate (path B — Python fast-close of #5).
+
+Spec deviation documented: original spec (#5) required hourly 2008-2022 data
+and a Go eval harness (#4). Path B uses daily 2019-2023 + Python harness to
+close the gate quickly, since val_vol_r2 > 0 already demonstrates SSL
+feasibility. The Go harness (#4) remains open for production experiments.
+
+Gate criteria (from #5):
+  - Silhouette > 0.20 on held-out 2022-2023 (binary HV label: top-33% RV days)
+  - TS-JEPA silhouette > PCA baseline silhouette
+  - Rerun x3 seeds within ±10% of mean silhouette
+  - PC1/HV correlation < 0.95 (sanity: not trivially memorising the label)
+
+    python scripts/phase0_gate.py
+"""
+import json
+import math
+import os
+import numpy as np
+import pandas as pd
+import torch
+import torch.nn as nn
+from sklearn.decomposition import PCA
+from sklearn.metrics import silhouette_score
+from sklearn.preprocessing import StandardScaler
+
+SEEDS     = [0, 1, 2]
+WINDOW    = 30
+PATCH_LEN = 5
+STRIDE    = 5
+D_MODEL   = 32
+DEPTH     = 2
+N_HEADS   = 4
+EPOCHS    = 400
+LR        = 3e-4
+SIGREG_LAM = 0.5
+HV_PERCENTILE = 67   # top-33% = "high volatility"
+
+dev = "cuda" if torch.cuda.is_available() else "cpu"
+
+
+# ── SIGReg ──────────────────────────────────────────────────────────────────
+
+def sigreg(tokens: torch.Tensor, knots: int = 17) -> torch.Tensor:
+    B, T, D = tokens.shape
+    z = tokens.reshape(B * T, D).float()
+    t = torch.linspace(0, 3, knots, device=z.device, dtype=z.dtype)
+    dt = 3.0 / (knots - 1)
+    w = torch.full((knots,), 2 * dt, device=z.device, dtype=z.dtype)
+    w[0] = dt; w[-1] = dt
+    phi = torch.exp(-t.square() / 2.0)
+    A = torch.randn(D, 256, device=z.device, dtype=z.dtype)
+    A = A / A.norm(p=2, dim=0)
+    x_t = (z @ A).unsqueeze(-1) * t
+    err = (x_t.cos().mean(0) - phi).square() + x_t.sin().mean(0).square()
+    return ((err @ (w * phi)) * z.shape[0]).mean()
+
+
+# ── Encoder ──────────────────────────────────────────────────────────────────
+
+class PatchEncoder(nn.Module):
+    def __init__(self, in_feats, patch_len, stride, d_model, depth, n_heads):
+        super().__init__()
+        self.patch_len = patch_len
+        self.stride    = stride
+        self.embed = nn.Linear(patch_len * in_feats, d_model)
+        layer = nn.TransformerEncoderLayer(d_model, n_heads, 2 * d_model,
+                                           dropout=0.0, batch_first=True)
+        self.tf = nn.TransformerEncoder(layer, num_layers=depth)
+        n_patches = (WINDOW - patch_len) // stride + 1
+        pos = torch.zeros(n_patches, d_model)
+        for p in range(n_patches):
+            for i in range(0, d_model, 2):
+                pos[p, i]     = math.sin(p / 10000 ** (i / d_model))
+                if i + 1 < d_model:
+                    pos[p, i+1] = math.cos(p / 10000 ** (i / d_model))
+        self.register_buffer("pos", pos)
+
+    def forward(self, x):
+        B, W, F = x.shape
+        n_patches = (W - self.patch_len) // self.stride + 1
+        patches = torch.stack([x[:, i*self.stride:i*self.stride+self.patch_len, :]
+                                .reshape(B, -1) for i in range(n_patches)], dim=1)
+        tokens = self.embed(patches) + self.pos[:n_patches]
+        return self.tf(tokens)
+
+
+# ── Data ─────────────────────────────────────────────────────────────────────
+
+def load_data():
+    df = pd.read_parquet("data/processed/eurusd_daily.parquet").reset_index(drop=True)
+    df["date"] = pd.to_datetime(df["date"])
+    train = df[df["date"].dt.year <= 2021].copy()
+    oos   = df[df["date"].dt.year >= 2022].copy()
+
+    feats_all  = df[["ret", "realized_vol"]].to_numpy(np.float32)
+    target_all = df["realized_vol"].to_numpy(np.float32)
+    dates_all  = df["date"].values
+
+    mu = feats_all[:len(train)].mean(0)
+    sd = feats_all[:len(train)].std(0) + 1e-8
+
+    def windows(df_subset, feats_norm, dates):
+        idx_start = df.index[df["date"].isin(df_subset["date"])][0]
+        X, oos_dates, oos_rv = [], [], []
+        for t in range(idx_start + WINDOW, idx_start + len(df_subset)):
+            X.append(feats_norm[t - WINDOW:t])
+            oos_dates.append(dates[t])
+            oos_rv.append(target_all[t])
+        return np.stack(X), np.array(oos_rv), np.array(oos_dates)
+
+    feats_norm = (feats_all - mu) / sd
+
+    Xtr, rvtr, _ = windows(train, feats_norm, dates_all)
+    Xte, rvte, te_dates = windows(oos, feats_norm, dates_all)
+
+    # binary HV label: top-33% realized vol days in OOS = "high volatility"
+    hv_threshold = np.percentile(rvte, HV_PERCENTILE)
+    hv_labels = (rvte >= hv_threshold).astype(int)
+
+    return Xtr, rvtr, Xte, rvte, hv_labels
+
+
+# ── Train + embed ─────────────────────────────────────────────────────────────
+
+def train_and_embed(Xtr, Xte, seed):
+    torch.manual_seed(seed)
+    np.random.seed(seed)
+    enc  = PatchEncoder(Xtr.shape[2], PATCH_LEN, STRIDE, D_MODEL, DEPTH, N_HEADS).to(dev)
+    pred = nn.Sequential(nn.Linear(D_MODEL, D_MODEL), nn.GELU(),
+                         nn.Linear(D_MODEL, D_MODEL)).to(dev)
+    opt  = torch.optim.AdamW(list(enc.parameters()) + list(pred.parameters()), lr=LR)
+
+    Xtr_t = torch.tensor(Xtr, device=dev)
+    n_patches = (WINDOW - PATCH_LEN) // STRIDE + 1
+    n_mask = max(1, int(0.30 * n_patches))
+
+    for ep in range(EPOCHS):
+        idx_mask = torch.randperm(n_patches)[:n_mask]
+        tokens_ctx    = enc(Xtr_t)
+        tokens_target = enc(Xtr_t).detach()
+        jepa_loss = ((pred(tokens_ctx[:, idx_mask, :]) -
+                      tokens_target[:, idx_mask, :]) ** 2).mean()
+        reg = sigreg(tokens_ctx)
+        loss = jepa_loss + SIGREG_LAM * reg
+        opt.zero_grad(); loss.backward(); opt.step()
+
+    enc.eval()
+    with torch.no_grad():
+        Ete = enc(torch.tensor(Xte, device=dev)).mean(1).cpu().numpy()
+
+    return Ete
+
+
+# ── Gate ─────────────────────────────────────────────────────────────────────
+
+def pca_baseline(Xte, hv_labels):
+    flat = Xte.reshape(len(Xte), -1)
+    sc = StandardScaler().fit(flat)
+    emb = PCA(n_components=8).fit_transform(sc.transform(flat))
+    return silhouette_score(emb, hv_labels), emb
+
+
+def main():
+    os.makedirs("results/summaries", exist_ok=True)
+    Xtr, rvtr, Xte, rvte, hv_labels = load_data()
+    print(f"train={len(Xtr)} OOS={len(Xte)} HV={hv_labels.sum()}/{len(hv_labels)}")
+
+    pca_sil, pca_emb = pca_baseline(Xte, hv_labels)
+    pc1 = pca_emb[:, 0]
+    pc1_hv_corr = abs(np.corrcoef(pc1, hv_labels)[0, 1])
+    print(f"PCA baseline silhouette = {pca_sil:.4f} | PC1/HV |r| = {pc1_hv_corr:.4f}")
+
+    sils = []
+    for seed in SEEDS:
+        emb = train_and_embed(Xtr, Xte, seed)
+        sc  = StandardScaler().fit(emb)
+        sil = silhouette_score(sc.transform(emb), hv_labels)
+        sils.append(sil)
+        print(f"  seed={seed}  silhouette={sil:.4f}")
+
+    mean_sil = np.mean(sils)
+    spread   = (max(sils) - min(sils)) / mean_sil if mean_sil != 0 else 99
+
+    # gate checks
+    g_sil    = mean_sil > 0.20
+    g_beats  = mean_sil > pca_sil
+    g_stable = spread < 0.10
+    g_corr   = pc1_hv_corr < 0.95
+    passed   = all([g_sil, g_beats, g_stable, g_corr])
+
+    label = "pass" if passed else "null"
+    print(f"\nsilhouette mean={mean_sil:.4f} spread={spread:.2%} PCA={pca_sil:.4f} PC1/HV={pc1_hv_corr:.4f}")
+    print(f"gate: sil>0.20={g_sil} beats_pca={g_beats} stable={g_stable} corr<0.95={g_corr}")
+    print(f"PHASE-0: {label.upper()}")
+
+    summary = f"""# Phase-0 SSL feasibility gate — {label}
+
+**Date:** 2026-06-24
+**Path B deviation:** Daily 2019-2023 (not hourly 2008-2022); Python harness
+(not Go #4); gate metric adapted from silhouette-on-embedding to match
+available data. Go harness (#4) remains open for production experiments.
+
+## Data
+- Train: EUR/USD daily 2019-2021 ({len(Xtr)} windows)
+- OOS: EUR/USD daily 2022-2023 ({len(Xte)} windows)
+- HV label: top-{100-HV_PERCENTILE}% realized-vol days = high-volatility ({hv_labels.sum()} days)
+
+## Results
+| | Value | Gate |
+|---|---|---|
+| TS-JEPA mean silhouette (3 seeds) | {mean_sil:.4f} | > 0.20 → **{g_sil}** |
+| Beats PCA baseline ({pca_sil:.4f}) | {mean_sil:.4f} | > PCA → **{g_beats}** |
+| Seed stability (spread) | {spread:.2%} | < 10% → **{g_stable}** |
+| PC1/HV correlation | {pc1_hv_corr:.4f} | < 0.95 → **{g_corr}** |
+
+Per-seed: {[f"{s:.4f}" for s in sils]}
+
+## Verdict: **{label.upper()}**
+
+{"All 4 gate criteria met. TS-JEPA embeddings separate HV regimes significantly above PCA baseline with stable reproducibility." if passed else "One or more gate criteria not met. See null result protocol in #5."}
+"""
+    path = f"results/summaries/phase-0-{label}.md"
+    with open(path, "w") as f:
+        f.write(summary)
+    print(f"Written: {path}")
+
+    result = {"label": label, "mean_sil": mean_sil, "pca_sil": pca_sil,
+              "spread": spread, "pc1_hv_corr": pc1_hv_corr,
+              "per_seed": sils, "passed": passed}
+    with open("results/summaries/phase-0-metrics.json", "w") as f:
+        json.dump(result, f, indent=2)
+    return 0 if passed else 1
+
+
+if __name__ == "__main__":
+    raise SystemExit(main())

scripts/prepare_data.py

@@ -0,0 +1,57 @@
+"""LOCKED data pipeline (toy) — agent must NOT edit (brain Phase-1 contract).
+
+Parses histdata EUR/USD M1 zips → daily series with realized volatility (the
+val_vol_r2 target = 1-day realized vol from intraday squared returns).
+Output: data/processed/eurusd_daily.parquet  [date, close, ret, realized_vol].
+"""
+import glob
+import os
+import zipfile
+
+import numpy as np
+import pandas as pd
+
+RAW = "data/raw"
+OUT = "data/processed/eurusd_daily.parquet"
+
+
+def load_m1() -> pd.DataFrame:
+    frames = []
+    for zp in sorted(glob.glob(os.path.join(RAW, "DAT_ASCII_EURUSD_M1_*.zip"))):
+        with zipfile.ZipFile(zp) as z:
+            csv = [n for n in z.namelist() if n.endswith(".csv")][0]
+            with z.open(csv) as f:
+                df = pd.read_csv(
+                    f, sep=";", header=None,
+                    names=["dt", "open", "high", "low", "close", "vol"],
+                )
+        df["ts"] = pd.to_datetime(df["dt"], format="%Y%m%d %H%M%S")
+        frames.append(df[["ts", "close"]])
+    out = pd.concat(frames).sort_values("ts").reset_index(drop=True)
+    return out
+
+
+def main():
+    m1 = load_m1()
+    m1["r"] = np.log(m1["close"]).diff()
+    m1["day"] = m1["ts"].dt.normalize()
+    daily = m1.groupby("day").agg(
+        close=("close", "last"),
+        realized_vol=("r", lambda x: np.sqrt(np.nansum(x.values ** 2))),
+        n_min=("r", "count"),
+    ).reset_index()
+    daily = daily[daily["n_min"] > 60]  # drop thin days (holidays)
+    daily["ret"] = np.log(daily["close"]).diff()
+    daily = daily.dropna().reset_index(drop=True)
+    os.makedirs(os.path.dirname(OUT), exist_ok=True)
+    daily[["day", "close", "ret", "realized_vol"]].rename(columns={"day": "date"}).to_parquet(OUT)
+    print("rows:", len(daily), "| dates:", daily["day"].min().date(), "→", daily["day"].max().date())
+    # sanity: the COVID crash (March 2020) must show a realized-vol spike
+    rv = daily.set_index("day")["realized_vol"]
+    mar20 = rv["2020-03-01":"2020-03-31"].max()
+    typ = rv["2019-01-01":"2019-12-31"].median()
+    print("median 2019 RV: %.5f | max Mar-2020 RV: %.5f | spike x%.1f" % (typ, mar20, mar20 / typ))
+
+
+if __name__ == "__main__":
+    main()

scripts/prepare_hourly.py

@@ -0,0 +1,108 @@
+"""Prepare EUR/USD hourly OHLCV + realized vol from histdata M1 zips.
+
+Aggregates all M1 bars in data/raw/DAT_ASCII_EURUSD_M1_*.zip to hourly.
+Realized vol per hour = sqrt(sum(log-return²)) over the constituent M1 bars.
+Weekend hours are naturally absent (FX market closed Sat/Sun); NO interpolation.
+Hours with fewer than MIN_BARS M1 bars are dropped (holidays, thin sessions).
+
+Output: data/processed/eurusd_hourly.parquet
+Columns: datetime (UTC, tz-naive), close, ret (log), realized_vol
+
+    python scripts/prepare_hourly.py
+    RAW=data/raw OUT=data/processed/eurusd_hourly.parquet python scripts/prepare_hourly.py
+"""
+import glob
+import os
+import zipfile
+
+import numpy as np
+import pandas as pd
+
+RAW_DEFAULT = "data/raw"
+OUT_DEFAULT = "data/processed/eurusd_hourly.parquet"
+MIN_BARS    = 30   # drop hours thinner than this (holidays, DST boundary artefacts)
+
+
+# ── Core transformation ──────────────────────────────────────────────────────
+
+def resample_to_hourly(m1: pd.DataFrame) -> pd.DataFrame:
+    """Aggregate M1 DataFrame to hourly bars.
+
+    Args:
+        m1: DataFrame with columns ['ts' (datetime), 'close' (float)]
+
+    Returns:
+        DataFrame with columns ['datetime', 'close', 'ret', 'realized_vol']
+        sorted by datetime; hours with fewer than MIN_BARS M1 ticks dropped.
+    """
+    m1 = m1.sort_values("ts").copy()
+    m1["log_r"] = np.log(m1["close"]).diff()
+    m1["hour"]  = m1["ts"].dt.floor("h")
+
+    agg = m1.groupby("hour").agg(
+        close       = ("close",  "last"),
+        realized_vol= ("log_r",  lambda x: np.sqrt(np.nansum(x.values ** 2))),
+        n_bars      = ("log_r",  "count"),
+    ).reset_index()
+
+    agg = agg[agg["n_bars"] >= MIN_BARS].copy()
+    agg["ret"] = np.log(agg["close"]).diff()
+    agg = agg.dropna(subset=["ret"]).reset_index(drop=True)
+    agg = agg.rename(columns={"hour": "datetime"})
+    return agg[["datetime", "close", "ret", "realized_vol"]]
+
+
+def load_m1_from_zips(raw_dir: str) -> pd.DataFrame:
+    """Load and concatenate all M1 zips from raw_dir (histdata format)."""
+    pattern = os.path.join(raw_dir, "DAT_ASCII_EURUSD_M1_*.zip")
+    zips = sorted(glob.glob(pattern))
+    if not zips:
+        raise FileNotFoundError(f"No M1 zips found at {pattern}")
+    frames = []
+    for zp in zips:
+        with zipfile.ZipFile(zp) as z:
+            csv = [n for n in z.namelist() if n.endswith(".csv")][0]
+            with z.open(csv) as f:
+                df = pd.read_csv(
+                    f, sep=";", header=None,
+                    names=["dt", "open", "high", "low", "close", "vol"],
+                )
+        df["ts"] = pd.to_datetime(df["dt"], format="%Y%m%d %H%M%S")
+        frames.append(df[["ts", "close"]])
+        print(f"  loaded {os.path.basename(zp)}: {len(df):,} rows")
+    return pd.concat(frames).sort_values("ts").reset_index(drop=True)
+
+
+def build_hourly_parquet(
+    raw_dir: str = RAW_DEFAULT,
+    out_path: str = OUT_DEFAULT,
+) -> pd.DataFrame:
+    """Full pipeline: load all M1 zips → hourly parquet. Returns the DataFrame."""
+    print(f"Loading M1 zips from {raw_dir}...")
+    m1 = load_m1_from_zips(raw_dir)
+    print(f"Total M1 bars: {len(m1):,}  ({m1['ts'].min().date()} → {m1['ts'].max().date()})")
+
+    print("Resampling to hourly...")
+    hourly = resample_to_hourly(m1)
+    print(f"Hourly rows: {len(hourly):,}  ({hourly['datetime'].min()} → {hourly['datetime'].max()})")
+
+    # Sanity: COVID crash (Mar 2020) should show realized vol spike if data covers it
+    if hourly["datetime"].dt.year.isin([2020]).any():
+        rv = hourly.set_index("datetime")["realized_vol"]
+        try:
+            mar20 = rv["2020-03-01":"2020-03-31"].max()
+            typ   = rv["2019-01-01":"2019-12-31"].median()
+            print(f"Sanity — median 2019 RV: {typ:.6f} | max Mar-2020 RV: {mar20:.6f} | spike ×{mar20/typ:.1f}")
+        except Exception:
+            pass
+
+    os.makedirs(os.path.dirname(os.path.abspath(out_path)), exist_ok=True)
+    hourly.to_parquet(out_path, index=False)
+    print(f"Written: {out_path}")
+    return hourly
+
+
+if __name__ == "__main__":
+    raw_dir  = os.environ.get("RAW",  RAW_DEFAULT)
+    out_path = os.environ.get("OUT",  OUT_DEFAULT)
+    build_hourly_parquet(raw_dir=raw_dir, out_path=out_path)

tests/test_hepa.py

@@ -0,0 +1,112 @@
+"""Failing tests for HEPA backbone in train.py.
+
+Run: cd ~/dev/AI/jepa-fx-risk && .venv/bin/python -m pytest tests/test_hepa.py -v
+These tests define what the new backbone must satisfy BEFORE implementation.
+"""
+import math
+import torch
+import torch.nn as nn
+import pytest
+
+# ── Tests import the classes from train.py ────────────────────────────────────
+# They will fail until train.py implements: CausalEncoder, HorizonPredictor, vicreg_loss
+
+
+def _import():
+    import importlib.util, sys
+    spec = importlib.util.spec_from_file_location("train", "train.py")
+    mod = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(mod)
+    return mod
+
+
+@pytest.fixture(scope="module")
+def train_mod():
+    return _import()
+
+
+# 1. CausalEncoder exists and has correct output shape
+def test_causal_encoder_shape(train_mod):
+    enc = train_mod.CausalEncoder(n_channels=2, patch_len=10, d_model=32, n_heads=4, depth=1)
+    x = torch.randn(4, 60, 2)
+    tokens = enc(x)          # should return all tokens (B, N, D) for JEPA pretraining
+    assert tokens.shape == (4, 6, 32), f"expected (4, 6, 32), got {tokens.shape}"
+
+
+# 2. CausalEncoder is actually causal: earlier token outputs don't change when later inputs change
+def test_causal_masking(train_mod):
+    enc = train_mod.CausalEncoder(n_channels=2, patch_len=10, d_model=32, n_heads=4, depth=2)
+    enc.eval()
+    torch.manual_seed(0)
+    x = torch.randn(1, 60, 2)
+    x_perturbed = x.clone()
+    # non-uniform noise (constant shift absorbed by per-patch LayerNorm; variance change is not)
+    torch.manual_seed(99)
+    x_perturbed[:, 30:, :] += torch.randn_like(x[:, 30:, :]) * 5.0
+
+    with torch.no_grad():
+        h1 = enc(x)
+        h2 = enc(x_perturbed)
+
+    # First 3 tokens must be identical (causal — don't see future patches)
+    assert torch.allclose(h1[:, :3, :], h2[:, :3, :], atol=1e-5), \
+        "causal masking broken: early tokens change when later input changes"
+    # Last token should differ (it can see the perturbed patches)
+    assert not torch.allclose(h1[:, -1, :], h2[:, -1, :], atol=1e-5), \
+        "last token should differ when later input changes"
+
+
+# 3. HorizonPredictor exists, takes (h, delta_t_float) → same shape as h
+def test_horizon_predictor_shape(train_mod):
+    pred = train_mod.HorizonPredictor(d_model=32)
+    h = torch.randn(4, 32)
+    dt = torch.tensor([1.0, 2.0, 3.0, 1.0])
+    out = pred(h, dt)
+    assert out.shape == (4, 32), f"expected (4, 32), got {out.shape}"
+
+
+# 4. vicreg_loss is a scalar and backward doesn't error
+def test_vicreg_loss_backward(train_mod):
+    h_pred = torch.randn(8, 32, requires_grad=True)
+    h_target = torch.randn(8, 32)
+    loss = train_mod.vicreg_loss(h_pred, h_target, alpha=0.1)
+    assert loss.shape == (), f"expected scalar, got {loss.shape}"
+    loss.backward()
+    assert h_pred.grad is not None
+
+
+# 5. Full JEPA step: encode context, predict future, compute loss, backward
+def test_jepa_step_end_to_end(train_mod):
+    enc = train_mod.CausalEncoder(n_channels=2, patch_len=10, d_model=32, n_heads=4, depth=1)
+    pred = train_mod.HorizonPredictor(d_model=32)
+    opt = torch.optim.SGD(list(enc.parameters()) + list(pred.parameters()), lr=1e-3)
+
+    x = torch.randn(4, 60, 2)
+    tokens = enc(x)          # (4, 6, 32)
+    c, dt = 2, 2             # context position 2, horizon 2
+    h_ctx = tokens[:, c, :]
+    h_tgt = tokens[:, c + dt, :].detach()
+    h_hat = pred(h_ctx, torch.full((4,), float(dt)))
+    loss = train_mod.vicreg_loss(h_hat, h_tgt, alpha=0.1)
+    opt.zero_grad(); loss.backward(); opt.step()
+    assert loss.item() < 100, "loss exploded"
+
+
+# 6. build() returns year-based OOS split (2022-2023); hourly gives many more windows
+def test_build_year_split(train_mod):
+    (Xtr, ytr), (Xte, yte) = train_mod.build()
+    assert Xtr.shape[1] == train_mod.WINDOW
+    assert Xte.shape[1] == train_mod.WINDOW
+    assert len(Xtr) > 0 and len(Xte) > 0
+    # OOS: daily ≈ 600; hourly ≈ 17,000 (2 years × ~8,500 trading hours/year)
+    assert len(Xte) > 400, f"OOS too small: {len(Xte)}"
+
+
+# 7. hourly build gives > 10× more training windows than daily
+def test_build_hourly_more_windows(train_mod):
+    import os
+    if not os.path.exists("data/processed/eurusd_hourly.parquet"):
+        pytest.skip("eurusd_hourly.parquet not present — run data:prepare:hourly first")
+    (Xtr, _), _ = train_mod.build()
+    # Daily had ~877 train windows; hourly with 2008-2021 should have > 50,000
+    assert len(Xtr) > 10_000, f"expected >10k hourly train windows, got {len(Xtr)}"

tests/test_prepare_hourly.py

@@ -0,0 +1,126 @@
+"""Failing tests for scripts/prepare_hourly.py.
+
+Tests the M1 → hourly aggregation logic using synthetic data before touching
+real downloads.
+
+Run: cd ~/dev/AI/jepa-fx-risk && .venv/bin/python -m pytest tests/test_prepare_hourly.py -v
+"""
+import numpy as np
+import pandas as pd
+import pytest
+import importlib.util, sys, os
+
+
+def _import():
+    spec = importlib.util.spec_from_file_location(
+        "prepare_hourly", "scripts/prepare_hourly.py"
+    )
+    mod = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(mod)
+    return mod
+
+
+@pytest.fixture(scope="module")
+def ph():
+    return _import()
+
+
+def _make_m1(n_days: int = 3, price: float = 1.1000, noise: float = 0.0005) -> pd.DataFrame:
+    """Synthetic M1 DataFrame starting 2020-01-06 (Monday), 390 ticks/day."""
+    rng = np.random.default_rng(42)
+    # generate full trading hours: Mon-Fri 00:00-23:59 (FX is 24h weekday)
+    start = pd.Timestamp("2020-01-06 00:00:00")  # Monday
+    periods = n_days * 24 * 60
+    ts = pd.date_range(start, periods=periods, freq="min")
+    # remove weekends
+    ts = ts[ts.day_of_week < 5]
+    prices = price + np.cumsum(rng.normal(0, noise, len(ts)))
+    return pd.DataFrame({"ts": ts, "close": prices})
+
+
+# 1. resample_to_hourly: DataFrame has correct columns
+def test_columns(ph):
+    m1 = _make_m1()
+    hourly = ph.resample_to_hourly(m1)
+    assert set(["datetime", "close", "ret", "realized_vol"]).issubset(hourly.columns), \
+        f"missing columns: {hourly.columns.tolist()}"
+
+
+# 2. No cross-weekend interpolation: gap between Friday 23:xx and Sunday/Monday must remain
+def test_no_weekend_interpolation(ph):
+    # Make 2 days: Friday + Monday (skip Saturday/Sunday)
+    fri = pd.date_range("2020-01-10 00:00", "2020-01-10 23:59", freq="min")  # Friday
+    mon = pd.date_range("2020-01-13 00:00", "2020-01-13 23:59", freq="min")  # Monday
+    ts = fri.append(mon)
+    prices = 1.1 + np.cumsum(np.random.default_rng(0).normal(0, 0.0001, len(ts)))
+    m1 = pd.DataFrame({"ts": ts, "close": prices})
+    hourly = ph.resample_to_hourly(m1)
+    dates = pd.DatetimeIndex(hourly["datetime"]).date
+    import datetime
+    sat = datetime.date(2020, 1, 11)
+    sun = datetime.date(2020, 1, 12)
+    assert sat not in dates and sun not in dates, "weekend rows found in hourly output"
+
+
+# 3. Realized vol = sqrt(sum(r²)) over minute returns in each hour
+def test_realized_vol_formula(ph):
+    # Two hours: anchor gives 10:00 a valid ret; measurement hour has one known log-return.
+    ts0 = pd.date_range("2020-01-06 09:00", periods=60, freq="min")
+    ts1 = pd.date_range("2020-01-06 10:00", periods=60, freq="min")
+    prices0 = np.ones(60) * 1.0
+    # price jumps at minute 1 and STAYS (no reversion) → one non-zero log-return
+    prices1 = np.full(60, np.exp(0.01))
+    prices1[0] = 1.0  # only first tick is at 1.0; jump happens at tick 1
+    m1 = pd.DataFrame({
+        "ts": np.concatenate([ts0, ts1]),
+        "close": np.concatenate([prices0, prices1]),
+    })
+    hourly = ph.resample_to_hourly(m1)
+    assert len(hourly) >= 1, "no rows after resample"
+    rv = hourly.iloc[-1]["realized_vol"]
+    expected = np.sqrt(0.01 ** 2)
+    assert abs(rv - expected) < 1e-6, f"realized_vol={rv:.8f}, expected≈{expected:.8f}"
+
+
+# 4. Only hours with ≥ 30 M1 bars are kept (thin hours dropped)
+def test_thin_hours_dropped(ph):
+    # 4 hours: pre-anchor gives 09:00 a valid ret; full survives; thin (11:00) is dropped.
+    # pre-anchor (08:00): gives 09:00 a valid ret
+    # anchor (09:00): 60 bars, valid ret  → kept
+    # full  (10:00): 60 bars, valid ret  → kept
+    # thin  (11:00): 10 bars             → dropped
+    # Result: 3 hourly candidates, first (pre-anchor) gets NaN ret → dropped → 2 rows
+    pre   = pd.date_range("2020-01-06 08:00", periods=60, freq="min")
+    anchor= pd.date_range("2020-01-06 09:00", periods=60, freq="min")
+    full  = pd.date_range("2020-01-06 10:00", periods=60, freq="min")
+    thin  = pd.date_range("2020-01-06 11:00", periods=10, freq="min")
+    ts = pre.append(anchor).append(full).append(thin)
+    m1 = pd.DataFrame({"ts": ts, "close": np.ones(len(ts)) * 1.1})
+    hourly = ph.resample_to_hourly(m1)
+    assert len(hourly) == 2, f"expected 2 rows (pre-anchor NaN ret dropped + thin dropped), got {len(hourly)}"
+
+
+# 5. Output parquet path and schema (integration — reads actual M1 zips if present)
+def test_output_schema_from_zips(ph, tmp_path):
+    # Build a minimal fake zip structure
+    import zipfile, io
+    # synthetic M1 CSV (histdata format: YYYYMMDD HHMMSS;O;H;L;C;V)
+    rows = []
+    for h in range(24):
+        for m in range(60):
+            rows.append(f"20200106 {h:02d}{m:02d}00;1.10000;1.10100;1.09900;1.10000;100")
+    csv_content = "\n".join(rows).encode()
+    zip_buf = io.BytesIO()
+    with zipfile.ZipFile(zip_buf, "w") as zf:
+        zf.writestr("DAT_ASCII_EURUSD_M1_2020.csv", csv_content)
+    zip_buf.seek(0)
+    raw_dir = tmp_path / "raw"
+    raw_dir.mkdir()
+    (raw_dir / "DAT_ASCII_EURUSD_M1_2020.zip").write_bytes(zip_buf.read())
+
+    out_path = str(tmp_path / "eurusd_hourly.parquet")
+    ph.build_hourly_parquet(raw_dir=str(raw_dir), out_path=out_path)
+    assert os.path.exists(out_path), "output parquet not created"
+    df = pd.read_parquet(out_path)
+    assert set(["datetime", "close", "ret", "realized_vol"]).issubset(df.columns)
+    assert len(df) > 0

train.py

@@ -0,0 +1,262 @@
+"""train.py — autoresearch agent file (only this may be edited).
+
+HEPA backbone (Petersen et al., arXiv:2605.11130, ICML 2026 Spotlight):
+Causal Transformer pre-trained via horizon-conditioned JEPA. Predictor
+maps (h_t, Δt) → predicted future embedding; loss = VICReg (L1 alignment
+on L2-normalised reps + variance-covariance regulariser, no stop-gradient).
+Probe: ridge regression on the last-token embedding (true OOS split).
+
+Agent may tune: encoder depth/width, patch geometry, ALPHA, DELTA_T_MAX,
+optimizer, LR. Do NOT touch prepare_data.py, loop.py, or the data pipeline.
+"""
+import json
+import math
+import numpy as np
+import pandas as pd
+import torch
+import torch.nn as nn
+import torch.nn.functional as F
+
+# --- agent-tunable knobs ---
+USE_HOURLY  = True   # prefer eurusd_hourly.parquet when available
+WINDOW      = 240    # hourly: 10 trading days; if USE_HOURLY=False reset to 60
+PATCH_LEN   = 24     # hourly: 1-day patches (10 tokens); if USE_HOURLY=False reset to 10
+D_MODEL     = 128
+DEPTH       = 2
+N_HEADS     = 4
+ALPHA       = 0.1   # VICReg mixing weight (fixed at 0.1 in HEPA paper)
+DELTA_T_MAX = 3     # max prediction horizon in patches (1..min(DELTA_T_MAX, N-1-c))
+BATCH_SIZE  = 512   # mini-batch per step (hourly dataset is too large for full-batch)
+EPOCHS      = 300
+LR          = 3e-4
+SEED        = 0
+# ---------------------------
+
+torch.manual_seed(SEED)
+np.random.seed(SEED)
+dev = "cuda" if torch.cuda.is_available() else "cpu"
+
+
+# ── VICReg pretraining loss ──────────────────────────────────────────────────
+
+def vicreg_loss(h_pred: torch.Tensor, h_target: torch.Tensor, alpha: float = 0.1) -> torch.Tensor:
+    """L = (1-α)·L1(normalize(ĥ), normalize(h*)) + α·(L_var + L_cov).
+
+    Both encoders receive gradients (joint training — no stop-grad on h_target).
+    Variance-covariance terms prevent embedding collapse.
+    """
+    pred_n = F.normalize(h_pred, dim=-1)
+    targ_n = F.normalize(h_target, dim=-1)
+    l1 = F.l1_loss(pred_n, targ_n)
+    # variance hinge: push each feature std toward ≥ 1
+    std = h_pred.std(dim=0) + 1e-4
+    l_var = F.relu(1.0 - std).mean()
+    # covariance penalty: decorrelate features
+    B, D = h_pred.shape
+    h_c = h_pred - h_pred.mean(dim=0, keepdim=True)
+    cov = (h_c.t() @ h_c) / max(B - 1, 1)
+    off = cov - torch.diag(torch.diag(cov))
+    l_cov = (off ** 2).sum() / D
+    return (1 - alpha) * l1 + alpha * (l_var + l_cov)
+
+
+# ── CausalEncoder ─────────────────────────────────────────────────────────────
+
+class CausalEncoder(nn.Module):
+    """Non-overlapping patches → per-patch LayerNorm → causal Transformer → all tokens (B, N, D).
+
+    Per-patch LayerNorm instead of full-window RevIN: each patch is normalised
+    using only its own timesteps, so no future statistics leak into past tokens.
+    Use [:, -1, :] for probing (last token sees full context).
+    Use [:, c, :] for JEPA pretraining (context-at-c).
+    """
+    def __init__(self, n_channels: int, patch_len: int, d_model: int,
+                 n_heads: int, depth: int):
+        super().__init__()
+        self.patch_len = patch_len
+        self.d_model   = d_model
+        patch_dim = patch_len * n_channels
+        self.patch_norm = nn.LayerNorm(patch_dim)   # applied per-patch, no future leakage
+        self.embed     = nn.Linear(patch_dim, d_model)
+        layer = nn.TransformerEncoderLayer(d_model, n_heads, 2 * d_model,
+                                           dropout=0.0, batch_first=True)
+        self.tf   = nn.TransformerEncoder(layer, num_layers=depth)
+        self.norm = nn.LayerNorm(d_model)
+
+    def forward(self, x: torch.Tensor) -> torch.Tensor:
+        B, W, F = x.shape
+        P = self.patch_len
+        N = W // P
+        tokens = x[:, :N * P, :].reshape(B, N, P * F)
+        tokens = self.embed(self.patch_norm(tokens))
+        # sinusoidal PE
+        pos = torch.arange(N, device=x.device).float()
+        div = torch.exp(torch.arange(0, self.d_model, 2, device=x.device).float()
+                        * -(math.log(10000.0) / self.d_model))
+        pe = torch.zeros(N, self.d_model, device=x.device)
+        pe[:, 0::2] = torch.sin(pos.unsqueeze(1) * div)
+        pe[:, 1::2] = torch.cos(pos.unsqueeze(1) * div)
+        tokens = tokens + pe
+        # causal mask
+        mask = nn.Transformer.generate_square_subsequent_mask(N, device=x.device)
+        return self.norm(self.tf(tokens, mask=mask, is_causal=True))
+
+
+# ── HorizonPredictor ─────────────────────────────────────────────────────────
+
+class HorizonPredictor(nn.Module):
+    """MLP(cat(h_t, Δt)) → predicted future embedding."""
+    def __init__(self, d_model: int):
+        super().__init__()
+        self.net = nn.Sequential(
+            nn.Linear(d_model + 1, d_model), nn.GELU(),
+            nn.Linear(d_model, d_model),     nn.GELU(),
+            nn.Linear(d_model, d_model),
+        )
+
+    def forward(self, h: torch.Tensor, delta_t: torch.Tensor) -> torch.Tensor:
+        dt = delta_t.float().unsqueeze(-1)
+        return self.net(torch.cat([h, dt], dim=-1))
+
+
+# ── Data ─────────────────────────────────────────────────────────────────────
+
+def build():
+    """Year-based split: encoder trains on ≤2021; probe evaluates on ≥2022 OOS.
+
+    Uses eurusd_hourly.parquet when USE_HOURLY=True and the file exists;
+    falls back to eurusd_daily.parquet otherwise.
+    """
+    import os
+    hourly_path = "data/processed/eurusd_hourly.parquet"
+    daily_path  = "data/processed/eurusd_daily.parquet"
+    if USE_HOURLY and os.path.exists(hourly_path):
+        df = pd.read_parquet(hourly_path).reset_index(drop=True)
+        df["date"] = pd.to_datetime(df["datetime"])
+    else:
+        df = pd.read_parquet(daily_path).reset_index(drop=True)
+        df["date"] = pd.to_datetime(df["date"])
+    feats  = df[["ret", "realized_vol"]].to_numpy(np.float32)
+    target = df["realized_vol"].to_numpy(np.float32)
+    tr_idx = df.index[df["date"].dt.year <= 2021].tolist()
+    te_idx = df.index[df["date"].dt.year >= 2022].tolist()
+    mu = feats[:tr_idx[-1]+1].mean(0)
+    sd = feats[:tr_idx[-1]+1].std(0) + 1e-8
+    fn = (feats - mu) / sd
+    def windows(idx):
+        X, y = [], []
+        for t in idx:
+            if t - WINDOW >= 0 and t + 1 < len(df):
+                X.append(fn[t - WINDOW:t]); y.append(target[t + 1])
+        return np.stack(X).astype(np.float32), np.array(y, np.float32)
+    return windows(tr_idx), windows(te_idx)
+
+
+# ── Training ──────────────────────────────────────────────────────────────────
+
+def main():
+    (Xtr, ytr), (Xte, yte) = build()
+    n_feats   = Xtr.shape[2]
+    n_patches = WINDOW // PATCH_LEN
+    N_tr      = len(Xtr)
+    bs        = min(BATCH_SIZE, N_tr)
+
+    enc  = CausalEncoder(n_feats, PATCH_LEN, D_MODEL, N_HEADS, DEPTH).to(dev)
+    pred = HorizonPredictor(D_MODEL).to(dev)
+    opt  = torch.optim.AdamW(list(enc.parameters()) + list(pred.parameters()), lr=LR)
+
+    for ep in range(EPOCHS):
+        # Random mini-batch (avoids OOM on large hourly dataset)
+        idx_b = torch.randperm(N_tr)[:bs]
+        Xb    = torch.tensor(Xtr[idx_b.numpy()], device=dev)
+
+        # Sample random context position and horizon
+        c  = torch.randint(0, n_patches - 1, ()).item()
+        dt = torch.randint(1, max(2, min(DELTA_T_MAX, n_patches - 1 - c) + 1), ()).item()
+
+        tokens = enc(Xb)                        # (bs, N, D)
+        h_ctx  = tokens[:, c, :]               # context embedding
+        h_tgt  = tokens[:, c + dt, :]          # target embedding (joint training)
+        h_hat  = pred(h_ctx, torch.full((bs,), float(dt), device=dev))
+        loss   = vicreg_loss(h_hat, h_tgt, alpha=ALPHA)
+        opt.zero_grad(); loss.backward(); opt.step()
+
+    enc.eval()
+    with torch.no_grad():
+        def embed(X_np):
+            chunks = []
+            for i in range(0, len(X_np), bs):
+                t = torch.tensor(X_np[i:i+bs], device=dev)
+                chunks.append(enc(t)[:, -1, :].cpu().numpy())
+            return np.concatenate(chunks, axis=0)
+
+        Etr = embed(Xtr)
+        Ete = embed(Xte)
+
+    # Ridge probe: fit on train, evaluate on OOS (true OOS R²)
+    mu_e = Etr.mean(0); sd_e = Etr.std(0) + 1e-8
+    Etr_n = (Etr - mu_e) / sd_e
+    Ete_n = (Ete - mu_e) / sd_e
+    A  = np.hstack([Etr_n, np.ones((len(Etr_n), 1))])
+    w  = np.linalg.solve(A.T @ A + 1e-3 * np.eye(A.shape[1]), A.T @ ytr)
+    pred_np = np.hstack([Ete_n, np.ones((len(Ete_n), 1))]) @ w
+    ss_res  = ((yte - pred_np) ** 2).sum()
+    ss_tot  = ((yte - yte.mean()) ** 2).sum()
+    val_vol_r2 = float(1 - ss_res / ss_tot)
+
+    json.dump({
+        "val_vol_r2": val_vol_r2, "n_test": len(yte),
+        "knobs": {"WINDOW": WINDOW, "PATCH_LEN": PATCH_LEN,
+                  "D_MODEL": D_MODEL, "DEPTH": DEPTH, "ALPHA": ALPHA,
+                  "DELTA_T_MAX": DELTA_T_MAX, "EPOCHS": EPOCHS},
+    }, open("metrics.json", "w"), indent=2)
+    print("val_vol_r2 = %.4f  (n_test=%d, dev=%s)" % (val_vol_r2, len(yte), dev))
+
+    # ── EXPORT BLOCK — do NOT edit (agent boundary) ──────────────────────────
+    # Set EXPORT_EMBEDDINGS=1 to write embeddings.json for the Go eval harness.
+    import os
+    if os.environ.get("EXPORT_EMBEDDINGS") == "1":
+        hourly_path2 = "data/processed/eurusd_hourly.parquet"
+        daily_path2  = "data/processed/eurusd_daily.parquet"
+        if USE_HOURLY and os.path.exists(hourly_path2):
+            df2 = pd.read_parquet(hourly_path2).reset_index(drop=True)
+            df2["date"] = pd.to_datetime(df2["datetime"])
+        else:
+            df2 = pd.read_parquet(daily_path2).reset_index(drop=True)
+            df2["date"] = pd.to_datetime(df2["date"])
+        tr_mask = df2["date"].dt.year <= 2021
+        feats2 = df2[["ret", "realized_vol"]].to_numpy(np.float32)
+        mu2 = feats2[tr_mask].mean(0); sd2 = feats2[tr_mask].std(0) + 1e-8
+        fn2 = (feats2 - mu2) / sd2
+        def _export_windows(year_mask):
+            idx = df2.index[year_mask].tolist()
+            Xs, dates, rvs = [], [], []
+            for t in idx:
+                if t - WINDOW >= 0 and t + 1 < len(df2):
+                    Xs.append(fn2[t - WINDOW:t])
+                    dates.append(str(df2["date"].iloc[t].date()))
+                    rvs.append(float(df2["realized_vol"].iloc[t + 1]))
+            if not Xs:
+                return [], [], []
+            Xa = np.stack(Xs)
+            chunks = []
+            with torch.no_grad():
+                for i in range(0, len(Xa), bs):
+                    chunks.append(enc(torch.tensor(Xa[i:i+bs], device=dev))[:, -1, :].cpu().numpy())
+            E = np.concatenate(chunks, axis=0).tolist()
+            return E, dates, rvs
+        Etr2, dates_tr, rv_tr = _export_windows(tr_mask)
+        Eoos, dates_oos, rv_oos = _export_windows(df2["date"].dt.year >= 2022)
+        hv_thr = float(np.percentile(rv_oos, 67))
+        hv_label = [1 if v >= hv_thr else 0 for v in rv_oos]
+        json.dump({"embeddings": Eoos, "dates": dates_oos,
+                   "realized_vol": rv_oos, "hv_label": hv_label,
+                   "train_embeddings": Etr2, "train_realized_vol": rv_tr},
+                  open("embeddings.json", "w"))
+        print("exported embeddings.json  train=%d  oos=%d  HV=%d/%d" % (
+            len(Etr2), len(Eoos), sum(hv_label), len(hv_label)))
+    # ── END EXPORT BLOCK ─────────────────────────────────────────────────────
+
+
+if __name__ == "__main__":
+    main()