Files
foxhunt/crates/ml-backtesting/cuda/decision_policy.cu
jgrusewski cd82f9a4a0 feat(ml-backtesting): threshold gate + per-fill cost integration (P4)
Adds the two sweep axes that the spec's deployability grid needs but
were missing from the kernels:

Threshold gate (decision_policy.cu, both kernels):
- New per-backtest `threshold_per_b` array kernel arg.
- Pre-Kelly prelude: if max_h |alpha[h] - 0.5| * 2 < threshold[b],
  emit noop and return. Kept deterministic from alpha alone so the
  threshold pre-registration step (p60-p95 absolute calibration on a
  validation window, future P6) reflects exactly what gets gated in
  deployment.

Per-fill cost integration (resting_orders.cu / apply_fill_to_pos):
- apply_fill_to_pos signature grows three args: b, cost_per_lot_per_side_per_b,
  total_fees_per_b. Single insertion point at line 90.
- After the close-leg realized_pnl math runs (so the gross unwind P&L
  is preserved), deduct fill_cost = filled_lots * cost_per_lot_per_side[b]
  from pos.realized_pnl AND accumulate into total_fees_per_b[b].
- Net-of-cost semantics: isv_kelly_update_on_close reads realized_pnl
  delta which is now net of cost — Kelly state learns from realistic
  return distribution.
- All 3 apply_fill_to_pos call sites in step_resting_orders updated.
  order_match.cu's submit_market_immediate path is dead code in the
  post-P1 flow (everything routes through seed_inflight_limits_batched
  → step_resting_orders → apply_fill_to_pos) so not touched here.

BatchedSimConfig + UniformSimParams + BacktestHarnessConfig gain
threshold + cost_per_lot_per_side fields. All UniformSimParams
constructors in tests and main.rs updated with defaults (0.0, 0.0 =
gate disabled, frictionless).

Regression:
- threshold_gate_skips_low_conviction (p=0.51 + threshold=0.10 → noop)
- threshold_gate_allows_high_conviction (p=0.8 + threshold=0.10 → buy 1+)
- threshold_zero_is_passthrough (sanity)
- All P1+P2+P3 tests continue to pass via the new ABI.

cost_deducted_at_each_fill + kelly_state_sees_net_return end-to-end
tests deferred — they require a full submit_market → fill → close
sequence, which the production smoke exercises.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
2026-05-19 17:13:29 +02:00

