Files
foxhunt/crates/ml-alpha/examples/phase1d_long_horizon.rs
jgrusewski e190ecfa61 feat(ml-alpha): backtest cost sweep + annualised Sharpe (Phase 1d.4)
Initial single-cost backtest at τ=0.25 cost=0.25 revealed the binding
constraint: mean_ret = -0.176, pre-cost EV ≈ +0.074, so the model has
a real directional edge but K=6000 price moves are too small to clear
1-tick round-trip cost. The cost-vs-edge balance is the real Phase 1d.4
verdict question, not whether the model has signal.

Two enhancements per the insight block in the previous run:

1. **Cost sweep**: GpuBacktest::run now takes `costs: &[f32]` and
   produces (C × T) rows instead of T. The smoke runs at five costs:
     - 0.0     frictionless upper bound (theoretical max Sharpe)
     - 0.0625  quarter-tick (very aggressive execution)
     - 0.125   half-tick (professional desk)
     - 0.25    1 tick = $12.50/contract (retail / pessimistic)
     - 0.50    2 ticks (very pessimistic)
   Tells us the break-even cost where Sharpe crosses zero.

2. **Annualised Sharpe**: per-trade Sharpe × sqrt(trades_per_year).
   trades_per_year = n_trades × (seconds_per_year / test_time_span_seconds).
   test_time_span_seconds derived from first/last test sequence end-bar
   timestamps via FxCacheReader::record_timestamp. Standard Sharpe-time-
   scaling assumption (trades roughly i.i.d.); imperfect when signals
   cluster in correlated regimes, but the right ballpark for comparison
   with industry benchmarks.

Output adds per-cost-band "best operating point" tables plus a clear
"REALISTIC VERDICT" line at cost=0.125 (half-tick — what a professional
desk would actually pay) with three gates:
- Sharpe_ann > 2.0 → deployable
- 0.5 < Sharpe_ann ≤ 2.0 → marginal
- Sharpe_ann ≤ 0.5 → fail at realistic cost

The "FRICTIONLESS UPPER BOUND" line reports the intrinsic edge — what
the model could theoretically achieve at zero cost. Even if realistic
Sharpe fails, this number tells us whether the model has anything to
optimise toward at deployment.

Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>
2026-05-15 09:44:55 +02:00

637 lines
32 KiB
Rust

