Files
resonance-engine/_analyze_overnight.py
T
2026-06-06 20:34:31 +07:00

289 lines
12 KiB
Python

"""_analyze_overnight.py — analyze overnight_v4 results.
Reads /mnt/d/Resonance_Engine/traj/overnight_<RUN_ID>/ and produces:
1. Per-variable mean/std of each channel during pulses, averaged over 3 reps.
2. CRITICAL: variable ranking PER rep. If the ranking is preserved across all
3 reps with different injection orderings -> real per-variable effect.
If ranking tracks injection position -> time-drift artifact.
3. Pairwise variable trajectory distances (Euclidean over resampled trace),
compared to within-variable self-distances across reps.
4. Decay shape: tail-decay (last 30s of pulses vs final 60s post) per variable.
Usage:
python3 _analyze_overnight.py <RUN_ID>
python3 _analyze_overnight.py latest
"""
from __future__ import annotations
import json, sys
from pathlib import Path
import numpy as np
import pandas as pd
from scipy.stats import spearmanr
ROOT = Path("/mnt/d/Resonance_Engine/traj")
CHANNELS = ["asymmetry", "coherence", "stress_xx", "stress_yy", "stress_xy",
"vorticity_mean", "vel_mean", "vel_max", "vel_var"]
VARS = ["mid_price", "signed_flow_usd", "trade_count",
"large_print_cnt", "wallet_entropy", "vwap_drift"]
def find_run(arg: str) -> Path:
if arg == "latest":
runs = sorted([p for p in ROOT.glob("overnight_*") if p.is_dir()])
if not runs:
sys.exit("no overnight_* runs in " + str(ROOT))
return runs[-1]
p = ROOT / f"overnight_{arg}" if not arg.startswith("overnight_") else ROOT / arg
if not p.exists():
sys.exit(f"missing: {p}")
return p
def load_rep(rep_dir: Path) -> dict[str, pd.DataFrame]:
out = {}
for v in VARS:
f = rep_dir / f"{v}.parquet"
if f.exists():
out[v] = pd.read_parquet(f)
return out
def pulse_window_stats(df: pd.DataFrame, rel_t_col: str = "rel_t") -> dict:
"""Stats during pulse phase (rel_t between 0 and end-of-pulses) per channel."""
# exclude pre (rel_t<0) and the very-far post-tail; keep pulse + first 5s post
if rel_t_col not in df.columns:
return {}
pulse_end = df[rel_t_col].max() - 55 # approx end of pulses
mask = (df[rel_t_col] >= 0) & (df[rel_t_col] <= pulse_end)
sub = df[mask]
out = {}
for ch in CHANNELS:
if ch in sub.columns:
out[ch + "_mean"] = float(sub[ch].mean())
out[ch + "_std"] = float(sub[ch].std())
out[ch + "_max"] = float(sub[ch].max())
out[ch + "_min"] = float(sub[ch].min())
return out
def decay_stats(df: pd.DataFrame, rel_t_col: str = "rel_t") -> dict:
"""Compare last 20s of pulse phase vs last 20s of post (decay magnitude)."""
if rel_t_col not in df.columns or len(df) < 10:
return {}
t_max = df[rel_t_col].max()
pulse_end = t_max - 55 # 60s post window approx
tail_pulse_mask = (df[rel_t_col] >= pulse_end - 20) & (df[rel_t_col] <= pulse_end)
tail_post_mask = (df[rel_t_col] >= t_max - 20) & (df[rel_t_col] <= t_max)
tp = df[tail_pulse_mask]
pp = df[tail_post_mask]
out = {}
for ch in ["asymmetry", "coherence", "stress_xy"]:
if ch in df.columns and len(tp) > 2 and len(pp) > 2:
out[ch + "_decay"] = float(tp[ch].mean() - pp[ch].mean())
return out
def resample_trajectory(df: pd.DataFrame, channel: str, n_pts: int = 200) -> np.ndarray:
"""Resample one channel on a uniform rel_t grid."""
if "rel_t" not in df.columns or channel not in df.columns:
return np.full(n_pts, np.nan)
d = df.dropna(subset=["rel_t", channel]).sort_values("rel_t")
if len(d) < 4:
return np.full(n_pts, np.nan)
t0, t1 = d["rel_t"].min(), d["rel_t"].max()
grid = np.linspace(t0, t1, n_pts)
return np.interp(grid, d["rel_t"].values, d[channel].values)
def main():
run_id_arg = sys.argv[1] if len(sys.argv) > 1 else "latest"
run_dir = find_run(run_id_arg)
print(f"=== analyzing {run_dir.name} ===\n")
meta_p = run_dir / "meta.json"
if meta_p.exists():
meta = json.loads(meta_p.read_text())
print(f"run_id={meta.get('run_id')} started={meta.get('started')}")
print(f"vars={meta.get('vars')}")
print(f"orderings_planned={meta.get('orderings', 'n/a')}")
print()
# baseline
base_stats_p = run_dir / "baseline_clean_stats.json"
if base_stats_p.exists():
bs = json.loads(base_stats_p.read_text())
print("=== verified-clean baseline sigma ===")
for k, v in bs.items():
if isinstance(v, float):
print(f" {k:30s} {v:.6f}")
else:
print(f" {k:30s} {v}")
print()
sigma = {k.replace("sigma_", ""): v for k, v in bs.items() if k.startswith("sigma_")}
else:
sigma = {}
print("(no baseline_clean_stats.json found)\n")
# ----- load all reps -----
rep_dirs = sorted([p for p in run_dir.iterdir() if p.is_dir() and p.name.startswith("rep_")])
if not rep_dirs:
sys.exit("no rep_* dirs found yet")
print(f"found {len(rep_dirs)} reps: {[p.name for p in rep_dirs]}\n")
# per-rep, per-var stats
rep_stats = {} # rep_name -> var -> stats dict
rep_order = {} # rep_name -> var -> int (order in this rep)
for rd in rep_dirs:
rep_stats[rd.name] = {}
# try to read order from a meta file
order_p = rd / "rep_meta.json"
if order_p.exists():
rm = json.loads(order_p.read_text())
for i, v in enumerate(rm.get("order", [])):
rep_order.setdefault(rd.name, {})[v] = i
reps = load_rep(rd)
for v, df in reps.items():
rep_stats[rd.name][v] = pulse_window_stats(df)
rep_stats[rd.name][v].update(decay_stats(df))
# ----- table 1: per-var per-rep asymmetry mean -----
print("=== ASYMMETRY mean during pulse phase (per rep, per var) ===")
header = f"{'var':<20}" + "".join(f"{rd.name:>12}" for rd in rep_dirs) + f"{'mean':>10}{'std':>8}"
print(header)
print("-" * len(header))
rows = []
for v in VARS:
vals = []
for rd in rep_dirs:
s = rep_stats[rd.name].get(v, {})
vals.append(s.get("asymmetry_mean", np.nan))
m = np.nanmean(vals)
sd = np.nanstd(vals)
rows.append((v, vals, m, sd))
line = f"{v:<20}" + "".join(f"{x:>12.3f}" if not np.isnan(x) else f"{'--':>12}" for x in vals)
line += f"{m:>10.3f}{sd:>8.3f}"
print(line)
print()
# ----- table 2: cross-rep ranking preservation (THE critical test) -----
print("=== RANKING preservation test (variable effect vs drift artifact) ===")
rankings = {}
for rd in rep_dirs:
valid = [(v, rep_stats[rd.name].get(v, {}).get("asymmetry_mean", np.nan)) for v in VARS]
valid = [vv for vv in valid if not np.isnan(vv[1])]
valid.sort(key=lambda x: x[1])
rankings[rd.name] = [vv[0] for vv in valid]
print(f" {rd.name}: low->high asym = {' < '.join(rankings[rd.name])}")
print()
# Spearman rank correlation between all rep pairs
if len(rep_dirs) >= 2:
print(" Spearman rank corr between rep pairs (1.0 = perfect agreement, 0 = random):")
all_vars_sorted = sorted(VARS)
for i in range(len(rep_dirs)):
for j in range(i + 1, len(rep_dirs)):
ri = [rep_stats[rep_dirs[i].name].get(v, {}).get("asymmetry_mean", np.nan) for v in all_vars_sorted]
rj = [rep_stats[rep_dirs[j].name].get(v, {}).get("asymmetry_mean", np.nan) for v in all_vars_sorted]
mask = ~(np.isnan(ri) | np.isnan(rj))
if mask.sum() >= 3:
corr, p = spearmanr(np.array(ri)[mask], np.array(rj)[mask])
print(f" {rep_dirs[i].name} vs {rep_dirs[j].name}: rho={corr:+.3f} p={p:.3f}")
print()
# Also test: does asym_mean correlate with INJECTION ORDER (drift artifact)?
if rep_order:
print(" ORDER-POSITION vs asymmetry_mean (high pos corr = time-drift artifact):")
for rd in rep_dirs:
if rd.name not in rep_order:
continue
pairs = []
for v in VARS:
pos = rep_order[rd.name].get(v)
val = rep_stats[rd.name].get(v, {}).get("asymmetry_mean", np.nan)
if pos is not None and not np.isnan(val):
pairs.append((pos, val))
if len(pairs) >= 3:
pos = np.array([p[0] for p in pairs])
val = np.array([p[1] for p in pairs])
r, pp = spearmanr(pos, val)
print(f" {rd.name}: rho(order_pos, asym_mean) = {r:+.3f} p={pp:.3f}")
print()
# ----- table 3: per-channel between/within variance ratio (across-rep means) -----
print("=== Between-var vs within-var std ratio (per channel, using rep-averaged values) ===")
print(" ratio > 1 means variables differ from each other more than reps of same var disagree")
print(f"{'channel':<22}{'between_std':>14}{'within_std':>14}{'ratio':>10}")
for ch in CHANNELS:
per_var_means = []
per_var_within = []
for v in VARS:
vals = [rep_stats[rd.name].get(v, {}).get(ch + "_mean", np.nan) for rd in rep_dirs]
vals = [x for x in vals if not np.isnan(x)]
if len(vals) >= 1:
per_var_means.append(np.mean(vals))
if len(vals) >= 2:
per_var_within.append(np.std(vals))
if len(per_var_means) >= 2 and per_var_within:
bs = float(np.std(per_var_means))
ws = float(np.mean(per_var_within))
ratio = bs / ws if ws > 1e-12 else float('inf')
print(f"{ch:<22}{bs:>14.4f}{ws:>14.4f}{ratio:>10.2f}")
print()
# ----- table 4: pairwise variable trajectory distances on asymmetry -----
print("=== Pairwise variable trajectory distance (asymmetry, rep-averaged) ===")
# Build rep-averaged trajectory per var
avg_traj = {}
for v in VARS:
traces = []
for rd in rep_dirs:
f = rd / f"{v}.parquet"
if f.exists():
df = pd.read_parquet(f)
t = resample_trajectory(df, "asymmetry", 200)
if not np.all(np.isnan(t)):
traces.append(t)
if traces:
avg_traj[v] = np.nanmean(np.stack(traces), axis=0)
keys = list(avg_traj.keys())
if len(keys) >= 2:
print(f"{'':<20}" + "".join(f"{k[:10]:>12}" for k in keys))
for i, ka in enumerate(keys):
row = f"{ka:<20}"
for kb in keys:
if ka == kb:
row += f"{'-':>12}"
else:
d = float(np.linalg.norm(avg_traj[ka] - avg_traj[kb]))
row += f"{d:>12.3f}"
print(row)
print()
# Self-distance: distance between rep-pair trajectories for same var
print(" Self-distance (cross-rep traj distance for same var, sanity):")
for v in VARS:
ts = []
for rd in rep_dirs:
f = rd / f"{v}.parquet"
if f.exists():
df = pd.read_parquet(f)
t = resample_trajectory(df, "asymmetry", 200)
if not np.all(np.isnan(t)):
ts.append(t)
if len(ts) >= 2:
ds = []
for i in range(len(ts)):
for j in range(i + 1, len(ts)):
ds.append(float(np.linalg.norm(ts[i] - ts[j])))
if ds:
print(f" {v:<20} mean_self_dist={np.mean(ds):.3f} n_pairs={len(ds)}")
print("\n=== DONE ===")
if __name__ == "__main__":
main()