Eigenstructure and Persistent Directions

Part I

Persistence is repeated reweighting

A matrix usually changes both the length and direction of a vector. An eigenvector is exceptional because the transformation changes only its scale, or its scale and sign. Gilbert Strang describes eigenvectors as the special vectors for which Ax stays parallel to x, and connects eigenstructure directly to dynamic systems. [1] The definition is compact:

The important word is not “special”. It is “invariant”. Once a component lies along an eigenvector, repeated application cannot leak it into another eigendirection in the diagonalizable idealisation. The component follows a scalar recurrence. After k steps, its multiplier is λᵢᵏ. A difficult multivariate evolution has been decomposed into independent one-dimensional growth laws.

Why invariant directions matter

Suppose the eigenvectors form a basis, so the starting state can be written as a weighted sum:

Now order the eigenvalues by magnitude, with |λ₁| > |λ₂| ≥ …. Divide the whole expression by λ₁ᵏ. Every competing coefficient contains a factor (λᵢ/λ₁)ᵏ. Those factors approach zero when their moduli are below one. Provided c₁ ≠ 0, the normalised state approaches the line through v₁.

This is why “persistent” is relative. In the press, both modes vanish in absolute magnitude. Yet v₁ persists relative to v₂. In a growing system, both may expand, but the larger modulus wins. In a Markov chain, the stationary mode has eigenvalue one while the others decay. The same algebra describes absolute decay, absolute growth and convergence to a stable composition.

The ledger also separates three questions that are often collapsed. The dominant eigenvector answers which shape survives. The dominant eigenvalue answers how that shape scales. The ratio between leading eigenvalue moduli answers how quickly other shapes become negligible. A stability statement that names only one of these is usually underspecified.

For two clean modes, the leading error behaves approximately like

Markov-chain theory gives this separation an operational name. For a reversible chain, the absolute spectral gap is one minus the largest non-stationary eigenvalue modulus, and its reciprocal defines a relaxation time. [5] The exact mixing bound needs more than one number, but the gap explains why two systems with the same stationary state can forget their starting conditions at very different speeds.

Part II

Power iteration turns persistence into an instrument

Once repeated multiplication selects a direction, the mechanism itself becomes an algorithm. Choose a starting vector, multiply by the matrix, normalise to avoid overflow or underflow, and repeat. The operation never asks for the characteristic polynomial. It only needs a matrix-vector product. That makes the method useful when the matrix is sparse, implicit or too large to manipulate as a dense object.

The power method appears in authoritative numerical-eigenvalue templates as the simplest single-vector iteration. [2] Modern libraries usually expose stronger general solvers. NumPy’s dense eig routine returns all right eigenvectors of a square array, while SciPy’s sparse eigs interface finds selected eigenpairs for an array, sparse matrix or matrix-as-operator. [3] [4] The conceptual value of power iteration remains unusually high because every line of the algorithm mirrors the causal explanation.

Power iteration is not searching over directions. It is letting the operator amplify its own preference. That distinction matters when reading the output. The result is a property of the operator, the start and the stopping rule. It is not automatically the most meaningful direction for the application, and it is not guaranteed to be the direction of greatest one-step amplification. The latter is a singular-vector question.

Positive case and negative control

Begin the press experiment at (3, 2). This vector contains both eigenmodes. Twelve normalised multiplications leave an angle of roughly 0.0004 degrees to the dominant line, and the Rayleigh quotient is essentially 0.9. The method succeeds because the dominant component is present and separated.

Now begin exactly at v₂. Every multiplication remains on v₂. Normalisation removes its shrinking magnitude, so the algorithm reports eigenvalue 0.4 forever. The dominant mode exists, but the experiment never excites it. In exact arithmetic, the method cannot invent a missing component.

