diff options
| -rw-r--r-- | experiments/verify_theory.py | 49 |
1 files changed, 49 insertions, 0 deletions
diff --git a/experiments/verify_theory.py b/experiments/verify_theory.py index d4325e9..21a7c06 100644 --- a/experiments/verify_theory.py +++ b/experiments/verify_theory.py @@ -239,6 +239,54 @@ def check_residual_coupling_instability(): assert too_negative_radius > 1.0 +def check_dynamic_neutral_projection(): + """A paired neutral fit nulls affine coupling across covariance scales.""" + rng = np.random.default_rng(571) + examples = 128 + cells = 64 + soma = rng.normal(size=(examples, cells)) + coupling = rng.normal(scale=0.04, size=cells) + offset = rng.normal(scale=0.01, size=cells) + neutral_residual = soma * coupling + offset + centered_soma = soma - soma.mean(axis=0) + centered_residual = neutral_residual - neutral_residual.mean(axis=0) + variance = np.square(centered_soma).mean(axis=0) + correction = (centered_soma * centered_residual).mean(axis=0) / variance + remainder = centered_residual - correction * centered_soma + remainder_ratio = (np.linalg.norm(remainder) + / np.linalg.norm(neutral_residual)) + post_slope = ((centered_soma * remainder).mean(axis=0) / variance) + maximum_post_slope = float(np.abs(post_slope).max()) + + eta = 0.1 + momentum = 0.9 + decay = 1e-4 + + def radius(k): + transition = np.asarray(( + (1.0 + eta * k, eta * momentum), + (k, momentum))) + return float(np.abs(np.linalg.eigvals(transition)).max()) + + # With zero residual coupling, covariance drops out and only decay remains. + null_radii = [radius(-decay) for _ in (1e-3, 1.0, 1e3, 1e6)] + # A fixed negative coefficient eventually crosses the lower Jury boundary + # as the nonnegative covariance eigenvalue grows. + fixed_coefficient = -0.03 + high_covariance_k = fixed_coefficient * 1e6 - decay + high_covariance_radius = radius(high_covariance_k) + print("\nDYNAMIC NEUTRAL PROJECTION") + print(f"remainder_ratio={remainder_ratio:.3e} " + f"post_slope={maximum_post_slope:.3e} " + f"rho_null={max(null_radii):.9f} " + f"rho_fixed_high_cov={high_covariance_radius:.3f}") + assert remainder_ratio < 2e-15 + assert maximum_post_slope < 2e-16 + assert max(null_radii) < 1.0 + assert high_covariance_k < -2.0 * (1.0 + momentum) / eta + assert high_covariance_radius > 1.0 + + def check_intermittent_feedback_tracking(): eta_m = 0.1 cadence = 16 @@ -308,6 +356,7 @@ def main(): check_predictor_timescale() check_innovation_identification() check_residual_coupling_instability() + check_dynamic_neutral_projection() check_intermittent_feedback_tracking() check_kolen_pollack_difference_dynamics() print("\nALL THEORY CHECKS PASSED") |
