Spectral density of $D + XX^T$

Question

Let $D$ be a fixed diagonal matrix with real entries, and $X$ a random $m\times n$ matrix. More precisely, the entries $X_{ij}$ are real independent and identically distributed. It can be shown that the eigenvalues of $XX^T$ have a density following the Marcenko-Pastur law with probability 1, in the limit that $n,m\rightarrow\infty$ with a fixed ratio $m/n$, and this result depends only on the variance of the entries $X_{ij}$.

Can we say something about the spectral density of $D + XX^T$, also in the limit of large $m,n$? Is an analytical expression (as in the Marcenko-Pastur law) known in this case? At least something numerically tractable?

Answer 1

The solution was obtained already by Marchkenko and Pastur, in terms of the Cauchy transform $g(z)$ of the spectral density of $D+XX^{\rm T}$, see for example equations 2 + 3 of Spectrum of deformed random matrices and free probability:

$$g(z)=\int\frac{1}{z-t-(m/n)(z-t)g(z)-1+m/n}\rho(t)\,dt$$

where $\rho(t)$ is the spectral density of $D$. This holds for any random Hermitian perturbation $D$, irrespective of whether it is diagonal or not.

---

Let me try an alternative approach; The Cauchy transform of the Marchenko-Pastur distribution, supported on $[(1-\sqrt{\lambda})^2,(1+\sqrt{\lambda})^2]$, is given by $$g_{\rm MP}(z)=\frac{-\sqrt{(-\lambda+z-1)^2-4 \lambda}-\lambda+z+1}{2 z}.$$ The perturbation $D$, with $\rho(t)=w\delta(t-1)+(1-w)\delta(t)$, has Cauchy transform $$g_D(z)=\frac{w}{z-1}+\frac{1-w}{z}.$$ The Cauchy transform $g_{D+{\rm MP}}(z)$ of $D+XX^{\rm T}$ is given by the composition of free probability, $$R_{D+{\rm MP}}=R_{\rm MP}+R_D,\;\;R(z)=z-1/g(z).$$ This gives the Cauchy transform $$g_{D+{\rm MP}}(z)=\frac{1}{z}\left[ -\frac{2}{\sqrt{(\lambda-z+1)^2-4 \lambda}+\lambda-z-1}+\frac{z-1}{w+z-1}-1\right]^{-1}.$$

Now the goal is to integrate $\langle\log |t|\rangle_t=\int \rho_{D+{\rm MP}}(t)\log t\,dt$, which might be obtainable directly from the Cauchy transform via $$\frac{d}{dz}\langle\log |z-t|\rangle_t = \int \frac{\rho_{D+{\rm MP}}(t)}{z-t} \, dt = g_{D+{\rm MP}}(z).$$ The integration constant can be obtained from the large-$z$ limit, which should be just $\log z$.

Answer 2

Let me add some remarks in addition to the answer of Carlo, in particular, concerning the numerical side. In free probability terms, the solution is given by the free convolution of the Marchenko-Pastur distribution with the limiting distribution of $D$. There are various ways of describing the Cauchy transform $g_{D+MP}$ of this free convolution. The original approach in free probability was to use the $R$-transform and the fact that this is additive for free variables. This is applied (for the special case where $D$ has only atoms at 0 and 1) in the alternative approach of Carlo. What I want to point out is that this usually not the best way for numerical approximations, in cases where analytical solutions are not obtainable. Then the subordination formulation of free convolution $$g_{D+MP}(z)=g_D(\omega(z))$$ is more useful (in this case this is related to the integral equation in Carlo's answer). Here the Cauchy transform of the free convolution is subordinated to the Cauchy transform of one of the summands, say $D$, and the main point is now to calculate this subordination function $\omega$ (which depends on both input distributions). The nice point about this approach is that $\omega(z)$ can be calculated as the fixed point of an equation and that actually this solution is given as the limit of iterates of the fixed point map; this is very useful for numerical calculations. For more concrete information and examples, see, e.g., the paper of Capitaine and Donati-Martin, which was mentioned in Carlo's answer, or Chapter 3 of my book Free probability and random matrices or also a video-recorded lecture of mine on this topic.