The negative control prevents a common verbal inflation: “power iteration finds the dominant eigenvector from any start.” The defensible claim requires a non-zero projection on the dominant eigenspace, a unique dominant modulus, adequate numerical precision and enough iterations. A random start makes exact orthogonality unlikely in many ordinary settings. It does not cure a narrow gap, poor conditioning, finite-precision loss or a mis-specified operator.

The executable experiment

The following NumPy programme reproduces the positive case, the missing-excitation negative case and the finite-horizon non-normal counterexample used later. It operates entirely on synthetic 2 × 2 matrices. Its assertions are part of the argument: a claim about persistence is accepted only when the expected mode appears, the negative control stays negative, and the stable-spectrum counterexample exhibits transient amplification.

import numpy as np

A = np.array([
    [11 / 15, 1 / 3],
    [1 / 6, 17 / 30],
], dtype=float)

V1 = np.array([2.0, 1.0])
V2 = np.array([1.0, -1.0])


def unit(vector: np.ndarray) -> np.ndarray:
    norm = np.linalg.norm(vector)
    if np.isclose(norm, 0.0):
        raise ValueError("The iteration produced the zero vector.")
    return vector / norm


def line_angle_degrees(left: np.ndarray, right: np.ndarray) -> float:
    """Angle between lines, so v and -v count as the same direction."""
    cosine = abs(float(unit(left) @ unit(right)))
    return float(np.degrees(np.arccos(np.clip(cosine, -1.0, 1.0))))


def power_iteration(matrix: np.ndarray, start: np.ndarray, steps: int = 12):
    vector = unit(start.astype(float))
    history = [vector.copy()]
    for _ in range(steps):
        vector = unit(matrix @ vector)
        history.append(vector.copy())
    eigenvalue = float(vector @ matrix @ vector)  # Rayleigh quotient
    return vector, eigenvalue, np.array(history)


# Positive case: the start contains a non-zero component along V1.
pos_vector, pos_value, pos_history = power_iteration(A, np.array([3.0, 2.0]))
pos_angle = line_angle_degrees(pos_vector, V1)

# Negative control: the start is exactly V2, so V1 is never excited.
neg_vector, neg_value, neg_history = power_iteration(A, V2)
neg_angle = line_angle_degrees(neg_vector, V1)

# Finite-horizon counterexample: same stable eigenvalues, different geometry.
normal = np.diag([0.8, 0.6])
non_normal = np.array([[0.8, 3.0], [0.0, 0.6]])
probe = np.array([0.0, 1.0])
normal_norms = [np.linalg.norm(np.linalg.matrix_power(normal, k) @ probe)
                for k in range(21)]
non_normal_norms = [np.linalg.norm(np.linalg.matrix_power(non_normal, k) @ probe)
                    for k in range(21)]

def finite_horizon_gain(matrix: np.ndarray, horizon: int):
    """Worst Euclidean amplification over every unit start through horizon."""
    gains = [np.linalg.svd(np.linalg.matrix_power(matrix, k),
                           compute_uv=False)[0]
             for k in range(horizon + 1)]
    return float(max(gains)), int(np.argmax(gains))

normal_gain, normal_gain_step = finite_horizon_gain(normal, horizon=8)
non_normal_gain, non_normal_gain_step = finite_horizon_gain(non_normal, horizon=8)

print(f"positive: eigenvalue={pos_value:.6f}, angle={pos_angle:.6f} degrees")
print(f"negative: eigenvalue={neg_value:.6f}, angle-to-v1={neg_angle:.6f} degrees")
print(f"normal peak={max(normal_norms):.6f}")
print(f"non-normal peak={max(non_normal_norms):.6f} at step {np.argmax(non_normal_norms)}")
print(f"normal worst-case gain={normal_gain:.6f} at step {normal_gain_step}")
print(f"non-normal worst-case gain={non_normal_gain:.6f} at step {non_normal_gain_step}")

assert abs(pos_value - 0.9) < 1e-5
assert pos_angle < 0.001
assert abs(neg_value - 0.4) < 1e-12
assert neg_angle > 70.0
assert max(normal_norms) <= 1.0
assert max(non_normal_norms) > 4.4
assert abs(normal_gain - 1.0) < 1e-12
assert non_normal_gain > 4.47

