LTEE B2 — Anderson Disorder Analogy for Fitness Dynamics

Paper: Wiser, Ribeck & Lenski (2013) “Long-Term Dynamics of Adaptation in Asexual Populations” Science 342, 1364–1367

LTEE GuideStone Queue: B2 (hotSpring assignment)

What we reproduce: The fitness trajectory of the E. coli Long-Term Evolution Experiment follows a power-law with no asymptote. We apply Anderson localization diagnostics (level spacing ratio, Lyapunov exponent) to fitness increment sequences, testing whether the dynamics exhibit localization–delocalization transitions analogous to the Anderson model.

Existing infrastructure:

• Papers 14–20: Full Anderson 1D/2D/3D localization + spectral theory
• control/spectral_theory/scripts/spectral_control.py: Python baseline
• barracuda::spectral: Rust parity with level_spacing_ratio, GOE/Poisson

Bridge: This notebook extends the physical Anderson model to biological fitness trajectories — same spectral diagnostics, different substrate.

Tier 1 Python baseline. Rust Tier 2 target: validate_ltee_anderson scenario. Feeds lithoSpore module 7 (anderson).

Physics–Biology Mapping

Key Prediction

If fitness dynamics are analogous to Anderson localization:

• Early LTEE (generations 0–10,000): fitness increments are correlated → $\langle r \rangle$ near GOE (extended/exploring)
• Late LTEE (generations 40,000+): increments become uncorrelated → $\langle r \rangle$ drifts toward Poisson (localized/diminishing returns)
• The power-law (no plateau) means full Poisson localization is never reached — the system is always at a critical point

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import json
from pathlib import Path

POISSON_R = 2.0 * np.log(2.0) - 1.0  # 0.3863
GOE_R = 0.531

C_PASS = '#2ecc71'
C_FAIL = '#e74c3c'
C_INFO = '#3498db'
C_RUST = '#1abc9c'
C_PYTHON = '#f39c12'
C_GPU = '#9b59b6'
C_BIO = '#e67e22'

rng = np.random.default_rng(42)

print(f"Poisson <r> = {POISSON_R:.4f}")
print(f"GOE <r> = {GOE_R:.3f}")

1. Wiser Power-Law Fitness Model

Wiser et al. fit the LTEE fitness data to a power-law:

$$w(t) = (1 + 2 \alpha t)^{\beta}$$

where $\alpha$ and $\beta$ are fitted from 12 replicate populations over 50,000 generations. The key finding: $\beta > 0$ with no asymptote, meaning adaptation continues indefinitely.

Published parameter estimates (Wiser et al. 2013, Table S2):

• $\alpha \approx 6.2 \times 10^{-4}$ per generation
• $\beta \approx 0.056$ (diminishing returns exponent)

# Published Wiser et al. parameters (Table S2, mean across populations)
ALPHA = 6.2e-4
BETA = 0.056

# LTEE sampling points (generations where fitness was measured)
GENERATIONS = np.array([
    500, 1000, 1500, 2000, 5000, 10000, 15000, 20000,
    30000, 40000, 50000
])

def wiser_fitness(t, alpha=ALPHA, beta=BETA):
    """Power-law fitness model from Wiser et al. 2013."""
    return (1.0 + 2.0 * alpha * t) ** beta

# Generate fine-grained trajectory for analysis
t_fine = np.linspace(500, 50000, 1000)
w_fine = wiser_fitness(t_fine)

# Fitness at measured time points
w_measured = wiser_fitness(GENERATIONS)

print(f"Fitness at gen 500: {wiser_fitness(500):.4f}")
print(f"Fitness at gen 50000: {wiser_fitness(50000):.4f}")
print(f"Relative gain: {wiser_fitness(50000) / wiser_fitness(500):.2f}x")

# Replicate populations with noise (12 populations, as in LTEE)
N_POPS = 12
NOISE_SIGMA = 0.02  # Measurement noise from competitive fitness assays

pop_trajectories = []
for pop in range(N_POPS):
    # Each population has slightly different alpha/beta + measurement noise
    a_pop = ALPHA * (1.0 + 0.15 * rng.standard_normal())
    b_pop = BETA * (1.0 + 0.10 * rng.standard_normal())
    w_pop = wiser_fitness(GENERATIONS, a_pop, b_pop)
    w_pop += NOISE_SIGMA * rng.standard_normal(len(GENERATIONS))
    pop_trajectories.append(w_pop)

