The way principal component analysis is usually computed is a one-shot recipe: form the covariance matrix \(\mathbf{C}\), diagonalize it, and read off the eigenvectors. That recipe is exact, but it assumes \(\mathbf{C}\) can actually be formed and diagonalized, which costs \(O(d^3)\) time and \(O(d^2)\) space for \(d\) variables. For a modest image of a few hundred pixels this is routine; for the internal activations of a large neural network, where \(d\) can run into the tens of millions, it is not even close to feasible — the covariance matrix alone would need more entries than there are atoms in a grain of sand. This raises a genuinely different question from the one the main chapter answers: can the same principal components be found by a process that never materializes the full covariance matrix at all, and that can be spread across many machines working at once?
One elegant answer, due to Gemp, McWilliams, Vernade, and Graepel, recasts principal component analysis as a competitive game: each principal component is controlled by its own player, every player tries to maximize its own private payoff, and the unique stable outcome of the game — its Nash equilibrium — turns out to be exactly the ordered set of principal components. The resulting method, EigenGame, needs only matrix-vector products, updates every component with the same simple rule, and lets each player learn from a different machine while only ever broadcasting a single vector. This chapter derives the game from the two properties of the eigenvalue equation the main chapter already established, states the equilibrium result, and works through the resulting algorithm on a small worked example (Gemp et al., 2021, EigenGame: PCA as a Nash Equilibrium, published at ICLR 2021).
Why Maximizing Variance Alone Is Not Enough
Write any \(d\) candidate directions as the columns of an orthonormal matrix \(\hat{\mathbf{V}}\) (\(\hat{\mathbf{V}}^T\hat{\mathbf{V}}=\mathbf{I}\)), and collect how they interact with the covariance matrix \(\mathbf{C}\) into a single small matrix,
\[\mathbf{R}:=\hat{\mathbf{V}}^T\mathbf{C}\hat{\mathbf{V}}.\]
The diagonal entries \(R_{ii}=\hat{\mathbf{v}}_i^T\mathbf{C}\hat{\mathbf{v}}_i\) are exactly the Rayleigh quotients from the main chapter, so maximizing their sum looks like a natural way to recover the principal directions. But if \(\hat{\mathbf{V}}\) is a full \(d\times d\) orthonormal matrix, this sum is constant:
\[\operatorname{tr}(\mathbf{R})=\operatorname{tr}(\hat{\mathbf{V}}^T\mathbf{C}\hat{\mathbf{V}}) =\operatorname{tr}(\mathbf{C}\hat{\mathbf{V}}\hat{\mathbf{V}}^T)=\operatorname{tr}(\mathbf{C}),\]
using the cyclic property of the trace and \(\hat{\mathbf{V}}\hat{\mathbf{V}}^T=\mathbf{I}\) for a square orthonormal matrix. Every orthonormal basis of the whole space gives the same total variance, whether or not its vectors have anything to do with the eigenvectors of \(\mathbf{C}\). Restricting to only the first \(k<d\) columns does not fix this: any orthonormal basis of the same \(k\)-dimensional subspace achieves the same trace, even one that is rotated arbitrarily within that subspace and is nowhere close to the true components. A quick numerical check makes this concrete. For a \(5\times5\) covariance matrix with top two eigenvalues \(15.509\) and \(13.347\), the true top-2 eigenvectors give \(\operatorname{tr}(\mathbf{R})=28.856\), and so does an arbitrary rotation of those same two eigenvectors within their span — the objective genuinely cannot tell them apart. Maximizing variance alone, in other words, only ever recovers the right subspace, never the individual directions within it.
A Complementary Objective: Making the Cross-Terms Vanish
The off-diagonal entries of \(\mathbf{R}\), \(R_{ij}=\hat{\mathbf{v}}_i^T\mathbf{C}\hat{\mathbf{v}}_j\), measure how much the candidate directions \(i\) and \(j\) agree under the metric induced by \(\mathbf{C}\). If \(\hat{\mathbf{v}}_i\) and \(\hat{\mathbf{v}}_j\) are genuine, distinct eigenvectors \(\mathbf{v}_i,\mathbf{v}_j\) of \(\mathbf{C}\) with \(i\neq j\), then
\[R_{ij}=\mathbf{v}_i^T\mathbf{C}\mathbf{v}_j=\lambda_j\mathbf{v}_i^T\mathbf{v}_j=0\]
by orthogonality of eigenvectors belonging to different eigenvalues, so \(\mathbf{R}\) is exactly diagonal only when the columns of \(\hat{\mathbf{V}}\) are eigenvectors. Continuing the numerical example from above, the same arbitrarily rotated top-2 subspace that tied the trace objective gives \(\mathbf{R}=\left(\begin{smallmatrix} 13.426 & 0.406\\0.406 & 15.430\end{smallmatrix}\right)\) — visibly off-diagonal — whereas the true eigenvectors give the diagonal \(\mathbf{R}=\left(\begin{smallmatrix}15.509 & 0\\0 & 13.347\end{smallmatrix} \right)\), recovering both eigenvalues on the diagonal. Minimizing the squared off-diagonal terms of \(\mathbf{R}\) is therefore a natural second objective, complementary to maximizing the trace: it pins down the specific directions inside a subspace that maximizing variance alone leaves undetermined.
Turning Two Objectives into a Hierarchy of Player Utilities
Simply adding the two objectives together, maximize the trace while minimizing every squared off-diagonal term, treats every pair of directions symmetrically. That is undesirable: the first principal component should be free to seek out the single direction of maximum variance without being penalized for what the fifth, tenth, or hundredth component happens to be doing. EigenGame breaks the symmetry by imposing an order on the components and letting each one only be judged against the components ranked above it, its parents. Player \(i\) controls a unit vector \(\hat{\mathbf{v}}_i\), observes the current guesses \(\hat{\mathbf{v}}_1,\ldots,\hat{\mathbf{v}}_{i-1}\) of its parents, and tries to maximize its own utility
\[\boxed{u_i(\hat{\mathbf{v}}_i\mid \hat{\mathbf{v}}_{j<i})=\hat{\mathbf{v}}_i^T\mathbf{C}\hat{\mathbf{v}}_i -\sum_{j<i}\frac{(\hat{\mathbf{v}}_i^T\mathbf{C}\hat{\mathbf{v}}_j)^2}{\hat{\mathbf{v}}_j^T\mathbf{C} \hat{\mathbf{v}}_j}}\qquad\text{subject to }\|\hat{\mathbf{v}}_i\|=1.\]
The first term is the familiar Rayleigh quotient, rewarding directions of high variance; the second term penalizes alignment with each higher-ranked parent, scaled by that parent's own variance so the two terms stay comparable in size. Player 1 has no parents, so its utility is simply the Rayleigh quotient, maximized exactly by the leading eigenvector as shown in the main chapter. At the point where every parent already sits exactly on its true eigenvector, \(\hat{\mathbf{v}}_j=\mathbf{v}_j\), the penalty term for a candidate \(\hat{\mathbf{v}}_i\) becomes \((\hat{\mathbf{v}}_i^T\mathbf{C}\mathbf{v}_j)^2/\lambda_j=\lambda_j(\hat{\mathbf{v}}_i^T\mathbf{v}_j)^2\), which vanishes only when \(\hat{\mathbf{v}}_i\) is orthogonal to \(\mathbf{v}_j\) — exactly the constraint needed to keep the players from collapsing onto each other.
Nash Equilibrium: Every Player Content With Its Own Choice
A set of strategies, one per player, is a Nash equilibrium of a game if no single player can improve its own utility by unilaterally switching to a different strategy while every other player keeps theirs fixed. Formally, \(\hat{\mathbf{v}}_1,\ldots,\hat{\mathbf{v}}_k\) is a (strict) Nash equilibrium of the game above if, for every player \(i\),
\[u_i(\hat{\mathbf{v}}_i\mid \hat{\mathbf{v}}_{j<i})>u_i(\mathbf{z}\mid \hat{\mathbf{v}}_{j<i})\qquad \text{for every other unit vector }\mathbf{z}\neq\hat{\mathbf{v}}_i.\]
Gemp et al. (2021) prove that, whenever the top \(k\) eigenvalues of \(\mathbf{C}\) are positive and distinct, this game has exactly one such equilibrium (up to the usual sign ambiguity already discussed in the main chapter), and it is precisely the ordered sequence of the true top-\(k\) eigenvectors of \(\mathbf{C}\). No player can do better by deviating once every player already holds its true principal component, and this is the only configuration with that property. This reframes an entirely classical linear-algebra problem in the language of game theory: principal component analysis is the unique equilibrium of a particular \(k\)-player, non-cooperative game played on the unit sphere.
Utility Gradients: Oja's Rule Meets Gram–Schmidt
Differentiating the utility with respect to \(\hat{\mathbf{v}}_i\) (holding the parents fixed) gives
\[\nabla_{\hat{\mathbf{v}}_i}u_i=2\left(\mathbf{C}\hat{\mathbf{v}}_i-\sum_{j<i} \frac{\hat{\mathbf{v}}_i^T\mathbf{C}\hat{\mathbf{v}}_j}{\hat{\mathbf{v}}_j^T\mathbf{C}\hat{\mathbf{v}}_j}\, \mathbf{C}\hat{\mathbf{v}}_j\right),\]
which was checked against a numerical finite-difference gradient on a random \(5\)-dimensional example and matched to within \(10^{-9}\). The two pieces of this gradient have familiar names. The term \(\mathbf{C}\hat{\mathbf{v}}_i\) alone is exactly the direction used by Oja's rule, a classical update for tracking the single dominant eigenvector by repeated multiplication with \(\mathbf{C}\) (Oja, 1982). The subtracted sum is a generalized Gram–Schmidt step: it removes from \(\mathbf{C}\hat{\mathbf{v}}_i\) exactly the component that would otherwise pull \(\hat{\mathbf{v}}_i\) back towards each of its already-settled parents. Every player runs the same rule; only the identity of its parents changes from player to player.
A Sequential Algorithm on the Unit Sphere
Because \(\|\hat{\mathbf{v}}_i\|=1\) is a hard constraint, the utility is not maximized over all of \(\mathbb{R}^d\) but only over the unit sphere, a curved surface rather than a flat space. Naively stepping in the direction of \(\nabla_{\hat{\mathbf{v}}_i}u_i\) would walk off the sphere, so the update first discards the component of the gradient that points radially away from the sphere,
\[\nabla^R_{\hat{\mathbf{v}}_i}u_i=\nabla_{\hat{\mathbf{v}}_i}u_i -\left(\hat{\mathbf{v}}_i^T\nabla_{\hat{\mathbf{v}}_i}u_i\right)\hat{\mathbf{v}}_i,\]
leaving only the component tangent to the sphere at \(\hat{\mathbf{v}}_i\), and then takes a small step and renormalizes back onto the sphere,
\[\boxed{\hat{\mathbf{v}}_i\leftarrow\frac{\hat{\mathbf{v}}_i+\alpha\,\nabla^R_{\hat{\mathbf{v}}_i}u_i} {\left\|\hat{\mathbf{v}}_i+\alpha\,\nabla^R_{\hat{\mathbf{v}}_i}u_i\right\|}}\]
for a small step size \(\alpha>0\). Repeating this update learns the parents first, then their children, so that by the time a player differentiates its own utility, the parents it depends on are already fixed at their converged values.
The covariance matrix behind this picture is \(\mathbf{C}=\left(\begin{smallmatrix}3.0 & 1.2\\1.2 & 1.5 \end{smallmatrix}\right)\), with true eigenvalues \(3.665\) and \(0.835\). Starting player 1 at an arbitrary initial guess \(\hat{\mathbf{v}}_1(0)\) and running the update above converges it, within a handful of steps, onto the true leading eigenvector \(\mathbf{v}_1^*\); starting player 2 from another arbitrary guess and running the same update, but now with player 1's converged vector fixed as its parent, sweeps it along the circle until it settles orthogonally onto \(\mathbf{v}_2^*\), exactly as the theorem predicts.
Decentralizing the Computation
The hierarchy of parents forms a directed chain: player \(i\) needs to know the current vectors of players \(1,\ldots,i-1\), but nothing about players \(i+1,\ldots,k\). This one-directional dependency is what makes the method naturally parallelizable: each player can run its own copy of the update rule on its own machine, broadcasting only its current unit vector to whichever players still need it, and never needs to see the full covariance matrix — only matrix-vector products with the data are required. In practice, updating every player simultaneously rather than strictly in sequence still converges reliably, because as parents approach their equilibrium they become nearly stationary targets for their children, even though the formal convergence guarantee assumes the strict sequential order. Gemp et al. (2021) demonstrate the resulting scalability by computing the top \(32\) principal components of the activations of a trained residual network on a large image dataset, with a dimensionality of roughly \(d\approx2\times10^7\); forming and diagonalizing the corresponding covariance matrix directly would require storing and decomposing a matrix with about \(4\times10^{14}\) entries, far beyond what any direct eigendecomposition could handle, while the decentralized game scales by simply adding one more player per additional component.
Relation to Classical Iterative Methods
None of the individual ingredients here are new on their own. Oja's rule has long been used to track a single dominant eigenvector without ever forming the covariance matrix, and extending this to several components classically means alternating an update step with an explicit re-orthonormalization of the whole set of vectors, typically via a QR factorization — an inherently centralized operation that has to see every vector at once before it can proceed. What the game formulation contributes is folding that orthonormalization directly into each player's own utility function as a smooth penalty term, so that no separate, centralized orthonormalization pass is ever required: every player's local gradient already points in a direction that keeps it orthogonal to its parents, without anyone having to coordinate the whole set of vectors at once.
References
- Gemp21 Gemp, I., McWilliams, B., Vernade, C., & Graepel, T. (2021). EigenGame: PCA as a Nash Equilibrium. International Conference on Learning Representations (ICLR 2021). arXiv:2010.00554.
- Oja82 Oja, E. (1982). A Simplified Neuron Model as a Principal Component Analyzer. Journal of Mathematical Biology, 15(3), 267–273.