path: root/THEORY.md

# 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.

### Perturbing the representable convolutional subspace

The failed ResNet-20 run exposed a separate projection problem.  Its
translation-shared feedback field for example `i`, channel `c`, and spatial
site `s` is

```text
p_ics = b_ic + t_ics v_ic,       t_ics = tanh(h_ics).
```

The legacy estimator perturbs every `h_ics` independently, constructs a full
unit target, and only then averages its delta-rule update into the shared
`b/v` parameterization.  At K1, almost all sampled energy can lie outside this
two-dimensional per-channel subspace.  SDIL-v2 instead draws independent
Rademacher coefficients `z^b_ic,z^v_ic` and intervenes along

```text
xi_ics = (z^b_ic + t_ics z^v_ic) / sqrt(2).
```

Let `q_ics=-dL/dh_ics`, let `S` be the number of spatial sites, and let `D` be
the centered directional loss derivative.  Then

```text
qhat^b_ic = -sqrt(2) D z^b_ic / S,
qhat^v_ic = -sqrt(2) D z^v_ic / S
```

are unbiased for the two negative-gradient moments

```text
qbar^b_ic = mean_s q_ics,
qbar^v_ic = mean_s (q_ics t_ics).
```

All other examples, channels, and layers contribute zero-mean cross terms.
For BatchNorm, `D` is the derivative of the summed minibatch loss, so the same
identity targets each example's contribution including the true cross-example
Jacobian coupling.

This does not mistake moments for expansion coefficients.  The predicted
moments of `p=b+t v` are

```text
pbar^b = b + mean(t) v,
pbar^v = mean(t) b + mean(t^2) v.
```

Consequently updating the two vectorizer pathways with
`qhat^b-pbar^b` and `qhat^v-pbar^v` is exactly the expected gradient of the
original full-field squared prediction error.  It uses the same two causal
loss queries per direction as unit perturbation, but avoids estimating
unrepresentable spatial components.  It remains a variance reduction, not a
guarantee of useful feedback or stable end-to-end optimization; those require
a separately frozen empirical gate.

### Perturbing the vectorizer parameter subspace

The channel-moment estimator still samples one independent base/gate
coefficient for every example and channel and only afterward regresses those
noisy coefficients on the output context.  The quantity the local delta rule
actually needs is much smaller.  For

```text
p_ics = (A c_i)_c + t_ics (G c_i)_c
```

define the normalized prediction objective

```text
J(A,G) = (1 / 2B) sum_i mean_s ||p_i - q_i||^2,
q_i = -d ell_i / dh_i.
```

Its causal target terms are matrices, not per-example coefficient fields:

```text
U_A = (1/B) sum_i mean_s(q_ics) c_i^T,
U_G = (1/B) sum_i mean_s(q_ics t_ics) c_i^T.
```

Draw independent Rademacher matrices `Xi_A,Xi_G` with the same shapes as
`A,G`, independently across layers, and intervene along

```text
d_ics = ((Xi_A c_i)_c + t_ics (Xi_G c_i)_c) / sqrt(2).
```

Let `D` be the centered derivative of the summed minibatch loss under all
layer interventions. Since `D=-sum_ics q_ics d_ics`, independence gives

```text
Uhat_A = -sqrt(2) D Xi_A / (B S),     E[Uhat_A] = U_A,
Uhat_G = -sqrt(2) D Xi_G / (B S),     E[Uhat_G] = U_G.
```

Cross-layer and A/G cross terms vanish in expectation. The prediction terms
`(1/B) sum mean(p)c^T` and `(1/B) sum mean(pt)c^T` are locally available
exactly, so subtracting them from `Uhat_A,Uhat_G` gives an unbiased stochastic
gradient of `J`. BatchNorm changes neither identity: `D` is explicitly the
summed minibatch derivative and therefore includes its true cross-example
Jacobian.

This estimator perturbs only the shared mapping that will survive
amortization. It uses the same two loss queries per direction as the earlier
estimators. It does not uniformly dominate them for every context covariance
or batch size, but at the frozen batch-128 synthetic audit its one-direction
MSE is `0.03381x` the coefficient-then-regress estimator's MSE. The executable
mean has cosine `0.96928` and norm ratio `1.04888` to the exact A/G target.
Those mechanics do not establish task utility; a separately frozen gate is
still required.

## 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

The regression residual has a stronger interpretation than the linear model
alone suggests. Let `F_h` be the information available to a soma-only
predictor in a neutral period, let `n` denote ordinary apical traffic, and set

```text
m(h) = E[n | F_h].
```

For every square-integrable soma-measurable predictor `g(h)`, conditional
expectation gives the orthogonality and Pythagorean identities

```text
E <n - m(h), g(h)> = 0,

E ||n - g(h)||^2
  = E ||n - m(h)||^2 + E ||m(h) - g(h)||^2.
```

Thus `n-m(h)` is the unique minimum-power innovation after subtracting any
function of the available somatic statistic. This is an `L2` projection
statement, not a Gaussian or linear assumption. If task-period apical activity
is `a=s+n` and the neutral traffic law is invariant, neutral fitting yields

```text
r = a - m(h) = s + n - m(h).
```

It removes the largest soma-predictable nuisance component without projecting
out any part of the task-only instruction `s`. The implemented diagonal affine
predictor is the corresponding optimal projection only within each cell's
restricted span `{1,h_i}`; it cannot remove traffic predictable only from
population context or unobserved state. This restriction is consistent with
the failed broad endogenous-traffic gate and prevents the theorem from being
used to claim arbitrary contextual denoising.

The norm-matched raw control isolates direction from magnitude. For every
positive per-example scalar `alpha`,

```text
cos(alpha a, -g) = cos(a, -g).
```

Choosing `alpha=||r||/||a||` gives raw feedback exactly the innovation's norm
but cannot rotate it toward the descent direction at the same network state.
Across training, rescaling can of course change the subsequent trajectory;
this is why the matched control is still trained end to end rather than being
replaced by the algebraic observation.

The finite-dimensional special case used for the main experiment writes
neutral traffic as a linear somatic component plus unpredictable noise.

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; the
conditional-projection Pythagorean identity and norm-matching direction
invariance; neutral predictor convergence across `eta_P` and update counts;
and the task-fit teaching-signal absorption predicted above.