pop_trajectories = np.array(pop_trajectories)

fig, ax = plt.subplots(figsize=(10, 6))
for i, traj in enumerate(pop_trajectories):
    ax.plot(GENERATIONS / 1000, traj, 'o-', alpha=0.4, markersize=3,
            label=f'Pop {i+1}' if i < 3 else None)
ax.plot(t_fine / 1000, w_fine, 'k-', lw=2, label='Mean power-law')
ax.set_xlabel('Generations (thousands)')
ax.set_ylabel('Relative fitness w(t)')
ax.set_title(f'Wiser et al. 2013: w(t) = (1 + 2αt)^β\n'
             f'α = {ALPHA:.1e}, β = {BETA:.3f} — no plateau')
ax.legend(ncol=2)
ax.grid(True, alpha=0.3)
fig.savefig('/tmp/hotspring_ltee_b2_fitness_trajectories.png', dpi=150,
            bbox_inches='tight')
plt.close(fig)
print('Saved fitness trajectories plot')

2. Fitness Increments as Disorder Potential

The Anderson analogy maps fitness increments $\Delta w_i$ to the on-site disorder potential $V_i$ in the Anderson Hamiltonian. If we construct a tight-binding Hamiltonian where the diagonal elements are the fitness increments, we can apply standard spectral diagnostics.

from scipy import sparse
from scipy.sparse.linalg import eigsh

def level_spacing_ratio(eigenvalues):
    """Compute mean level spacing ratio <r> for sorted eigenvalues.
    
    GOE: <r> ~ 0.531 (extended/correlated)
    Poisson: <r> ~ 0.386 (localized/uncorrelated)
    """
    evals = np.sort(eigenvalues)
    spacings = np.diff(evals)
    spacings = spacings[spacings > 1e-14]
    if len(spacings) < 3:
        return float('nan')
    ratios = []
    for i in range(len(spacings) - 1):
        s_n = spacings[i]
        s_next = spacings[i + 1]
        ratios.append(min(s_n, s_next) / max(s_n, s_next))
    return np.mean(ratios)


def fitness_anderson_hamiltonian(fitness_increments, hopping=1.0):
    """Build Anderson-like Hamiltonian from fitness increments.
    
    H = -t * (sum |i><i+1| + h.c.) + sum V_i |i><i|
    where V_i = fitness_increments[i] (normalized to unit variance)
    """
    n = len(fitness_increments)
    # Normalize increments to control effective disorder
    v = fitness_increments.copy()
    if np.std(v) > 1e-14:
        v = (v - np.mean(v)) / np.std(v)
    
    diag = v
    off_diag = -hopping * np.ones(n - 1)
    H = np.diag(diag) + np.diag(off_diag, 1) + np.diag(off_diag, -1)
    return H


# Compute fitness increments from the mean trajectory
increments = np.diff(w_fine)
print(f"Fitness increments: {len(increments)} points")
print(f"Mean increment: {np.mean(increments):.6f}")
print(f"Std increment: {np.std(increments):.6f}")
print(f"Increment ratio (last/first): {increments[-1]/increments[0]:.4f}")
print(f"  → Diminishing returns: each step contributes less")

# Build Anderson Hamiltonian from fitness increments and analyze
H_fitness = fitness_anderson_hamiltonian(increments)
evals_fitness = np.linalg.eigvalsh(H_fitness)

r_fitness = level_spacing_ratio(evals_fitness)
print(f"Level spacing ratio for fitness trajectory:")
print(f"  <r> = {r_fitness:.4f}")
print(f"  GOE reference: {GOE_R:.4f}")
print(f"  Poisson reference: {POISSON_R:.4f}")
print()

if r_fitness > (GOE_R + POISSON_R) / 2:
    phase = 'EXTENDED (GOE-like)'
elif r_fitness < POISSON_R + 0.03:
    phase = 'LOCALIZED (Poisson-like)'
else:
    phase = 'CRITICAL (intermediate)'
print(f"  Phase: {phase}")

3. Sliding-Window Localization Analysis

The key test: does $\langle r \rangle$ evolve over the trajectory? If fitness dynamics undergo a localization transition, we expect $\langle r \rangle$ to decrease from GOE toward Poisson as the system moves from rapid adaptation to diminishing returns.