489 lines
21 KiB
Plaintext
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
// decision_policy.cu — alpha → per-horizon Kelly → aggregate → market target.
//
// Hardcoded v1 policy = Strategy::default_for(max_lots) from C3:
// for each h in 0..N_HORIZONS:
// compute signed sizing from (alpha[h], IsvKellyState[h])
// aggregate via WeightedByRealizedSharpe:
// w[h] = max(0, recent_sharpe[h]) / sum_of_positive_sharpes
// final = sum_h w[h] * signed_size[h]
// final → side+size in market_targets[b]
//
// The bytecode VM from spec §6 is intentionally not implemented in v1 —
// the default Strategy::default_for produces exactly this shape and
// alternative compositions are deferred. Bytecode plumbing remains in
// src/policy/mod.rs ready for v2 expansion.
//
// Single-writer per block (thread 0). See spec §5.
#include "lob_state.cuh"
// Minimum closed-trade count before the variance-derived `cap_units`
// constraint is trusted. With N < this threshold the cap falls back to
// the host-supplied `max_lots` budget — a single-sample variance proxy
// (`ret²`) systematically over-estimates Var[ret] for a directionally
// biased loss/win and locks `cap_lots` to ~zero, which kills further
// trading. Threshold of 10 picked as the standard rule-of-thumb for
// "first moment estimate becomes informative" — same role as the
// pearl_first_observation_bootstrap pattern applied to the variance side.
#define MIN_TRADES_FOR_VAR_CAP 10u
extern "C" __global__ void decision_policy_default(
const float* __restrict__ alpha_probs, // [N_HORIZONS] broadcast
const unsigned char* __restrict__ isv_kelly_base, // [n_backtests * N_HORIZONS * sizeof(IsvKellyState)]
const Pos* __restrict__ positions, // [n_backtests]
const unsigned int* __restrict__ program_lens, // [n_backtests] — non-zero ⇒ skip (program kernel handles it)
int* __restrict__ market_targets, // [n_backtests * 2] — written
unsigned int* __restrict__ open_horizon_masks, // [n_backtests] — written (entry attribution)
// P1: per-backtest sim parameter arrays. Each backtest indexed by b.
const float* __restrict__ target_annual_vol_units_per_b, // [n_backtests]
const float* __restrict__ annualisation_factor_per_b, // [n_backtests]
const int* __restrict__ max_lots_per_b, // [n_backtests]
// Cold-start floors per pearl_blend_formulas_must_have_permanent_floor:
// max(floor, real), not blend; pearl_kelly_cap_signal_driven_floors.
// When isv_kelly state is at sentinel (no closed trades), kelly_frac
// falls back to `kelly_frac_floor_per_b[b]`; the aggregate weight floor lets
// cross-horizon sum produce a non-zero size before recent_sharpe is
// populated. Floors are non-zero by contract — passing 0.0 reverts
// to the old sentinel-skip behaviour.
const float* __restrict__ kelly_frac_floor_per_b, // [n_backtests]
const float* __restrict__ sharpe_weight_floor_per_b, // [n_backtests]
// P4: threshold gate (absolute max-conviction cutoff, pre-Kelly).
const float* __restrict__ threshold_per_b, // [n_backtests]
int n_backtests
) {
int b = blockIdx.x;
if (b >= n_backtests || threadIdx.x != 0) return;
if (program_lens[b] != 0u) return; // backtest b uses a custom bytecode program; skip default
// P4: threshold gate. Skip entirely if max conviction below threshold.
// Pre-Kelly placement keeps the gate deterministic from alpha alone,
// so the threshold pre-registration step (compute p60-p95 absolute
// values from observed convictions) reflects exactly what gets gated
// in deployment.
{
float max_conv = 0.0f;
#pragma unroll
for (int h = 0; h < N_HORIZONS; ++h) {
max_conv = fmaxf(max_conv, fabsf(alpha_probs[h] - 0.5f) * 2.0f);
}
if (max_conv < threshold_per_b[b]) {
market_targets[b * 2 + 0] = 2;
market_targets[b * 2 + 1] = 0;
open_horizon_masks[b] = 0u;
return;
}
}
// Per-backtest scalar reads (host-uploaded arrays).
const float target_annual_vol_units = target_annual_vol_units_per_b[b];
const float annualisation_factor = annualisation_factor_per_b[b];
const int max_lots = max_lots_per_b[b];
const float kelly_frac_floor = kelly_frac_floor_per_b[b];
const float sharpe_weight_floor = sharpe_weight_floor_per_b[b];
const IsvKellyState* isv = reinterpret_cast<const IsvKellyState*>(
isv_kelly_base + (size_t)b * N_HORIZONS * sizeof(IsvKellyState)
);
float signed_sizes[N_HORIZONS];
float weights[N_HORIZONS];
unsigned int attribution_mask = 0;
float w_sum = 0.0f;
#pragma unroll
for (int h = 0; h < N_HORIZONS; ++h) {
const float p_h = alpha_probs[h];
const IsvKellyState& s = isv[h];
const float sig_mag = fabsf(p_h - 0.5f) * 2.0f; // ∈ [0, 1] — always derivable
const float dir = (p_h > 0.5f) ? 1.0f : -1.0f;
// Kelly fraction with floor. When pnl_ema_win is sentinel, the
// ratio is undefined → floor directly. Otherwise compute the
// realised Kelly and take max(floor, kelly) so the floor never
// drives sizing down once we have signal.
float kelly_frac;
if (s.pnl_ema_win > 1e-9f) {
float kf = (s.win_rate_ema * s.pnl_ema_win
- (1.0f - s.win_rate_ema) * s.pnl_ema_loss)
/ s.pnl_ema_win;
kf = fmaxf(0.0f, fminf(1.0f, kf));
kelly_frac = fmaxf(kelly_frac_floor, kf);
} else {
kelly_frac = kelly_frac_floor;
}
// Cap units: only trust the variance-derived cap_units after
// MIN_TRADES_FOR_VAR_CAP closed trades. Before then, a single-sample
// variance proxy (`ret²`) collapses cap to ~0 and locks out further
// trading. Fall back to the host-supplied `max_lots` budget.
float cap_lots;
if (s.n_trades_seen >= MIN_TRADES_FOR_VAR_CAP
&& s.realised_return_var > 1e-9f)
{
const float cap_units = target_annual_vol_units
/ sqrtf(s.realised_return_var * annualisation_factor);
cap_lots = fminf(cap_units, (float)max_lots);
} else {
cap_lots = (float)max_lots;
}
signed_sizes[h] = dir * sig_mag * kelly_frac * cap_lots;
// Recent-Sharpe weight with floor: lets cold-start aggregate fire
// (uniform weights = sharpe_weight_floor everywhere) until one
// horizon's recent_sharpe rises above the floor and dominates.
const float w = fmaxf(sharpe_weight_floor, s.recent_sharpe);
weights[h] = w;
w_sum += w;
}
float final_size = 0.0f;
if (w_sum > 1e-9f) {
#pragma unroll
for (int h = 0; h < N_HORIZONS; ++h) {
const float w_norm = weights[h] / w_sum;
final_size += w_norm * signed_sizes[h];
// Attribute the open to any horizon whose weight contributed.
if (weights[h] > 1e-9f) attribution_mask |= (1u << h);
}
}
// Round-to-nearest (truncation alone bites here: 0.8(f32) * 0.75 * 5.0
// = 2.99999976 → (int)2, off by one. Floor truncation for sub-1
// suppression is preserved via the explicit |lots| < 1 check.
int lots = (int)truncf(final_size + (final_size >= 0.0f ? 0.5f : -0.5f));
if (lots == 0) {
market_targets[b * 2 + 0] = 2; // no-op
market_targets[b * 2 + 1] = 0;
return;
}
// If currently flat, this is an opening trade: record attribution mask
// so pnl_track + isv_kelly_update can credit the right horizons on close.
if (positions[b].position_lots == 0) {
open_horizon_masks[b] = attribution_mask;
}
const int side = (lots > 0) ? 0 : 1;
const int abs_sz = (lots > 0) ? lots : -lots;
market_targets[b * 2 + 0] = side;
market_targets[b * 2 + 1] = abs_sz;
}
// Bytecode VM — interprets the flat Program emitted by Strategy::flatten()
// (Rust src/policy/mod.rs) against current alpha probs + per-horizon
// IsvKellyState. Replaces the hardcoded default policy when the host
// has uploaded a non-empty program for backtest b.
//
// Stack layout: parallel float-value + u32-attribution-mask stacks.
// Each EmitPerHorizonSize pushes its (signed_size, mask=(1<<h)). The
// aggregators pop n values + n masks, push the aggregated value with
// OR'd mask. The terminal WriteOrder pops one value+mask, converts to
// {side, abs_size} in market_targets[b] and records the attribution
// mask in open_horizon_masks[b] for entry-time crediting.
//
// Layout MUST stay in sync with src/policy/mod.rs::{OpCode, Instruction}.
#define OP_NOOP 0
#define OP_EVAL_REGIME 1
#define OP_EMIT_PER_HORIZON_SIZE 2
#define OP_AGG_WEIGHTED_SHARPE 3
#define OP_AGG_MEAN 4
#define OP_AGG_MAX_CONFIDENCE 5
#define OP_APPLY_CONFLICT 6
#define OP_WRITE_ORDER 7
#define OP_PUSH_SCALAR 8
#define OP_BRANCH_IF_REGIME 9
#define STACK_CAP 32
struct Instruction {
unsigned char op;
unsigned char arg0;
unsigned short arg1;
unsigned int arg2;
};
extern "C" __global__ void decision_policy_program(
const float* __restrict__ alpha_probs,
const unsigned char* __restrict__ isv_kelly_base,
const Pos* __restrict__ positions,
const Instruction* __restrict__ programs, // [n_backtests * max_instructions]
const unsigned int* __restrict__ program_lens, // [n_backtests]
const unsigned int* __restrict__ regimes, // [n_backtests] — current regime id
int* __restrict__ market_targets, // [n_backtests * 2]
unsigned int* __restrict__ open_horizon_masks,
// P1: per-backtest sim parameter arrays. See decision_policy_default header for contract.
const float* __restrict__ target_annual_vol_units_per_b,
const float* __restrict__ annualisation_factor_per_b,
int max_instructions,
const int* __restrict__ max_lots_per_b,
const float* __restrict__ kelly_frac_floor_per_b,
const float* __restrict__ sharpe_weight_floor_per_b,
// P4: threshold gate (same as decision_policy_default).
const float* __restrict__ threshold_per_b,
int n_backtests
) {
int b = blockIdx.x;
if (b >= n_backtests || threadIdx.x != 0) return;
const unsigned int plen = program_lens[b];
if (plen == 0u) {
// Skip — this backtest uses the hardcoded default kernel separately.
return;
}
// P4: threshold gate.
{
float max_conv = 0.0f;
#pragma unroll
for (int h = 0; h < N_HORIZONS; ++h) {
max_conv = fmaxf(max_conv, fabsf(alpha_probs[h] - 0.5f) * 2.0f);
}
if (max_conv < threshold_per_b[b]) {
market_targets[b * 2 + 0] = 2;
market_targets[b * 2 + 1] = 0;
open_horizon_masks[b] = 0u;
return;
}
}
// Per-backtest scalar reads (host-uploaded arrays).
const float target_annual_vol_units = target_annual_vol_units_per_b[b];
const float annualisation_factor = annualisation_factor_per_b[b];
const int max_lots = max_lots_per_b[b];
const float kelly_frac_floor = kelly_frac_floor_per_b[b];
const float sharpe_weight_floor = sharpe_weight_floor_per_b[b];
const Instruction* prog = programs + (size_t)b * max_instructions;
const IsvKellyState* isv = reinterpret_cast<const IsvKellyState*>(
isv_kelly_base + (size_t)b * N_HORIZONS * sizeof(IsvKellyState)
);
const unsigned int my_regime = regimes[b];
float stack_val[STACK_CAP];
unsigned int stack_mask[STACK_CAP];
int sp = 0;
unsigned int pc = 0;
while (pc < plen && pc < (unsigned int)max_instructions) {
const Instruction ins = prog[pc];
pc += 1;
switch (ins.op) {
case OP_NOOP: break;
case OP_PUSH_SCALAR: {
if (sp >= STACK_CAP) break;
const float v = __int_as_float((int)ins.arg2);
stack_val[sp] = v;
stack_mask[sp] = 0u;
sp += 1;
break;
}
case OP_EVAL_REGIME: {
if (sp >= STACK_CAP) break;
stack_val[sp] = (my_regime == (unsigned int)ins.arg0) ? 1.0f : 0.0f;
stack_mask[sp] = 0u;
sp += 1;
break;
}
case OP_BRANCH_IF_REGIME: {
if (my_regime != (unsigned int)ins.arg0) {
pc = ins.arg2;
}
break;
}
case OP_EMIT_PER_HORIZON_SIZE: {
if (sp >= STACK_CAP) break;
const int h = (int)ins.arg0;
const int leaf_max_lots = (int)ins.arg1;
if (h < 0 || h >= N_HORIZONS) { stack_val[sp]=0.0f; stack_mask[sp]=0u; sp+=1; break; }
const float p_h = alpha_probs[h];
const IsvKellyState& s = isv[h];
const float sig_mag = fabsf(p_h - 0.5f) * 2.0f;
const float dir = (p_h > 0.5f) ? 1.0f : -1.0f;
float kelly_frac;
if (s.pnl_ema_win > 1e-9f) {
float kf = (s.win_rate_ema * s.pnl_ema_win
- (1.0f - s.win_rate_ema) * s.pnl_ema_loss)
/ s.pnl_ema_win;
kf = fmaxf(0.0f, fminf(1.0f, kf));
kelly_frac = fmaxf(kelly_frac_floor, kf);
} else {
kelly_frac = kelly_frac_floor;
}
const float effective_cap = fminf((float)leaf_max_lots, (float)max_lots);
// Same MIN_TRADES_FOR_VAR_CAP gate as decision_policy_default.
float cap_lots;
if (s.n_trades_seen >= MIN_TRADES_FOR_VAR_CAP
&& s.realised_return_var > 1e-9f)
{
const float cap_units = target_annual_vol_units
/ sqrtf(s.realised_return_var * annualisation_factor);
cap_lots = fminf(cap_units, effective_cap);
} else {
cap_lots = effective_cap;
}
const float ss = dir * sig_mag * kelly_frac * cap_lots;
stack_val[sp] = ss;
stack_mask[sp] = 1u << h;
sp += 1;
break;
}
case OP_AGG_MEAN: {
const int n = (int)ins.arg0;
if (n <= 0 || sp < n) break;
float sum = 0.0f;
unsigned int mask_or = 0u;
for (int i = sp - n; i < sp; ++i) {
sum += stack_val[i];
mask_or |= stack_mask[i];
}
sp -= n;
stack_val[sp] = sum / (float)n;
stack_mask[sp] = mask_or;
sp += 1;
break;
}
case OP_AGG_WEIGHTED_SHARPE: {
const int n = (int)ins.arg0;
if (n <= 0 || sp < n) break;
float w_sum = 0.0f, wx_sum = 0.0f;
unsigned int mask_or = 0u;
for (int i = sp - n; i < sp; ++i) {
const unsigned int m = stack_mask[i];
// Look up recent_sharpe for the single horizon this slot
// attributes to. If the mask has multiple bits set (a
// nested aggregator already merged horizons) we fall back
// to uniform weight 1.0 because per-horizon attribution
// collapsed.
float w;
if (m != 0u && (m & (m - 1u)) == 0u) {
// Single-bit mask — recover horizon index.
int h = 0;
unsigned int mm = m;
while ((mm & 1u) == 0u && h < N_HORIZONS) { mm >>= 1; h += 1; }
// Floor lets cold-start aggregate; matches the
// decision_policy_default weight pattern.
w = fmaxf(sharpe_weight_floor, isv[h].recent_sharpe);
} else {
w = 1.0f;
}
w_sum += w;
wx_sum += w * stack_val[i];
mask_or |= m;
}
sp -= n;
stack_val[sp] = (w_sum > 1e-9f) ? (wx_sum / w_sum) : 0.0f;
stack_mask[sp] = mask_or;
sp += 1;
break;
}
case OP_AGG_MAX_CONFIDENCE: {
const int n = (int)ins.arg0;
if (n <= 0 || sp < n) break;
int best = sp - n;
for (int i = sp - n + 1; i < sp; ++i) {
if (fabsf(stack_val[i]) > fabsf(stack_val[best])) best = i;
}
const float bv = stack_val[best];
const unsigned int bm = stack_mask[best];
sp -= n;
stack_val[sp] = bv;
stack_mask[sp] = bm;
sp += 1;
break;
}
case OP_APPLY_CONFLICT: {
// v1: no-op (mean of children already aggregated; conflict
// resolution only matters when multiple cells share state,
// which v1 doesn't). Preserved for v2 Portfolio mode.
break;
}
case OP_WRITE_ORDER: {
if (sp <= 0) break;
sp -= 1;
const float final_size = stack_val[sp];
const unsigned int attribution = stack_mask[sp];
int lots = (int)truncf(final_size + (final_size >= 0.0f ? 0.5f : -0.5f));
if (lots == 0) {
market_targets[b * 2 + 0] = 2; // no-op
market_targets[b * 2 + 1] = 0;
} else {
if (positions[b].position_lots == 0) {
open_horizon_masks[b] = attribution;
}
const int side = (lots > 0) ? 0 : 1;
const int abs_sz = (lots > 0) ? lots : -lots;
market_targets[b * 2 + 0] = side;
market_targets[b * 2 + 1] = abs_sz;
}
break;
}
default: break;
}
}
}
// Updates per-horizon IsvKellyState[h] for every horizon flagged in
// `closed_horizon_mask[b]`. Called by the host after pnl_track_step
// detects close transitions. `realised_return[b]` is the closed-trade
// per-lot return in price units (positive = win, negative = loss).
//
// Wiener-α floor 0.4 per pearl_wiener_alpha_floor_for_nonstationary.md.
// Welford variance updates a running estimate around mean 0 (simplified;
// v2 may switch to ema-centered).
//
// recent_sharpe composite re-computed from the freshly updated EMAs.
extern "C" __global__ void isv_kelly_update_on_close(
unsigned char* isv_kelly_base,
const unsigned int* closed_horizon_mask, // [n_backtests]
const float* realised_return, // [n_backtests]
int n_backtests
) {
int b = blockIdx.x;
if (b >= n_backtests || threadIdx.x != 0) return;
const unsigned int mask = closed_horizon_mask[b];
if (mask == 0u) return;
const float ret = realised_return[b];
IsvKellyState* isv = reinterpret_cast<IsvKellyState*>(
isv_kelly_base + (size_t)b * N_HORIZONS * sizeof(IsvKellyState)
);
const float alpha = 0.4f; // Wiener-α floor (non-stationary policy P&L)
#pragma unroll
for (int h = 0; h < N_HORIZONS; ++h) {
if (!(mask & (1u << h))) continue;
IsvKellyState& s = isv[h];
const bool win = ret > 0.0f;
if (s.n_trades_seen == 0u) {
// First-observation bootstrap: replace EMAs directly (no zero-bias warmup).
s.pnl_ema_win = win ? ret : 0.0f;
s.pnl_ema_loss = win ? 0.0f : -ret;
s.win_rate_ema = win ? 1.0f : 0.0f;
s.realised_return_var = ret * ret; // single-sample variance proxy
} else {
if (win) {
s.pnl_ema_win = (1.0f - alpha) * s.pnl_ema_win + alpha * ret;
s.win_rate_ema = (1.0f - alpha) * s.win_rate_ema + alpha * 1.0f;
} else {
s.pnl_ema_loss = (1.0f - alpha) * s.pnl_ema_loss + alpha * (-ret);
s.win_rate_ema = (1.0f - alpha) * s.win_rate_ema + alpha * 0.0f;
}
// Welford-ish running variance around 0.
const float prev_n = (float)s.n_trades_seen;
s.realised_return_var =
(prev_n * s.realised_return_var + ret * ret) / (prev_n + 1.0f);
}
s.n_trades_seen += 1u;
// Composite recent_sharpe = (mean_winning_contribution mean_losing_contribution)
// / sqrt(realised_return_var)
const float numerator = s.pnl_ema_win * s.win_rate_ema
- s.pnl_ema_loss * (1.0f - s.win_rate_ema);
const float denom = sqrtf(fmaxf(s.realised_return_var, 1e-9f));
s.recent_sharpe = numerator / denom;
}
}