//! Phase 1d.2 — Multi-minute label smoke (K=6000 snapshots ≈ 1-5 min forward).
//!
//! THE decisive gate for the two-head architecture. The K-sweep (Phase 1c,
//! commit `db874b184`) showed alpha at stateless single-snapshot resolution
//! decays from K=50 peak to gone by K=500. Mamba2's job is to amplify
//! short-horizon evidence into a long-horizon prediction via SSM state
//! accumulation across `seq_len` snapshots.
//!
//! This test: feed a trailing window of 32 snapshots into Mamba2; predict
//! the binary direction `sign(mid[t+6000] - mid[t])`. Mamba's state must
//! integrate the per-step microstructure signal into a multi-minute call.
//!
//! Gate (per implementation plan):
//! - AUC > 0.55 at K=6000 → multi-minute alpha confirmed; design works.
//! - AUC < 0.52 → DECISIVE FAIL. The two-head architecture cannot turn
//! short-window snapshot context into multi-minute prediction with
//! the current model; need different inputs, model class, or both.
//!
//! Why this is decisive: if we *can* predict long horizons from short
//! tick context plus sequence state, the existing 81-dim snapshot stack
//! is the right foundation and the rest of FoxhuntQ-Δ proceeds. If we
//! can't, the whole architecture is wrong and needs a redesign before
//! any production work.
use anyhow::{Context, Result};
use clap::Parser;
use cudarc::driver::CudaContext;
use rand::{Rng, SeedableRng};
use rand_chacha::ChaCha8Rng;
use std::sync::Arc;
use tracing_subscriber::EnvFilter;
use ml_alpha::calibration::{Calibrator, IsotonicCalibrator, PlattScaler};
use ml_alpha::eval::{accuracy_from_logits, auc_from_logits};
use ml_alpha::fxcache_reader::{FxCacheReader, COL_RAW_CLOSE, FEAT_DIM};
use ml_alpha::mamba2_block::{
Mamba2AdamW, Mamba2AdamWConfig, Mamba2Block, Mamba2BlockConfig,
};
use ml_alpha::metrics_detail::{brier_score, log_loss, stratified_accuracy};
use ml_alpha::mlp::{MlpConfig, MlpModel};
use ml_alpha::backtest::GpuBacktest;
use ml_alpha::multi_horizon_labels::generate_labels;
use ml_core::cuda_autograd::gpu_tensor::GpuTensor;
#[derive(Parser, Debug)]
#[command(name = "phase1d_long_horizon", about = "FoxhuntQ-Δ Phase 1d.2 multi-minute Mamba smoke")]
struct Cli {
#[arg(long)] fxcache_path: String,
/// Label horizon in snapshots. K=6000 ≈ 1-5 min forward at typical rates.
#[arg(long, default_value_t = 6000)] horizon: usize,
/// Sequence length per training example (≤ kernel cap of 32).
#[arg(long, default_value_t = 32)] seq_len: usize,
#[arg(long, default_value_t = 64)] hidden_dim: usize,
#[arg(long, default_value_t = 16)] state_dim: usize,
#[arg(long, default_value_t = 5)] epochs: usize,
#[arg(long, default_value_t = 128)] batch_size: usize,
#[arg(long, default_value_t = 3e-3)] lr: f32,
/// 80/20 train/val split (purged via the embargo).
#[arg(long, default_value_t = 0.8)] train_frac: f32,
/// Embargo bars between train and val ranges (default = horizon for safety).
#[arg(long)] embargo: Option<usize>,
/// Subsample stride for train sequences (1 = use every starting position).
#[arg(long, default_value_t = 4)] train_stride: usize,
#[arg(long, default_value_t = 42)] seed: u64,
/// Fraction of val to use as calibration set (rest is held-out test).
/// 0 disables Platt/isotonic post-hoc calibration.
#[arg(long, default_value_t = 0.5)] cal_frac: f32,
/// Stacker hidden dim (small — input is 7-dim, plenty of capacity).
#[arg(long, default_value_t = 32)] stacker_hidden: usize,
/// Stacker training epochs.
#[arg(long, default_value_t = 20)] stacker_epochs: usize,
/// Stacker learning rate.
#[arg(long, default_value_t = 1e-2)] stacker_lr: f32,
/// Stacker batch size.
#[arg(long, default_value_t = 1024)] stacker_batch: usize,
/// Round-trip transaction cost in price units (ES.FUT: 0.25 = 1 tick = $12.50/contract).
#[arg(long, default_value_t = 0.25)] cost_per_trade: f32,
}
fn main() -> Result<()> {
tracing_subscriber::fmt()
.with_env_filter(EnvFilter::try_from_default_env().unwrap_or_else(|_| EnvFilter::new("warn")))
.init();
let cli = Cli::parse();
let ctx = CudaContext::new(0).context("init CUDA")?;
let stream = ctx.default_stream();
// ── Load fxcache + extract feature matrix and prices ──────────────
let reader = FxCacheReader::open(&cli.fxcache_path)?;
let alpha_dim = reader.alpha_feature_dim()
.ok_or_else(|| anyhow::anyhow!("fxcache lacks alpha column"))?;
let n_bars = reader.bar_count();
println!("fxcache: {} rows, alpha_dim={}, horizon K={}", n_bars, alpha_dim, cli.horizon);
// Flat feature matrix (CPU-resident) and prices.
let mut feature_matrix: Vec<f32> = Vec::with_capacity(n_bars * alpha_dim);
let mut prices: Vec<f32> = Vec::with_capacity(n_bars);
for i in 0..n_bars {
let row = reader.alpha_features(i).context("alpha row missing")?;
feature_matrix.extend_from_slice(row);
let rec = reader.record(i);
prices.push(rec.targets[COL_RAW_CLOSE - FEAT_DIM]);
}
// ── Generate K=6000 labels ────────────────────────────────────────
let labels = generate_labels(&prices, cli.horizon);
println!(
"labels: kept={} dropped_edge={} dropped_invalid={} up_frac={:.4}",
labels.labels.len(), labels.n_dropped_edge, labels.n_dropped_invalid,
labels.labels.iter().sum::<f32>() / labels.labels.len().max(1) as f32
);
if labels.labels.is_empty() {
anyhow::bail!("no valid labels at K={} — fxcache too small or all-tied", cli.horizon);
}
// ── Purged train/val split over labels.valid_indices ──────────────
// We split the VALID label list (not the raw bar list) at train_frac,
// then apply an embargo so train sequences ending near the boundary
// don't share forward-window prices with val.
let embargo = cli.embargo.unwrap_or(cli.horizon);
let n_valid = labels.valid_indices.len();
let n_train_target = (n_valid as f32 * cli.train_frac) as usize;
let train_split_bar = labels.valid_indices[n_train_target - 1];
let val_start_bar = train_split_bar + embargo;
let train_label_range: Vec<usize> = (0..n_valid)
.filter(|&p| labels.valid_indices[p] < train_split_bar.saturating_sub(embargo))
.collect();
let val_label_range: Vec<usize> = (0..n_valid)
.filter(|&p| labels.valid_indices[p] >= val_start_bar)
.collect();
println!(
"split: n_train_labels={} n_val_labels={} (embargo={} bars, val_start_bar={})",
train_label_range.len(), val_label_range.len(), embargo, val_start_bar
);
// Sequence start: bar `valid_indices[p] - seq_len + 1` to bar `valid_indices[p]`.
// Filter out any p where the sequence would underflow.
let make_seq_label_pos = |label_positions: &[usize], stride: usize| -> Vec<usize> {
label_positions.iter()
.filter(|&&p| labels.valid_indices[p] + 1 >= cli.seq_len)
.copied()
.step_by(stride)
.collect()
};
let train_pos = make_seq_label_pos(&train_label_range, cli.train_stride);
let val_pos = make_seq_label_pos(&val_label_range, 1);
println!("train sequences (stride={}): {} val sequences: {}",
cli.train_stride, train_pos.len(), val_pos.len());
// ── Build Mamba2 + AdamW ──────────────────────────────────────────
let block_cfg = Mamba2BlockConfig {
in_dim: alpha_dim,
hidden_dim: cli.hidden_dim,
state_dim: cli.state_dim,
seq_len: cli.seq_len,
};
let mut block = Mamba2Block::new(block_cfg.clone(), Arc::clone(&stream))?;
let mut opt = Mamba2AdamW::new(&block, Mamba2AdamWConfig {
lr: cli.lr, ..Default::default()
})?;
println!("Mamba2Block params: {}", block.param_count());
let gather = |sel: &[usize]| -> Result<(GpuTensor, Vec<f32>)> {
let b = sel.len();
let mut host = Vec::with_capacity(b * cli.seq_len * alpha_dim);
let mut ys = Vec::with_capacity(b);
for &p in sel {
let end_bar = labels.valid_indices[p];
let start_bar = end_bar + 1 - cli.seq_len;
for t in 0..cli.seq_len {
let off = (start_bar + t) * alpha_dim;
host.extend_from_slice(&feature_matrix[off..off + alpha_dim]);
}
ys.push(labels.labels[p]);
}
let dev = stream.clone_htod(&host)?;
let tensor = GpuTensor::new(dev, vec![b, cli.seq_len, alpha_dim])
.map_err(|e| anyhow::anyhow!("batch tensor: {e}"))?;
Ok((tensor, ys))
};
let bce_loss = |logits: &[f32], ys: &[f32]| -> f32 {
let mut s = 0.0_f32;
let eps = 1e-7_f32;
for (&z, &y) in logits.iter().zip(ys.iter()) {
let z = z.clamp(-50.0, 50.0);
let p = (1.0 / (1.0 + (-z).exp())).clamp(eps, 1.0 - eps);
s += -(y * p.ln() + (1.0 - y) * (1.0 - p).ln());
}
s / logits.len() as f32
};
// ── Training loop ─────────────────────────────────────────────────
let mut rng = ChaCha8Rng::seed_from_u64(cli.seed);
let train_len = train_pos.len();
if train_len < cli.batch_size {
anyhow::bail!("n_train ({train_len}) < batch_size ({})", cli.batch_size);
}
let batches_per_epoch = train_len / cli.batch_size;
for epoch in 0..cli.epochs {
let mut perm: Vec<usize> = (0..train_len).collect();
for i in (1..train_len).rev() {
let j = rng.gen_range(0..=i);
perm.swap(i, j);
}
let mut loss_sum = 0.0_f32;
for batch_idx in 0..batches_per_epoch {
let raw_sel = &perm[batch_idx * cli.batch_size .. (batch_idx + 1) * cli.batch_size];
let sel: Vec<usize> = raw_sel.iter().map(|&i| train_pos[i]).collect();
let (input, ys) = gather(&sel)?;
let (logit, cache) = block.forward_train(&input)?;
let logit_host = logit.to_host(&stream)?;
loss_sum += bce_loss(&logit_host, &ys);
let n_b = ys.len() as f32;
let d_logit_host: Vec<f32> = logit_host.iter().zip(ys.iter())
.map(|(&z, &y)| {
let p = 1.0 / (1.0 + (-z.clamp(-50.0, 50.0)).exp());
(p - y) / n_b
})
.collect();
let d_logit_dev = stream.clone_htod(&d_logit_host)?;
let d_logit = GpuTensor::new(d_logit_dev, vec![ys.len(), 1])?;
let grads = block.backward(&cache, &d_logit)?;
opt.step(&mut block, &grads)?;
}
println!("epoch {epoch:2} mean_train_bce={:.5} n_batches={batches_per_epoch}",
loss_sum / batches_per_epoch as f32);
}
// ── Validation ────────────────────────────────────────────────────
let mut val_logits: Vec<f32> = Vec::with_capacity(val_pos.len());
let mut val_ys: Vec<f32> = Vec::with_capacity(val_pos.len());
let mut i = 0;
while i < val_pos.len() {
let this = cli.batch_size.min(val_pos.len() - i);
let sel: Vec<usize> = val_pos[i..i + this].to_vec();
let (input, ys) = gather(&sel)?;
let logit = block.forward(&input)?;
let logit_host = logit.to_host(&stream)?;
val_logits.extend_from_slice(&logit_host);
val_ys.extend_from_slice(&ys);
i += this;
}
let labels_u8: Vec<u8> = val_ys.iter().map(|&y| if y > 0.5 { 1 } else { 0 }).collect();
let acc = accuracy_from_logits(&val_logits, &labels_u8);
let auc = auc_from_logits(&val_logits, &labels_u8);
let uncal_brier = brier_score(&val_logits, &val_ys);
let uncal_logl = log_loss(&val_logits, &val_ys);
let up_frac = val_ys.iter().sum::<f32>() / val_ys.len() as f32;
let chance_brier = up_frac * (1.0 - up_frac);
let chance_logl = -(up_frac.ln() * up_frac + (1.0 - up_frac).ln() * (1.0 - up_frac));
println!();
println!("===========================================================");
println!("PHASE 1d.2 — MAMBA2 MULTI-MINUTE SMOKE (K={})", cli.horizon);
println!("===========================================================");
println!("val sequences: {}", val_logits.len());
println!("up_fraction: {:.4} (chance baselines: Brier={:.5} log-loss={:.5})",
up_frac, chance_brier, chance_logl);
println!("--- Uncalibrated ---");
println!(" accuracy: {:.4}", acc);
println!(" AUC: {:.4}", auc);
println!(" Brier: {:.5}", uncal_brier);
println!(" log-loss: {:.5}", uncal_logl);
// ── Post-hoc Platt + Isotonic calibration on a held-out cal split ──
if cli.cal_frac > 0.0 && cli.cal_frac < 1.0 {
let n_val = val_logits.len();
let n_cal = (n_val as f32 * cli.cal_frac) as usize;
let (cal_l, test_l) = val_logits.split_at(n_cal);
let (cal_y, test_y) = val_ys.split_at(n_cal);
let test_labels_u8: Vec<u8> = test_y.iter().map(|&y| if y > 0.5 { 1 } else { 0 }).collect();
// Headline on test half BEFORE calibration (baseline for comparison).
let test_acc_raw = accuracy_from_logits(test_l, &test_labels_u8);
let test_auc_raw = auc_from_logits(test_l, &test_labels_u8);
let test_brier_raw = brier_score(test_l, test_y);
let test_logl_raw = log_loss(test_l, test_y);
// Fit + apply Platt.
let platt = PlattScaler::fit(cal_l, cal_y, 2000, 1e-2)
.map_err(|e| anyhow::anyhow!("platt fit: {e}"))?;
let platt_probs = platt.transform(test_l);
// accuracy_from_logits expects logits; convert calibrated probs back to logits.
let to_logits = |probs: &[f32]| -> Vec<f32> {
probs.iter().map(|&p| {
let pc = p.clamp(1e-7, 1.0 - 1e-7);
(pc / (1.0 - pc)).ln()
}).collect()
};
let platt_logits = to_logits(&platt_probs);
let platt_acc = accuracy_from_logits(&platt_logits, &test_labels_u8);
let platt_auc = auc_from_logits(&platt_logits, &test_labels_u8);
let platt_brier = brier_score(&platt_logits, test_y);
let platt_logl = log_loss(&platt_logits, test_y);
// Fit + apply Isotonic.
let iso = IsotonicCalibrator::fit(cal_l, cal_y)
.map_err(|e| anyhow::anyhow!("iso fit: {e}"))?;
let iso_probs = iso.transform(test_l);
let iso_logits = to_logits(&iso_probs);
let iso_acc = accuracy_from_logits(&iso_logits, &test_labels_u8);
let iso_auc = auc_from_logits(&iso_logits, &test_labels_u8);
let iso_brier = brier_score(&iso_logits, test_y);
let iso_logl = log_loss(&iso_logits, test_y);
println!("--- Held-out test (cal_frac={:.2}, n_cal={}, n_test={}) ---",
cli.cal_frac, n_cal, test_l.len());
println!(" Uncalibrated: acc={:.4} AUC={:.4} Brier={:.5} log-loss={:.5}",
test_acc_raw, test_auc_raw, test_brier_raw, test_logl_raw);
println!(" Platt: acc={:.4} AUC={:.4} Brier={:.5} log-loss={:.5} (a={:.4} b={:.4})",
platt_acc, platt_auc, platt_brier, platt_logl, platt.a, platt.b);
println!(" Isotonic: acc={:.4} AUC={:.4} Brier={:.5} log-loss={:.5}",
iso_acc, iso_auc, iso_brier, iso_logl);
// ── Stacked regime head ──────────────────────────────────────
// GPU-native MLP that takes [mamba_logit, 6 Block-S features] as
// input and predicts the same direction label. Trains on the cal
// half, evaluates on the test half. The cal half is used for both
// Platt and the stacker so all comparisons are on the same
// held-out test bars.
{
const N_BLOCK_S: usize = 6;
const STACKER_IN_DIM: usize = 1 + N_BLOCK_S; // mamba_logit + Block-S
// Block-S column offsets within the snapshot row.
const BS_OFFSETS: [usize; N_BLOCK_S] = [75, 76, 77, 78, 79, 80];
// Build a flat [N, STACKER_IN_DIM] feature matrix where N = n_val.
// Each row: [mamba_logit, time_since_trade, time_since_snap,
// book_event_rate, spread_bps, L1_imbalance, micro_mid_drift].
// Block-S values come from the END BAR of each val sequence.
let mut stacker_inputs: Vec<f32> = Vec::with_capacity(val_pos.len() * STACKER_IN_DIM);
for (i, &p) in val_pos.iter().enumerate() {
let end_bar = labels.valid_indices[p];
stacker_inputs.push(val_logits[i]);
for &off in &BS_OFFSETS {
stacker_inputs.push(feature_matrix[end_bar * alpha_dim + off]);
}
}
// Z-score normalise the Block-S columns on the cal half so the
// MLP's Xavier init is on a sensible scale. The mamba_logit is
// already close to standard-normal scale, so we leave it raw.
let mut col_mean = [0.0_f64; STACKER_IN_DIM];
let mut col_var = [0.0_f64; STACKER_IN_DIM];
let n_cal_rows = n_cal;
for r in 0..n_cal_rows {
for c in 1..STACKER_IN_DIM {
col_mean[c] += stacker_inputs[r * STACKER_IN_DIM + c] as f64;
}
}
for c in 1..STACKER_IN_DIM {
col_mean[c] /= n_cal_rows as f64;
}
for r in 0..n_cal_rows {
for c in 1..STACKER_IN_DIM {
let d = stacker_inputs[r * STACKER_IN_DIM + c] as f64 - col_mean[c];
col_var[c] += d * d;
}
}
let mut col_std = [1.0_f32; STACKER_IN_DIM];
for c in 1..STACKER_IN_DIM {
let s = (col_var[c] / n_cal_rows as f64).sqrt().max(1e-6);
col_std[c] = s as f32;
}
// Apply normalisation across the FULL val matrix using cal stats.
for r in 0..val_pos.len() {
for c in 1..STACKER_IN_DIM {
let v = stacker_inputs[r * STACKER_IN_DIM + c];
stacker_inputs[r * STACKER_IN_DIM + c] =
(v - col_mean[c] as f32) / col_std[c];
}
}
// Train/test split on the SAME boundary as Platt (first n_cal rows = train,
// rest = test). Temporal split (val_pos is in ascending end-bar order).
let stacker_train_n = n_cal;
let stacker_test_n = val_pos.len() - stacker_train_n;
// Build and train the stacker MLP.
let stacker_cfg = MlpConfig {
in_dim: STACKER_IN_DIM,
hidden_dim: cli.stacker_hidden,
out_dim: 1,
};
let mut stacker = MlpModel::new(stacker_cfg, Arc::clone(&stream))?;
stacker.set_learning_rate(cli.stacker_lr);
println!("--- Stacker MLP ---");
println!(" in_dim={} hidden={} params={}",
STACKER_IN_DIM, cli.stacker_hidden, stacker.param_count());
// Permutation for shuffling each epoch.
let mut srng = ChaCha8Rng::seed_from_u64(cli.seed.wrapping_add(7));
let batch = cli.stacker_batch.min(stacker_train_n);
let batches_per_epoch = stacker_train_n / batch;
for sep in 0..cli.stacker_epochs {
let mut perm: Vec<usize> = (0..stacker_train_n).collect();
for i in (1..stacker_train_n).rev() {
let j = srng.gen_range(0..=i);
perm.swap(i, j);
}
let mut loss_sum = 0.0_f32;
for bi in 0..batches_per_epoch {
let sel = &perm[bi * batch..(bi + 1) * batch];
let mut feats = Vec::with_capacity(batch * STACKER_IN_DIM);
let mut ys = Vec::with_capacity(batch);
for &r in sel {
let off = r * STACKER_IN_DIM;
feats.extend_from_slice(&stacker_inputs[off..off + STACKER_IN_DIM]);
ys.push(val_ys[r]);
}
let loss = stacker.train_step(&feats, &ys, batch)?;
loss_sum += loss;
}
if sep == 0 || (sep + 1) % 5 == 0 {
println!(" epoch {sep:2} mean_bce={:.5}", loss_sum / batches_per_epoch.max(1) as f32);
}
}
// Evaluate stacker on test half.
let mut stacker_logits: Vec<f32> = Vec::with_capacity(stacker_test_n);
let eval_batch = cli.stacker_batch;
let mut i = stacker_train_n;
while i < val_pos.len() {
let this = eval_batch.min(val_pos.len() - i);
let mut feats = Vec::with_capacity(this * STACKER_IN_DIM);
for r in i..i + this {
let off = r * STACKER_IN_DIM;
feats.extend_from_slice(&stacker_inputs[off..off + STACKER_IN_DIM]);
}
let out = stacker.forward_infer(&feats, this)?;
stacker_logits.extend_from_slice(&out);
i += this;
}
let stacker_test_ys: &[f32] = &val_ys[stacker_train_n..];
let stacker_labels_u8: Vec<u8> =
stacker_test_ys.iter().map(|&y| if y > 0.5 { 1 } else { 0 }).collect();
let s_acc = accuracy_from_logits(&stacker_logits, &stacker_labels_u8);
let s_auc = auc_from_logits(&stacker_logits, &stacker_labels_u8);
let s_brier = brier_score(&stacker_logits, stacker_test_ys);
let s_logl = log_loss(&stacker_logits, stacker_test_ys);
println!(" Stacker: acc={:.4} AUC={:.4} Brier={:.5} log-loss={:.5} (n_test={})",
s_acc, s_auc, s_brier, s_logl, stacker_test_n);
// Stratified accuracy of the STACKER on Block-S features —
// tells us whether the stacker absorbed the regime conditioning
// (uniform accuracy across quintiles) or just learned a sharper
// threshold gate (same Q4-elevated pattern as the raw Mamba).
// ── Phase 1d.4 GPU backtest ──────────────────────────────
// Use the stacker's test-half predictions as the trading
// signal. Apply confidence threshold sweep, compute per-trade
// PnL and Sharpe on GPU.
//
// Trade rule: if `|stacker_prob - 0.5| > τ`, take
// direction = sign(prob - 0.5), enter at end_bar, exit at
// end_bar + horizon. Cost in price units per round-trip.
//
// Test sequences span the SECOND HALF of val_pos; for each,
// we need price[end_bar] and price[end_bar + horizon].
let stacker_probs: Vec<f32> = stacker_logits.iter()
.map(|&z| 1.0_f32 / (1.0 + (-z.clamp(-50.0, 50.0)).exp()))
.collect();
let mut prices_t_host: Vec<f32> = Vec::with_capacity(stacker_test_n);
let mut prices_kt_host: Vec<f32> = Vec::with_capacity(stacker_test_n);
for r in stacker_train_n..val_pos.len() {
let end_bar = labels.valid_indices[val_pos[r]];
let kt_bar = end_bar + cli.horizon;
prices_t_host.push(prices[end_bar]);
prices_kt_host.push(prices[kt_bar]);
}
let probs_dev = stream.clone_htod(&stacker_probs)?;
let prices_t_dev = stream.clone_htod(&prices_t_host)?;
let prices_kt_dev = stream.clone_htod(&prices_kt_host)?;
let bt = GpuBacktest::from_block(&block)?;
// Compute the test-set wall-clock time span for annualisation.
// Use the first and last end-bar timestamps from the val sequences.
// FxCacheReader stores timestamps in nanoseconds.
let first_test_end_bar =
labels.valid_indices[val_pos[stacker_train_n]];
let last_test_end_bar =
labels.valid_indices[val_pos[val_pos.len() - 1]];
let first_ts_ns = reader.record_timestamp(first_test_end_bar);
let last_ts_ns = reader.record_timestamp(last_test_end_bar);
let test_time_span_s = ((last_ts_ns - first_ts_ns) as f64) * 1e-9;
println!();
println!("--- GPU backtest sweep (n_test={}, span={:.2}s ≈ {:.2}h) ---",
stacker_test_n, test_time_span_s, test_time_span_s / 3600.0);
// Cost sweep: frictionless / quarter-tick / half-tick / 1-tick / 2-tick.
// 0.25 = 1 tick = $12.50/contract round-trip on ES.FUT.
let thresholds: Vec<f32> = vec![0.00, 0.05, 0.10, 0.15, 0.20, 0.25];
let costs: Vec<f32> = vec![0.0, 0.0625, 0.125, 0.25, 0.50];
let bt_stats = bt.run(
&probs_dev, &prices_t_dev, &prices_kt_dev,
&thresholds, &costs, test_time_span_s,
)?;
println!(" {:>6} {:>5} {:>9} {:>10} {:>10} {:>9} {:>9} {:>10} {:>11} {:>11}",
"cost", "τ", "n_trades", "mean_ret", "std_ret",
"Sharpe", "hit_rate", "trades/yr", "Sharpe_ann", "total_pnl");
for s in &bt_stats {
println!(" {:>6.4} {:>5.2} {:>9} {:>10.5} {:>10.5} {:>9.4} {:>9.4} {:>10.0} {:>11.4} {:>11.2}",
s.cost, s.threshold, s.n_trades, s.mean_ret, s.std_ret,
s.sharpe_per_trade, s.hit_rate, s.trades_per_year,
s.sharpe_annualised, s.total_pnl);
}
// Identify best annualised Sharpe per cost band (≥ 100 trades for stat power).
println!();
println!("BEST OPERATING POINT per cost band (annualised-Sharpe-maximising, n_trades ≥ 100):");
println!(" {:>6} {:>5} {:>9} {:>11} {:>9} {:>11}",
"cost", "τ", "n_trades", "Sharpe_ann", "hit_rate", "total_pnl");
for &c in &costs {
let best = bt_stats.iter()
.filter(|s| (s.cost - c).abs() < 1e-6 && s.n_trades >= 100)
.max_by(|a, b| a.sharpe_annualised
.partial_cmp(&b.sharpe_annualised)
.unwrap_or(std::cmp::Ordering::Equal));
if let Some(b) = best {
println!(" {:>6.4} {:>5.2} {:>9} {:>11.4} {:>9.4} {:>11.2}",
b.cost, b.threshold, b.n_trades,
b.sharpe_annualised, b.hit_rate, b.total_pnl);
}
}
// Overall verdict at the most realistic professional cost: 0.125 (half-tick).
let realistic_best = bt_stats.iter()
.filter(|s| (s.cost - 0.125).abs() < 1e-6 && s.n_trades >= 100)
.max_by(|a, b| a.sharpe_annualised
.partial_cmp(&b.sharpe_annualised)
.unwrap_or(std::cmp::Ordering::Equal));
println!();
if let Some(best) = realistic_best {
println!("REALISTIC VERDICT (cost=0.125, half-tick — professional execution):");
println!(" threshold={:.2} n_trades={} Sharpe_ann={:.4} hit_rate={:.4}",
best.threshold, best.n_trades, best.sharpe_annualised, best.hit_rate);
if best.sharpe_annualised > 2.0 {
println!(" BACKTEST GATE PASS: annualised Sharpe > 2.0 at half-tick cost — deployable.");
} else if best.sharpe_annualised > 0.5 {
println!(" BACKTEST MARGINAL: 0.5 < Sharpe_ann ≤ 2.0 — improve execution or thresholds.");
} else {
println!(" BACKTEST FAIL at realistic cost: signal exists but doesn't survive half-tick friction.");
}
}
// Frictionless upper bound as a diagnostic.
let frictionless = bt_stats.iter()
.filter(|s| (s.cost - 0.0).abs() < 1e-6 && s.n_trades >= 100)
.max_by(|a, b| a.sharpe_annualised
.partial_cmp(&b.sharpe_annualised)
.unwrap_or(std::cmp::Ordering::Equal));
if let Some(f) = frictionless {
println!("FRICTIONLESS UPPER BOUND (cost=0.0): Sharpe_ann={:.4} at τ={:.2} ({} trades) — model's intrinsic edge",
f.sharpe_annualised, f.threshold, f.n_trades);
}
println!();
println!("--- Stacker stratified accuracy (test half, by Block-S feature) ---");
const BLOCK_S_NAMES_S: [&str; 6] = [
"time_since_trade_s", "time_since_snap_s", "book_event_rate_per_s",
"spread_bps", "L1_imbalance", "micro_mid_drift",
];
let stacker_labels_f32: Vec<f32> = stacker_labels_u8.iter().map(|&y| y as f32).collect();
for (col_off, name) in BLOCK_S_NAMES_S.iter().enumerate() {
let global_col = BS_OFFSETS[col_off];
let feat_test: Vec<f32> = (stacker_train_n..val_pos.len())
.map(|r| {
let end_bar = labels.valid_indices[val_pos[r]];
feature_matrix[end_bar * alpha_dim + global_col]
})
.collect();
let strats = stratified_accuracy(&stacker_logits, &stacker_labels_f32, &feat_test, 5);
let q4 = strats.last().unwrap();
let q0 = &strats[0];
println!(" {name:>24}: Q0_acc={:.4} Q4_acc={:.4} (Q0_n={}, Q4_n={})",
q0.accuracy, q4.accuracy, q0.n, q4.n);
}
}
}
// ── Stratified accuracy on Block-S features (regime diagnostic) ──
//
// Block-S sits at cols 75..81 of the 81-dim snapshot row. For each
// val sequence, gather the END-BAR's feature value, then stratify
// val accuracy across quintiles of that feature. This tells us
// whether Mamba's K=6000 alpha is regime-conditional (concentrated
// in spread-Q4 / micro-drift-Q4 like the Phase 1c stateless MLP) or
// uniform — which decides whether we need an explicit regime head
// before the backtest.
const BLOCK_S_OFF: usize = 75;
const BLOCK_S_NAMES: [&str; 6] = [
"time_since_trade_s",
"time_since_snap_s",
"book_event_rate_per_s",
"spread_bps",
"L1_imbalance",
"micro_mid_drift",
];
let labels_u8_full: Vec<u8> =
val_ys.iter().map(|&y| if y > 0.5 { 1 } else { 0 }).collect();
let val_labels_f32: Vec<f32> = labels_u8_full.iter().map(|&y| y as f32).collect();
println!();
println!("--- Stratified val accuracy across Block-S features (5 quintiles) ---");
for (col_off, name) in BLOCK_S_NAMES.iter().enumerate() {
let global_col = BLOCK_S_OFF + col_off;
let feat: Vec<f32> = val_pos.iter()
.map(|&p| {
let end_bar = labels.valid_indices[p];
feature_matrix[end_bar * alpha_dim + global_col]
})
.collect();
let strats = stratified_accuracy(&val_logits, &val_labels_f32, &feat, 5);
println!();
println!(" ▸ feature: {}", name);
println!(" {:>4} {:>14} {:>14} {:>10} {:>10} {:>12}",
"k", "feat_lo", "feat_hi", "n", "accuracy", "observed");
for s in &strats {
println!(
" {:>4} {:>14.4} {:>14.4} {:>10} {:>10.4} {:>12.4}",
s.stratum, s.feature_lo, s.feature_hi, s.n, s.accuracy, s.observed_pos_rate
);
}
}
println!();
if auc > 0.55 {
println!("GATE PASS: AUC > 0.55 at K={}; multi-minute alpha confirmed.", cli.horizon);
println!(" Two-head FoxhuntQ-Δ architecture validated.");
} else if auc < 0.52 {
println!("GATE FAIL (decisive): AUC < 0.52; design dead in current form.");
println!(" Need different inputs / model class / context window.");
} else {
println!("GATE MARGINAL: 0.52 ≤ AUC ≤ 0.55; tune hyperparameters or extend seq_len.");
}
Ok(())
}