My July 2025 introduction to SVDA focused on understanding attention. Here I want to explore a practical question: once we can inspect a model’s internal directions, how do we decide which ones to keep—and how do we check that an implementation change does what we intended?
Research context
Our paper, Prescriptive SVD-Inspired Attention via Spectral Energy Retention, by Vasileios Arampatzakis, Vasileios Sevetlidis and George Pavlidis, was first posted to arXiv on 21 September 2026 and identifies its venue as TMLR 2026. It introduces interventions on SVDA’s learned diagonal spectrum. Its energy-retention prescription selects directions within each head. The paper distinguishes masking from dimension reduction: preserving the original query/key normalization gives equivalence, whereas normalizing after selection changes the operator. Its reported compression results concern the score-forming pathway; they should not be read as equivalent reductions in total model size or measured inference time.
The example below is an independent, synthetic calculation. It does not reproduce the paper’s trained models or benchmarks.
Energy is a budget, not a coefficient cutoff
Take four coefficients:
Their squared magnitudes are 16, 4, 1 and 0.25, with total 21.25. Keeping the first two preserves 20/21.25, or about 94.12% of the energy. Keeping only the first preserves about 75.29%. A 90% budget therefore requires two coordinates in this example.
Notice what the percentage means. It does not say “remove coefficients smaller than 0.9,” and it does not promise “keep 90% of predictive performance.” Multiplying all four coefficients by 100 leaves these energy fractions unchanged. It changes the scale of the scores, however, so equal relative budgets need not produce equal changes in model behavior.
Here is a minimal selection function for this nonzero vector:
import numpy as np
sigma = np.array([4., 2., 1., .5])
rho = .90
energy = (sigma / np.max(np.abs(sigma))) ** 2
order = np.argsort(-energy, kind="stable")
cumulative = np.cumsum(energy[order])
width = np.searchsorted(cumulative, rho * cumulative[-1]) + 1
keep = np.sort(order[:width])
print(keep) # [0 1]
The accompanying complete example validates the input and handles an all-zero spectrum separately: a retained-energy fraction is undefined when its denominator is zero.
A small implementation choice changes the experiment
For a matrix calculation, let Q and K each contain six unit-length rows with four coordinates. Compare these two operations:
- Set the discarded coefficients to zero and multiply using the original rows.
- Select the retained coordinates from those same normalized rows and multiply the smaller arrays.
With identical scaling, both calculate the same sum. Selecting coordinates simply omits terms whose coefficients are already zero. But the selected rows generally have length less than one. Normalizing them again introduces a new transformation.
masked = (q * masked_sigma) @ k.T / np.sqrt(4)
selected = (q[:, keep] * sigma[keep]) @ k[:, keep].T / np.sqrt(4)
assert np.allclose(masked, selected)
Here q and k were normalized before selecting coordinates. Replacing the denominator with the square root of the retained width also changes the computation. These are useful alternative experiments, provided they are named and evaluated separately.
What to measure next
The executable example checks score differences and the resulting softmax distributions. A trained model needs additional checks: held-out task performance, variation across seeds, parameter counts, and actual latency on the intended hardware. A dense tensor full of zeros can produce exactly the desired mathematical intervention while providing little practical speed benefit.
For an implementation review, I would retain a reference output from the unmodified checkpoint, verify the masked calculation, and then compare the smaller implementation against that reference on the same inputs. This catches normalization mistakes before they become explanations for an apparent accuracy improvement or regression.
A useful experiment ends with evidence about both behavior and implementation cost. A percentage of discarded coordinates alone answers neither question.
Complete reproducible example
These files contain synthetic teaching data and calculations. They do not reproduce a trained model or a field deployment.
spectral_retention.py
"""Synthetic score-path experiment; not a benchmark reproduction."""
import json
import numpy as np
def retained_indices(sigma, rho):
sigma = np.asarray(sigma, dtype=float)
if sigma.ndim != 1 or sigma.size == 0 or not np.isfinite(sigma).all():
raise ValueError('sigma must be a nonempty finite vector')
if not 0 < rho <= 1:
raise ValueError('rho must be in (0, 1]')
# Rescale before squaring to avoid overflow; ordering and ratios unchanged.
scale = np.max(np.abs(sigma))
if scale == 0:
return np.array([], dtype=int)
energy = (sigma / scale) ** 2
order = np.argsort(-energy, kind='stable')
cumulative = np.cumsum(energy[order])
target = rho * cumulative[-1]
width = min(np.searchsorted(cumulative, target, side='left') + 1, sigma.size)
return np.sort(order[:width])
def normalize_rows(x):
norms = np.linalg.norm(x, axis=1, keepdims=True)
if np.any(norms == 0):
raise ValueError('zero rows cannot be normalized')
return x / norms
def softmax(x):
e = np.exp(x - x.max(axis=1, keepdims=True))
return e / e.sum(axis=1, keepdims=True)
def run():
rng = np.random.default_rng(7)
q, k = normalize_rows(rng.normal(size=(6, 4))), normalize_rows(rng.normal(size=(6, 4)))
sigma = np.array([4., 2., 1., .5])
keep = retained_indices(sigma, .90)
masked_sigma = np.zeros_like(sigma)
masked_sigma[keep] = sigma[keep]
original = (q * sigma) @ k.T / np.sqrt(4)
masked = (q * masked_sigma) @ k.T / np.sqrt(4)
selected = (q[:, keep] * sigma[keep]) @ k[:, keep].T / np.sqrt(4)
renormalized = (normalize_rows(q[:, keep]) * sigma[keep]) @ normalize_rows(k[:, keep]).T / np.sqrt(4)
exact_error = float(np.max(np.abs(masked-selected)))
score_error = float(np.max(np.abs(original-masked)))
bound = float(np.max(np.abs(sigma-masked_sigma)) / np.sqrt(4))
attention_error = float(np.max(np.sum(np.abs(softmax(original)-softmax(masked)), axis=1)))
assert keep.tolist() == [0, 1]
assert exact_error < 1e-12
assert score_error <= bound + 1e-12
assert attention_error <= score_error + 1e-12
assert retained_indices([0, 0], .9).size == 0
assert retained_indices([1, -1, 1], 1.).size == 3
assert np.array_equal(retained_indices(sigma, .9), retained_indices(sigma*100, .9))
return {'retained_indices':keep.tolist(), 'retained_energy_fraction':float(np.sum(sigma[keep]**2)/np.sum(sigma**2)),
'masked_vs_selected_max_error':exact_error, 'renormalized_vs_masked_max_error':float(np.max(np.abs(renormalized-masked))),
'score_error':score_error, 'score_bound':bound, 'attention_row_l1_error':attention_error}
if __name__ == '__main__':
print(json.dumps(run(), indent=2))