Expected output

positive: eigenvalue=0.900001, angle=0.000408 degrees
negative: eigenvalue=0.400000, angle-to-v1=71.565051 degrees
normal peak=1.000000
non-normal peak=4.445251 at step 3
normal worst-case gain=1.000000 at step 0
non-normal worst-case gain=4.474571 at step 3

What to inspect while the iteration runs

A robust experiment records more than the final vector. Track the residual ‖Aq − λ̂q‖, which measures how close the candidate is to an eigenvector. Track the angle or correlation between successive directions, while remembering that slow change can also mean stagnation. Repeat from several starts. Compare the leading pair of eigenvalues when possible. Perturb the matrix within plausible measurement error and see whether the direction remains stable.

Stopping because two consecutive vectors look similar is unsafe when the leading modes are close. They may rotate slowly within an almost invariant subspace. Stopping because the residual is small proves approximate invariance, but not that the eigenpair is the one the application needs. Convergence evidence and relevance evidence are separate ledgers.

Part III

When a workflow settles into a mode

Let each column describe where a case in one state goes at the next transition. The columns sum to one, so the matrix preserves total probability:

The dominant eigenvalue is one. Its eigenvector, normalised to sum to one, is approximately (0.439, 0.364, 0.197). Repeated transitions drive several different starting mixtures towards those proportions. The eigenvector is stationary because multiplying it by P returns the same mixture.

This is a useful persistent direction because the scale is fixed by probability conservation. In many other eigenvector problems, scaling is arbitrary: v and 7v represent the same direction. Here the “sum to one” constraint turns the direction into a unique distribution.

Yet it is easy to overread. The stationary proportions do not tell us how many cases exist, how long a transition takes, whether the queue is stable under new arrivals, whether specialists have enough capacity, or whether the transition matrix remains valid after a policy change. The eigenvector answers the question encoded by the operator, not the question the operator omitted.

Same stationary mixture, different operational memory

Create a stickier variant by mixing seventy per cent identity with thirty per cent of the original transition matrix:

Starting with every case waiting for evidence, the baseline comes within five percentage points of the stationary mixture after three transitions. The sticky version needs fourteen. Nothing about the destination reveals that difference. The subdominant modes carry the workflow’s memory of where it began.

In a production-shaped analysis, this pattern can become an early diagnostic. A large stationary share in “waiting” may indicate a structural accumulation only if the model also represents entry and exit correctly. A subdominant eigenvector can reveal which contrast decays slowly, such as “waiting versus active review”. The associated eigenvalue estimates persistence per transition. Operations teams can then test the specific transition probabilities that sustain that mode instead of treating the whole matrix as an opaque score.

The approach is also sensitive to perturbation. Ng, Zheng and Jordan showed in link-analysis settings that a small change to a matrix can cause a large change in its principal eigenvector when the relevant eigengap is small. [6] That result is a direct warning for governed use: if estimates come from sparse observations, or if policy changes alter transitions, report eigenvector stability under plausible perturbations. A precise-looking stationary vector can be an unstable consequence of uncertain inputs.

Why eigenvectors recur in machine intelligence

The same mechanism appears under different names. Principal-component methods seek directions of persistent variance, commonly through eigenvectors of a covariance matrix or singular vectors of a data matrix. Oja’s early neural rule showed how a simplified neuron can adapt towards a principal component. [9] Link-analysis algorithms turn repeated movement or endorsement into a stationary score. Recurrent systems repeatedly multiply or locally linearise state and gradient transformations, making spectral behaviour relevant to vanishing and exploding signals.

These are not interchangeable applications. A covariance matrix is symmetric positive semidefinite; a Markov operator preserves probability; a recurrent Jacobian can be non-normal and change at every step. Their eigenvectors live in different semantic spaces. The transferable idea is selective persistence under repeated transformation. The assumptions, norm and interpretation must be rebuilt for each operator.