# Sliding window analysis: <r> as a function of evolutionary time
WINDOW = 100  # Eigenvalues per window
STRIDE = 50

window_centers = []
window_r_values = []

for start in range(0, len(increments) - WINDOW, STRIDE):
    window = increments[start:start + WINDOW]
    H_w = fitness_anderson_hamiltonian(window)
    evals_w = np.linalg.eigvalsh(H_w)
    r_w = level_spacing_ratio(evals_w)
    
    # Map window position back to generation number
    gen_center = t_fine[start + WINDOW // 2]
    window_centers.append(gen_center)
    window_r_values.append(r_w)

window_centers = np.array(window_centers)
window_r_values = np.array(window_r_values)

fig, ax = plt.subplots(figsize=(10, 5))
ax.plot(window_centers / 1000, window_r_values, 'o-', color=C_BIO,
        markersize=4, label='Fitness <r>')
ax.axhline(GOE_R, color=C_RUST, ls='--', label=f'GOE = {GOE_R:.3f}')
ax.axhline(POISSON_R, color=C_PYTHON, ls='--', label=f'Poisson = {POISSON_R:.4f}')
ax.set_xlabel('Generation (thousands)')
ax.set_ylabel('Level spacing ratio <r>')
ax.set_title('Anderson Localization Diagnostic on Fitness Trajectory')
ax.legend()
ax.grid(True, alpha=0.3)
ax.set_ylim(0.3, 0.6)
fig.savefig('/tmp/hotspring_ltee_b2_sliding_r.png', dpi=150,
            bbox_inches='tight')
plt.close(fig)
print('Saved sliding-window <r> plot')

4. Anderson Reference: Physical System Comparison

Compare the fitness-derived Hamiltonian against a true Anderson model with the same dimensionality. This validates that our diagnostics detect known transitions and provides a calibration baseline.

def anderson_1d(n, W):
    """Standard 1D Anderson Hamiltonian with uniform disorder [-W/2, W/2]."""
    diag = W * (rng.random(n) - 0.5)
    off = -np.ones(n - 1)
    return np.diag(diag) + np.diag(off, 1) + np.diag(off, -1)

N = 200
disorder_strengths = [0.5, 1.0, 2.0, 4.0, 8.0, 16.0]
r_anderson = []

for W in disorder_strengths:
    rs = []
    for trial in range(20):
        H_a = anderson_1d(N, W)
        evals_a = np.linalg.eigvalsh(H_a)
        rs.append(level_spacing_ratio(evals_a))
    r_anderson.append((W, np.mean(rs), np.std(rs)))
    print(f"W = {W:5.1f}  <r> = {np.mean(rs):.4f} ± {np.std(rs):.4f}")

print(f"\nFitness trajectory <r> = {r_fitness:.4f}")

# Determine effective disorder strength
closest_W = min(r_anderson, key=lambda x: abs(x[1] - r_fitness))
print(f"Closest Anderson W: {closest_W[0]:.1f} (<r> = {closest_W[1]:.4f})")

# Comparison plot: Anderson <r> vs W with fitness overlay
fig, ax = plt.subplots(figsize=(10, 5))

Ws = [x[0] for x in r_anderson]
rs = [x[1] for x in r_anderson]
errs = [x[2] for x in r_anderson]

ax.errorbar(Ws, rs, yerr=errs, fmt='s-', color=C_INFO, label='1D Anderson',
            capsize=3, markersize=6)
ax.axhline(GOE_R, color=C_RUST, ls='--', alpha=0.5, label=f'GOE = {GOE_R:.3f}')
ax.axhline(POISSON_R, color=C_PYTHON, ls='--', alpha=0.5,
           label=f'Poisson = {POISSON_R:.4f}')
ax.axhline(r_fitness, color=C_BIO, ls='-', lw=2, alpha=0.8,
           label=f'Fitness <r> = {r_fitness:.4f}')

ax.set_xlabel('Disorder strength W')
ax.set_ylabel('<r>')
ax.set_title('Level Spacing Ratio: Anderson 1D vs LTEE Fitness Trajectory')
ax.set_xscale('log')
ax.legend()
ax.grid(True, alpha=0.3)
fig.savefig('/tmp/hotspring_ltee_b2_anderson_comparison.png', dpi=150,
            bbox_inches='tight')
plt.close(fig)
print('Saved Anderson comparison plot')

5. Population-Level Variance Analysis

Wiser et al. showed that replicate trajectories diverge. In the Anderson analogy, this corresponds to sensitivity to the specific disorder realization — each population samples a different random potential.

# Analyze each population's level spacing ratio
pop_r_values = []
for i, traj in enumerate(pop_trajectories):
    inc = np.diff(traj)
    H_p = fitness_anderson_hamiltonian(inc)
    evals_p = np.linalg.eigvalsh(H_p)
    r_p = level_spacing_ratio(evals_p)
    pop_r_values.append(r_p)
    print(f"Pop {i+1:2d}: <r> = {r_p:.4f}")

pop_r_values = np.array(pop_r_values)
print(f"\nMean <r> across populations: {np.mean(pop_r_values):.4f} ± {np.std(pop_r_values):.4f}")
print(f"Range: [{np.min(pop_r_values):.4f}, {np.max(pop_r_values):.4f}]")

6. Expected Values for lithoSpore

Produce the frozen JSON that lithoSpore module 7 (anderson) will absorb as validation targets.

expected = {
    "paper": "Wiser et al. 2013",
    "ltee_queue_id": "B2",
    "spring": "hotSpring",
    "model": {
        "type": "power_law",
        "formula": "w(t) = (1 + 2*alpha*t)^beta",
        "alpha": ALPHA,
        "beta": BETA,
        "source": "Wiser et al. 2013 Table S2"
    },
    "fitness_values": {
        "gen_500": float(wiser_fitness(500)),
        "gen_5000": float(wiser_fitness(5000)),
        "gen_10000": float(wiser_fitness(10000)),
        "gen_50000": float(wiser_fitness(50000))
    },
    "anderson_diagnostics": {
        "full_trajectory_r": float(r_fitness),
        "goe_reference": GOE_R,
        "poisson_reference": float(POISSON_R),
        "population_mean_r": float(np.mean(pop_r_values)),
        "population_std_r": float(np.std(pop_r_values)),
        "effective_disorder_W": float(closest_W[0]),
        "sliding_window": {
            "window_size": WINDOW,
            "stride": STRIDE,
            "n_windows": len(window_r_values),
            "r_early": float(np.mean(window_r_values[:3])),
            "r_late": float(np.mean(window_r_values[-3:]))
        }
    },
    "validation_checks": [
        {"name": "power_law_no_plateau", "check": "fitness at 50k > fitness at 10k",
         "expected": True},
        {"name": "diminishing_returns", "check": "increment ratio last/first < 1",
         "expected_bound": 1.0},
        {"name": "r_between_goe_poisson", "check": "Poisson < <r> < GOE",
         "expected_range": [float(POISSON_R), GOE_R]},
        {"name": "population_variance", "check": "std(<r>) > 0",
         "expected": True},
        {"name": "n_populations", "check": "12 replicate populations analyzed",
         "expected": N_POPS}
    ],
    "provenance": {
        "notebook": "notebooks/papers/13-ltee-anderson-fitness.ipynb",
        "date": "2026-05-11",
        "python": "numpy + scipy + matplotlib"
    }
}

# Write expected values JSON
out_dir = Path('..') / '..' / 'experiments' / 'results' / 'ltee'
out_dir.mkdir(parents=True, exist_ok=True)
out_path = out_dir / 'ltee_b2_anderson_expected.json'
with open(out_path, 'w') as f:
    json.dump(expected, f, indent=2)
print(f'Wrote expected values to {out_path}')
print(json.dumps(expected, indent=2))

Summary

Tier 1 status: Python baseline COMPLETE.

Next steps (Tier 2):

• Implement validate_ltee_anderson Rust scenario using barracuda::spectral
• Wire expected values into foundation Thread 7 (anderson localization)
• Compare against LTEE public data (when available via lithoSpore ingest)

Provenance: This notebook uses the power-law model from Wiser et al. 2013. The Anderson analogy is our contribution — applying the spectral diagnostics from Papers 14–20 to biological fitness trajectories.

hotSpring LTEE B2 Tier 1 — Anderson Disorder Analogy for Fitness Dynamics
Feeds: lithoSpore module 7 (anderson), foundation Thread 7