diff options
| -rw-r--r-- | RESULTS.md | 18 | ||||
| -rw-r--r-- | ROADMAP.md | 6 | ||||
| -rw-r--r-- | THEORY.md | 222 | ||||
| -rwxr-xr-x | experiments/finalize_claims.sh | 2 | ||||
| -rw-r--r-- | experiments/verify_theory.py | 157 |
5 files changed, 405 insertions, 0 deletions
@@ -429,6 +429,23 @@ the frozen 90% d4 threshold; **no prefix is eligible and task seeds 6--8 remain teaching directions are often strongly aligned even in mediocre runs, pointing next to regression gain/tracking over the moving trajectory rather than simply adding vectorizer capacity. +## Theory package and executable predictions + +`THEORY.md` proves the exact linearized variance of simultaneous Rademacher perturbation, +including its cross-layer term; bounds centered finite-sigma bias; gives a smooth-loss descent +condition that separates angle, gain, and curvature; proves the neutral-period innovation/SNR +result and task-fit absorption failure; and tabulates query/work/memory assumptions across the +main comparison set. + +The deterministic numerical verifier covers every frozen C5 axis. Across depth 2/4/8, width 8/32, +and K1/4/16, empirical perturbation MSE is within 1% of the exact prediction in the reported cases. +The finite-difference bias has log-log slope `2.000000` in sigma. A quadratic example changes from +descent to ascent on opposite sides of the predicted learning-rate threshold. Neutral predictor +residual power matches `(1-eta_P)^(2t)` across three rates and four update counts. Finally, +neutral residualization raises the teaching squared cosine from `0.7492` to `0.9656`, whereas +task-period fitting lowers it to `0.5395` by removing soma-correlated teaching. **The frozen C5 +theory package and numerical checks pass.** + ## Paper-faithful non-backprop baselines `PETITE` in the initial request is treated as **PEPITA**. Implementations and hyperparameters were @@ -496,6 +513,7 @@ canonical d2 comparison: `experiments/ep_depth2_sweep.sh`; near-parameter compar `experiments/ep_matched_sweep.sh`. C2 direct causal diagnosis: `experiments/c2_nodepert_validation.sh` followed by `python experiments/analyze_c2_nodepert_validation.py`. +Theory checks: `python experiments/verify_theory.py`. Audited aggregate tables: `experiments/analyze_verified.py`. The claim-locked finalizer is `experiments/finalize_claims.sh`; it rejects missing/dirty inputs, regenerates both main figures, their captions, source-hash manifest, and `results/audited_tables.md`, then runs all smoke checks. @@ -250,6 +250,12 @@ The minimum theory package contains: Numerical simulations must test the predicted scaling with depth, width, directions, perturbation scale, and predictor timescale. +**Status: passed on 2026-07-22.** `THEORY.md` gives the exact simultaneous-Rademacher MSE with its +cross-layer term, the Gaussian layerwise comparison, an `O(sigma^2)` antithetic bias bound, the +smooth-loss angle/gain/curvature descent condition, the neutral innovation/SNR theorem and +task-fit absorption result, and a phase/query/work/memory table. `experiments/verify_theory.py` +checks depth, width, K, sigma, and predictor timescale deterministically; all assertions pass. + ## Protocol discipline - Hyperparameters are selected on a validation split created only from the training set. diff --git a/THEORY.md b/THEORY.md new file mode 100644 index 0000000..25d27ef --- /dev/null +++ b/THEORY.md @@ -0,0 +1,222 @@ +# SDIL theory and cost audit + +This note states the minimum claims that can be proved for the implemented +estimators and innovation rule. It deliberately separates an unbiased neural +teaching direction from a useful parameter update: the completed experiments +show that a high hidden-state cosine can coexist with poor optimization. + +## 1. Simultaneous multi-layer perturbation + +Consider one example and collect the hidden-state gradients +`g_l = dL/dh_l` for hidden blocks `l=1,...,H`, where `g_l` has dimension +`d_l`. For one simultaneous Rademacher intervention, draw independent +`xi_l in {-1,+1}^{d_l}` and linearize the centered loss difference as + +```text +D = sum_j <g_j, xi_j>. +``` + +The gradient and descent estimators for layer `l` are + +```text +g_hat_l = xi_l D, q_l = -g_hat_l. +``` + +Independence and `E[xi_l xi_l^T] = I` give + +```text +E[g_hat_l] = g_l. +``` + +Because a Rademacher vector has deterministic squared norm `d_l`, the exact +linearized mean-squared error is + +```text +E ||g_hat_l - g_l||^2 + = (d_l - 1) ||g_l||^2 + + d_l sum_{j != l} ||g_j||^2. +``` + +Averaging `K` independent directions divides this expression by `K`. The first +term is ordinary directional-estimator variance; the second is cross-layer +interference. Thus simultaneous perturbation is unbiased but its relative +variance grows with layer width and with the total gradient energy injected at +other layers. This predicts why K1 can work in an easy/depth-flat regime yet +fail to train a moving context vectorizer on a compositional task. + +For comparison, layerwise Rademacher perturbation removes the cross-layer term. +The current layerwise implementation uses Gaussian directions, for which + +```text +E ||g_hat_l - g_l||^2 = (d_l + 1) ||g_l||^2 / K. +``` + +Layerwise targets therefore buy variance reduction with depth-dependent replay +cost; they are not a free scalable replacement for simultaneous calibration. + +### Finite perturbation bias + +Suppose the third derivative of the loss has operator norm at most `M` in the +perturbation neighborhood, and let `D_total = sum_l d_l`. A centered antithetic +difference has directional error at most + +```text +M sigma^2 ||xi||^3 / 6. +``` + +For simultaneous Rademacher directions this gives the conservative layerwise +bias bound + +```text +||E[g_hat_l] - g_l|| + <= (M sigma^2 / 6) sqrt(d_l) D_total^(3/2). +``` + +The important dependence is `O(sigma^2)`, with a dimension-dependent constant. +Reducing K removes variance but not this bias; increasing width/depth without a +sigma audit is therefore unsafe. + +## 2. What is sufficient for a descent step + +Flatten all forward parameters into `theta`, let `p = grad L(theta)`, and let +`v` be the complete local update before multiplying by the learning rate: + +```text +theta_next = theta + eta v. +``` + +If the loss is `beta`-smooth and + +```text +c = <-p, v> / (||p|| ||v||) > 0, +``` + +then the descent lemma gives + +```text +L(theta_next) - L(theta) + <= -eta c ||p|| ||v|| + (beta eta^2 / 2) ||v||^2. +``` + +Consequently the update is guaranteed to descend only when + +```text +0 < eta < 2 c ||p|| / (beta ||v||). +``` + +This separates three quantities that a neural-state cosine conflates: + +1. direction (`c`); +2. gain (`||v|| / ||p||`); +3. curvature (`beta`). + +For a local layer, multiplying a hidden teaching vector by the activation gain +and presynaptic activity maps it into parameter space. A positive +`cos(r_l, -g_l)` does not by itself imply a positive cosine after aggregating +examples, parameters, and layers, and even perfect parameter alignment can +increase loss if its gain makes the step exceed the bound. This is why future +diagnostics must include parameter-update cosine, norm ratio, and realized +single-step loss change rather than reporting hidden cosine alone. + +## 3. Somato-dendritic innovation as conditional nuisance removal + +Write apical activity during a task period as + +```text +a = s + M h + epsilon, +``` + +where `s` is the instructional component, `M h` is ordinary traffic predictable +from soma, and `epsilon` is unpredictable noise. In neutral periods `s=0`. +Assume `Sigma_h = E[h h^T]` is nonsingular and `E[epsilon h^T]=0`. Population +least squares on neutral periods yields + +```text +P_neutral = E[a h^T | neutral] Sigma_h^(-1) = M. +``` + +The task-period innovation is then + +```text +r_neutral = a - P_neutral h = s + epsilon. +``` + +If `s`, `M h`, and `epsilon` are mutually uncorrelated, a squared-cosine/SNR +proxy improves from + +```text +raw: S / (S + T + E) +innovation:S / (S + E), +``` + +where `S`, `T`, and `E` are their expected squared norms. The gain is strict +whenever predictable traffic has nonzero power. + +Task-period fitting has a different estimand. Let + +```text +C = E[s h^T] Sigma_h^(-1). +``` + +Then `P_task = M + C` and + +```text +r_task = s - C h + epsilon. +``` + +It removes the soma-predictable projection of the teaching signal itself. This +proves the identification role of neutral periods and also defines the failure +boundary: residualization cannot remove contextual traffic that remains +unpredictable from the chosen per-cell soma statistic, and it can erase credit +when the predictor is fit on task periods. + +For the diagonal normalized-LMS predictor used in the code, an ideal noiseless +cell with neutral relation `a=m h` follows + +```text +p_(t+1) = p_t + eta_P (m - p_t), +traffic residual power fraction = (1 - eta_P)^(2t). +``` + +The relevant timescale is therefore the product of learning rate and number of +neutral updates. Too few neutral updates leave traffic; task-period updates +converge to the wrong coefficient regardless of speed. + +## 4. Hardware-independent complexity + +Let `H` be the number of hidden layers, `K` the perturbation directions per +event, `e` the number of ordinary minibatches between events, `B` the batch +size, and `F_l` the affine work of layer `l`. Define `F=sum_l F_l` and suffix +work `S_l=sum_{j>l} F_j`. + +| method | activity phases / passes | transport | scalar task-loss observations | extra affine work | main transient memory | +|:--|:--|:--|:--|:--|:--| +| BP | forward + reverse | exact transposed Jacobians | ordinary supervised loss | about one reverse pass | saved activations / graph | +| FA | forward + reverse-like serial feedback | fixed random feedback | ordinary loss | `O(F)` feedback work | activations + feedback states | +| DFA | forward + direct feedback | fixed random output maps | ordinary loss | `O(sum_l d_l d_out)` | activations + feedback states | +| learned NP, simultaneous | ordinary forward; one duplicate clean plus `2K` noisy full forwards every `e` batches | learned direct maps; no weight transport | `2KB/e` per ordinary batch | `(1+2K)F/e` | up to `2KB` noisy examples in the vectorized implementation | +| learned NP, layerwise | ordinary forward; clean baseline plus antithetic suffix replays | learned direct maps | `(1+2KH)B/e` | `(F + 2K sum_l S_l)/e`, quadratic in depth for balanced layers | serial suffix replay, no `2K` batch expansion | +| direct NP | layerwise or simultaneous estimator every update | none beyond causal intervention | same formulas with `e=1` | same formulas with `e=1` | estimator dependent | +| canonical EP | free relaxation then nudged relaxation | symmetric energy weights | label injected in nudged phase | `T_free + T_nudged` recurrent steps | persistent/free and nudged states | +| Dual Prop `DP^T` VGG | one inference phase, 17 forward-only layer passes in the author CIFAR protocol | dyadic forward/top-down states | ordinary target loss | 17 layer passes before update | two states per dyadic unit | +| BurstCCN | one burst/event-probability phase with plastic W/Y/Q pathways | separate plastic feedback pathways | ordinary supervised target | architecture-dependent local W/Y/Q updates | event and burst states plus feedback weights | + +The project ledger reports the exact realized loss observations, +forward-equivalent examples, peak allocated memory, and wall time. Asymptotic +labels never replace those measurements. In particular, the frozen direct-NP +C2 diagnosis costs `68.4x` ordinary work, whereas the frozen CIFAR d20 K1/e4 +protocol uses 16x fewer logical queries, 11x less calibration work, and 5.3x +less total work than K16/e4 while retaining its gain over DFA. + +## 5. Executable checks + +Run: + +```bash +python experiments/verify_theory.py +``` + +The verifier checks the exact Rademacher MSE across depth, width, and K; the +quadratic finite-sigma bias law; the smooth-descent step threshold; neutral +predictor convergence across `eta_P` and update counts; and the task-fit +teaching-signal absorption predicted above. diff --git a/experiments/finalize_claims.sh b/experiments/finalize_claims.sh index 104d87c..047e86c 100755 --- a/experiments/finalize_claims.sh +++ b/experiments/finalize_claims.sh @@ -40,7 +40,9 @@ echo "predeclared Pareto/scaling claim gates passed" "$TRAIN_PY" experiments/baseline_smoke.py "$TRAIN_PY" experiments/synthetic_smoke.py +"$TRAIN_PY" experiments/protocol_smoke.py "$TRAIN_PY" experiments/smoke.py +"$TRAIN_PY" experiments/verify_theory.py summary_tmp=$(mktemp /tmp/sdil_audited_tables.XXXXXX) trap 'rm -f "$summary_tmp"' EXIT diff --git a/experiments/verify_theory.py b/experiments/verify_theory.py new file mode 100644 index 0000000..7ee87bc --- /dev/null +++ b/experiments/verify_theory.py @@ -0,0 +1,157 @@ +"""Deterministic numerical checks for the claims in THEORY.md.""" +import itertools +import math + +import numpy as np + + +def simultaneous_mse(rng, depth, width, directions, trials=12000): + gradients = rng.normal(size=(depth, width)) + gradients /= np.linalg.norm(gradients, axis=1, keepdims=True) + estimate = np.zeros((trials, width), dtype=np.float64) + for _ in range(directions): + xi = rng.integers(0, 2, size=(trials, depth, width), dtype=np.int8) + xi = 2.0 * xi - 1.0 + directional = np.einsum("tdw,dw->t", xi, gradients) + estimate += xi[:, 0, :] * directional[:, None] / directions + empirical = np.square(estimate - gradients[0]).sum(axis=1).mean() + total_energy = np.square(gradients).sum() + theoretical = (width * total_energy - np.square(gradients[0]).sum()) / directions + return empirical, theoretical + + +def check_simultaneous_variance(): + print("SIMULTANEOUS RADEMACHER MSE") + print("depth width K empirical theory ratio") + cases = ((2, 8, 1), (4, 8, 1), (8, 8, 1), + (4, 32, 1), (4, 32, 4), (4, 32, 16)) + for index, (depth, width, directions) in enumerate(cases): + empirical, theoretical = simultaneous_mse( + np.random.default_rng(100 + index), depth, width, directions) + ratio = empirical / theoretical + print(f"{depth:5d} {width:5d} {directions:2d} " + f"{empirical:9.4f} {theoretical:9.4f} {ratio:6.3f}") + assert abs(ratio - 1.0) < 0.06 + + +def cubic_loss(z, linear, cubic): + return linear @ z + cubic * np.power(z, 3).sum() / 6.0 + + +def check_sigma_bias(): + # Exhaustive directions remove Monte Carlo error, exposing only the + # centered-difference bias. For this separable cubic it is exactly O(sigma^2). + dim = 8 + z = np.linspace(-0.4, 0.5, dim) + linear = np.linspace(0.2, 0.9, dim) + cubic = 1.7 + true_gradient = linear + 0.5 * cubic * np.square(z) + directions = np.asarray(list(itertools.product((-1.0, 1.0), repeat=dim))) + sigmas = np.asarray((0.005, 0.01, 0.02, 0.04)) + biases = [] + print("\nCENTERED FINITE-DIFFERENCE BIAS") + print("sigma bias_norm") + for sigma in sigmas: + estimates = [] + for xi in directions: + coefficient = (cubic_loss(z + sigma * xi, linear, cubic) + - cubic_loss(z - sigma * xi, linear, cubic)) / (2.0 * sigma) + estimates.append(xi * coefficient) + bias = np.linalg.norm(np.mean(estimates, axis=0) - true_gradient) + biases.append(bias) + print(f"{sigma:5.3f} {bias:.9e}") + slope = np.polyfit(np.log(sigmas), np.log(biases), 1)[0] + print(f"log-log slope={slope:.6f}") + assert abs(slope - 2.0) < 1e-6 + + +def check_descent_threshold(): + rng = np.random.default_rng(301) + dim = 20 + gradient = rng.normal(size=dim) + orthogonal = rng.normal(size=dim) + orthogonal -= orthogonal.dot(gradient) * gradient / gradient.dot(gradient) + orthogonal /= np.linalg.norm(orthogonal) + neg_gradient = -gradient / np.linalg.norm(gradient) + cosine = 0.35 + update = cosine * neg_gradient + math.sqrt(1.0 - cosine ** 2) * orthogonal + beta = 2.3 + threshold = 2.0 * cosine * np.linalg.norm(gradient) / (beta * np.linalg.norm(update)) + + def quadratic_change(step): + return step * gradient.dot(update) + 0.5 * beta * step ** 2 * update.dot(update) + + below = quadratic_change(0.9 * threshold) + above = quadratic_change(1.1 * threshold) + print("\nSMOOTH DESCENT THRESHOLD") + print(f"cos={cosine:.3f} eta_max={threshold:.6f} " + f"delta_below={below:+.6f} delta_above={above:+.6f}") + assert below < 0 < above + + +def check_predictor_timescale(): + rng = np.random.default_rng(401) + examples = 50000 + cells = 12 + h = rng.normal(size=(examples, cells)) + slopes = np.exp(rng.normal(scale=0.25, size=cells)) + traffic = h * slopes + print("\nNEUTRAL PREDICTOR TIMESCALE") + print("eta steps residual_power empirical theory") + for eta in (0.002, 0.01, 0.05): + for steps in (0, 20, 100, 500): + predictor = slopes * (1.0 - (1.0 - eta) ** steps) + residual = traffic - h * predictor + empirical = np.square(residual).mean() / np.square(traffic).mean() + theoretical = (1.0 - eta) ** (2 * steps) + print(f"{eta:5.3f} {steps:5d} {empirical:14.8f} {theoretical:14.8f}") + assert abs(empirical - theoretical) < 1e-12 + + +def squared_cosine(x, y): + return float((x.ravel() @ y.ravel()) ** 2 + / ((x.ravel() @ x.ravel()) * (y.ravel() @ y.ravel()))) + + +def check_innovation_identification(): + rng = np.random.default_rng(501) + examples = 100000 + cells = 8 + h = rng.normal(size=(examples, cells)) + traffic_slopes = np.linspace(0.8, 1.5, cells) + teaching_slopes = np.linspace(0.3, 0.7, cells) + innovation = rng.normal(scale=0.6, size=(examples, cells)) + noise = rng.normal(scale=0.15, size=(examples, cells)) + teaching = h * teaching_slopes + innovation + traffic = h * traffic_slopes + apical = teaching + traffic + noise + + neutral_residual = apical - traffic + task_coeff = (apical * h).mean(axis=0) / np.square(h).mean(axis=0) + task_residual = apical - h * task_coeff + raw_alignment = squared_cosine(apical, teaching) + neutral_alignment = squared_cosine(neutral_residual, teaching) + task_alignment = squared_cosine(task_residual, teaching) + predictable_fraction = np.square(h * teaching_slopes).mean() / np.square(teaching).mean() + retained_fraction = np.square(teaching - h * teaching_slopes).mean() / np.square(teaching).mean() + print("\nINNOVATION IDENTIFICATION") + print(f"squared cosine raw={raw_alignment:.4f} neutral={neutral_alignment:.4f} " + f"task_fit={task_alignment:.4f}") + print(f"teaching predictable={predictable_fraction:.4f} " + f"retained_after_task_fit={retained_fraction:.4f}") + assert neutral_alignment > raw_alignment + 0.15 + assert task_alignment < neutral_alignment - 0.15 + assert abs(retained_fraction - (1.0 - predictable_fraction)) < 0.01 + + +def main(): + check_simultaneous_variance() + check_sigma_bias() + check_descent_threshold() + check_predictor_timescale() + check_innovation_identification() + print("\nALL THEORY CHECKS PASSED") + + +if __name__ == "__main__": + main() |