Part IV

The spectrum is the destination, not the whole journey

The clean story assumes a fixed, diagonalizable operator with a uniquely dominant eigenvalue modulus, a starting state that excites the dominant mode, and a decision that cares about sufficiently long horizons. Each clause can fail. The right response is not to abandon eigenstructure. It is to stop asking it to certify more than it contains.

Ties, signs, phases and missing bases

If two eigenvalues share the largest modulus, one direction need not win. The state may remain in a multi-dimensional invariant subspace. With eigenvalues 0.9 and −0.9, the relative sign alternates. With a complex-conjugate pair, a real state can rotate while its envelope grows or decays. NumPy’s documentation explicitly notes that real matrices may return real eigenvalues or complex-conjugate pairs. [3]

If the matrix is defective, there may not be enough independent eigenvectors to form a basis. Jordan chains introduce polynomial factors such as kλᵏ. The eigenvalue still controls the exponential envelope, but the polynomial can dominate finite horizons. Near-defective matrices also make eigenvectors highly sensitive. A report that prints eigenvalues without conditioning or invariant-subspace checks can be numerically precise and scientifically brittle.

Trefethen et al. documented the physical importance of this distinction in hydrodynamic stability. They showed that small perturbations can be linearly amplified even when every eigenmode decays, and used pseudospectra to reconcile that behaviour with the stable spectrum. [8] The production lesson is broader than fluid mechanics: when consequences occur before asymptopia, measure finite-horizon gain directly.

For a discrete linear operator, inspect max₀≤k≤H ‖Aᵏ‖ over the decision horizon H. Singular values identify the input direction with greatest amplification at a chosen step. Pseudospectra assess sensitivity and possible amplification around non-normal operators. These quantities answer different questions from the spectral radius. They should not be collapsed into a single “stability score”.

From a witness trajectory to a release test

Figure 7 follows one chosen start. That is enough to prove that amplification can occur, but it does not bound how large the excursion could be. A favourable start can miss the dangerous combination of coordinates entirely. At each step, the largest singular value of A^k gives the greatest Euclidean amplification over every unit-length starting direction. Taking the maximum across the relevant steps converts the counterexample into a finite-horizon test:

Apply this test to the two matrices for eight steps. The diagonal matrix has G₈ = 1.000 at step zero, and every transformed state contracts. The non-normal matrix reaches G₈ = 4.475 at step three. The plotted start (0, 1) reaches 4.445, close to the maximum, while the start (1, 0) shrinks monotonically. A test using only that favourable start would report no hazard. The matrices still share eigenvalues 0.8 and 0.6; only the probe changes.

This calculation is meaningful only after the state norm is meaningful. If one coordinate is seconds and another is thousands of cases, raw Euclidean length makes the larger numerical scale dominate. A production analysis should define a weighting matrix W from tolerances or consequence, then inspect ‖WA^kW⁻¹‖₂. The weighting is part of the claim: changing it changes which perturbation counts as large. Report it with the operator rather than hiding it in preprocessing.

A production-shaped matched test

Run the analysis as a matched intervention, not as an eigenvalue printout. Hold the operator, horizon, norm and threshold fixed while changing only the initial direction. Then hold the eigenvalues fixed while changing the eigenvector geometry, as above. Finally perturb uncertain matrix entries within their justified estimation ranges. The following release record keeps those comparisons tied to an action.

The observation that would weaken the persistent-direction account is now explicit. If the fitted operator cannot predict held-out one-step changes, if its leading subspace changes under plausible re-estimation, or if observed excursions fall outside its finite-horizon envelope, the matrix is not carrying the relevant mechanism. The response is to revise the operator or adopt a switching, forced or nonlinear model. More digits on the same eigenvalues do not repair a failed system boundary.

