# Efficient computation of leading eigenvector of a matrix product of the form $ADA^T$, where $D$ is diagonal

Let $$A=[A_1|\ldots|A_m] \in \mathbb R^{n \times m}$$ with $$n \gg m \gg 1$$ and $$D=\text{diag}(d_1,\ldots,d_m)$$ where $$d_1,\ldots,d_m > 0$$, and consider the $$n\times n$$ positive-definite matrix $$X=\sum_{i=1}^m d_i A_iA_i^T=ADA^T$$.

# Question

What is an efficient way to compute the leading eigenvector of $$X$$ without forming the product $$ADA^T$$ ?

Assuming that multiplication by a diagonal matrix is computationally easy, we can write it as: $$ADA^T = (A\sqrt{D})(\sqrt{D}^TA^T) = BB^T$$ where $$B=A\sqrt{D}$$. Motivated by the fact that the left-singular vectors of $$B$$ are a set of orthonormal eigenvectors of $$BB^T$$, we can further proceed: $$B = U\Sigma V^T$$ The first columns $$m$$ columns of $$U$$ (corresponding to $$m$$ highest singular values) should correspond to the first $$m$$ eigenvectors of $$ADA^T$$. Since you only want the leading eigenvector, this will suffice.

Here is a short MATLAB code to realize the idea:

m = 5; d = 3;
A = randn(m,d); % A is a random matrix - note that it can have negatives
D = diag(rand(d,1)); % just a random diagonal that is positive
[V, ~] = eigs(A*D*A', d); % what we actually like to have
[U, ~, ~] = svd(A*sqrt(D),'eco'); % the variant avoiding the product


A sample run of the code above produces the following $$U$$ and $$V$$:

$$U = \begin{bmatrix} -0.4276 & 0.0746 & -0.6573 \\ -0.3001 & 0.2736 & 0.6797 \\ -0.5469 & 0.5083 & -0.1503 \\ -0.4068 & 0.0579 & 0.2730 \\ -0.5124 & -0.8111 & 0.0941 \end{bmatrix}\quad V = \begin{bmatrix} -0.4276 & -0.0746 & 0.6573 \\ -0.3001 & -0.2736 & -0.6797 \\ -0.5469 & -0.5083 & 0.1503 \\ -0.4068 & -0.0579 & -0.2730 \\ -0.5124 & 0.8111 & -0.0941 \end{bmatrix}$$

Note the sign ambiguity. Signs of the vectors are arbitrary for the eigen-decomposition anyway. If you like to get the single leading eigenvector, another option is to directly replace the last line with:

[U, ~, ~] = svds(A*sqrt(D),1);

• Great, thanks. I also got to a similar idea based on power iterations. This will be super efficient since applying the matrix $X=ADA^T$ to a vector $v$ gives $Xv=A(d\circ (Av))$ which can be done efficiently (more precisely, in $\mathcal O(mn)$ flops) without ever forming the matrix product $ADA^T$ – dohmatob Oct 18 '19 at 12:59