A release decision needs both an asymptotic explanation and a finite-horizon receipt. The receipt should name the operator version, horizon, norm, threshold, peak gain, perturbation range, observed readback and resulting action. This makes a future disagreement diagnosable: reviewers can see whether the operator changed, the dangerous direction was missed, the threshold moved or the linear model stopped matching the system.

Recurrent systems expose the model boundary

In a linear recurrent system with a fixed recurrent matrix, powers of that matrix govern how earlier state components survive. Pascanu, Mikolov and Bengio use power-iteration reasoning to analyse when long-term gradient components vanish or explode in the linearised case. [7] This supports the mechanism, but it also marks its boundary. A practical recurrent network produces products of Jacobians that depend on time, input and state:

Applying the eigenvalues of one average Jacobian as a universal explanation can therefore fail. The relevant objects may be singular values of products, Lyapunov exponents, invariant bundles or empirical perturbation growth. Local eigenanalysis can still help near a fixed point, where a nonlinear map is well approximated by one Jacobian. It ceases to be decisive when the trajectory moves through substantially different regimes.

The strongest objection: the operator may be the invention

An eigenvector is exact relative to its matrix. The matrix may be an estimate, a convenient linearisation, a choice of variables, a discretisation, a policy snapshot or an arbitrary similarity representation. Eigenvalues survive a change of basis, but component meanings, Euclidean angles and norms do not automatically survive changes in units or semantics. A dramatic mode can be an artefact of encoding.

That objection is not external philosophy. It changes the experiment. Ask whether one application of the matrix corresponds to an observable transition. Check whether linear superposition is plausible in the regime. Scale variables according to physical or decision meaning. Compare predicted and observed trajectories. Re-estimate across periods. Perturb uncertain entries. If the mode is used to authorise an intervention, identify which measurable mechanism creates the matrix entries.

Conclusion

Design for the mode, test for the transient

A repeated linear transformation contains a selection mechanism. Decompose a state into invariant directions, and iteration raises each eigenvalue to a power. The directions do not compete through deliberation or optimisation. They persist unequally. The dominant mode is what remains after the operator has repeatedly expressed that inequality.

That mechanism earns several useful decisions. Use a dominant eigenvector to identify a stationary composition or asymptotic shape. Use the leading eigenvalue to describe its scaling. Use the runner-up and the eigengap to estimate how long the system remembers its start. Use power iteration when repeated matrix-vector products match the computational boundary, and use block or Krylov methods when one direction is not enough.

The argument also changes the release gate. A spectral radius below one is an asymptotic statement, not a finite-horizon safety case. Before relying on it, establish that the mode is excited, the separation is adequate, the eigenvectors or invariant subspace are well-conditioned, and transient amplification stays below the consequence threshold. Then verify that one fixed linear operator represents the horizon at all.

A defensible review should therefore trace the whole claim: from the definition of state, to the construction of the operator, to the estimated mode, to the horizon over which that mode changes a decision. An eigenvalue table alone stops the evidence chain too early. The claimed direction should predict an observable trajectory, remain recognisable under plausible perturbations and survive comparison with a serious finite-horizon baseline.

Persistent directions can also point towards interventions. To shorten unwanted memory, identify and change the transitions that sustain a slow subdominant mode. To preserve a useful signal, avoid transformations that repeatedly suppress its direction. To compute a clustered leading space, retain several vectors rather than forcing a single winner. These are intervention hypotheses until matrix entries are tied to controllable mechanisms and the altered system is tested again.

Negative results remain informative. A tied leading pair says that the durable object is a subspace. A missing projection says that the probe cannot see the mode. A large transient says that containment must be designed for the journey rather than the destination. A drifting operator says that the fixed-matrix question was posed at the wrong level. Each result narrows the next experiment instead of being treated as a failed eigendecomposition.

The practical decision is therefore precise: approve eigenstructure as the explanation only when mode, gap, excitation, transient and model boundary all survive their own tests. Otherwise retain the intuition of persistent directions, but move the evidence to the richer object the system actually requires.