[proofplan] We diagonalize the real symmetric matrix $A_d$ by an orthogonal change of basis and shift the diagonalization by the non-real scalar $z$. Since every eigenvalue of $A_d$ is real, none of the diagonal entries $\lambda_j(A_d)-z$ vanishes, so the shifted matrix is invertible and its inverse has diagonal entries $(\lambda_j(A_d)-z)^{-1}$ in the eigenbasis. Taking the trace is invariant under the orthogonal conjugation, and the resulting average is exactly the integral of $x \mapsto (x-z)^{-1}$ against the empirical spectral measure. [/proofplan]

[step:Diagonalize the symmetric matrix and record the empirical measure] By the finite-dimensional spectral theorem for real symmetric matrices (citing a result not yet in the wiki: Finite-Dimensional Spectral Theorem for Real Symmetric Matrices), there exist an orthogonal matrix $Q \in \mathbb{R}^{d \times d}$ and a diagonal matrix $\Lambda \in \mathbb{R}^{d \times d}$ such that \begin{align*} A_d = Q\Lambda Q^\top, \qquad \Lambda_{jj} = \lambda_j(A_d) \quad \text{for } j \in \{1,\dots,d\}. \end{align*} Here $Q^\top Q = QQ^\top = I_d$. By definition of the empirical spectral measure, \begin{align*} \int_{\mathbb{R}} \frac{1}{x-z}\,d\mu_{A_d}(x) = \frac{1}{d}\sum_{j=1}^{d}\frac{1}{\lambda_j(A_d)-z}, \end{align*} because $\mu_{A_d}$ is the finite measure assigning mass $1/d$ to each eigenvalue counted with multiplicity. [/step]

[step:Invert the shifted diagonal matrix] Since $z \in \mathbb{C}\setminus \mathbb{R}$ and each $\lambda_j(A_d) \in \mathbb{R}$, we have $\lambda_j(A_d)-z \neq 0$ for every $j \in \{1,\dots,d\}$. Hence the diagonal matrix $\Lambda-zI_d \in \mathbb{C}^{d \times d}$ is invertible, with inverse given by the diagonal matrix \begin{align*} (\Lambda-zI_d)^{-1}_{jj} = \frac{1}{\lambda_j(A_d)-z}, \qquad (\Lambda-zI_d)^{-1}_{ij} = 0 \quad \text{for } i \neq j. \end{align*} Using $A_d-zI_d = Q(\Lambda-zI_d)Q^\top$ and $Q^\top Q = QQ^\top = I_d$, we obtain \begin{align*} (A_d-zI_d)^{-1} = Q(\Lambda-zI_d)^{-1}Q^\top. \end{align*} [/step]

[step:Take the trace and identify the result with the Stieltjes transform]We compute the trace of the conjugated diagonal inverse. For any two matrices $B,C \in \mathbb{C}^{d \times d}$, the trace of $BC$ is \begin{align*} \operatorname{tr}(BC)=\sum_{i=1}^{d}\sum_{k=1}^{d} B_{ik}C_{ki}. \end{align*} Interchanging the two finite sums gives \begin{align*} \sum_{i=1}^{d}\sum_{k=1}^{d} B_{ik}C_{ki}=\sum_{k=1}^{d}\sum_{i=1}^{d} C_{ki}B_{ik}=\operatorname{tr}(CB). \end{align*} Thus trace is cyclic for products of two finite matrices. Applying this identity first to $Q$ and $(\Lambda-zI_d)^{-1}Q^\top$ gives \begin{align*} \operatorname{tr}\left((A_d-zI_d)^{-1}\right)=\operatorname{tr}\left((\Lambda-zI_d)^{-1}Q^\top Q\right). \end{align*} Since $Q^\top Q=I_d$, this becomes \begin{align*} \operatorname{tr}\left((A_d-zI_d)^{-1}\right)=\operatorname{tr}\left((\Lambda-zI_d)^{-1}\right). \end{align*} Because $(\Lambda-zI_d)^{-1}$ is diagonal, its trace is \begin{align*} \operatorname{tr}\left((\Lambda-zI_d)^{-1}\right)=\sum_{j=1}^{d}\frac{1}{\lambda_j(A_d)-z}. \end{align*} Dividing by $d$ and comparing with the formula for the integral against $\mu_{A_d}$ yields \begin{align*} \frac{1}{d}\operatorname{tr}\left((A_d-zI_d)^{-1}\right) = \frac{1}{d}\sum_{j=1}^{d}\frac{1}{\lambda_j(A_d)-z} = \int_{\mathbb{R}} \frac{1}{x-z}\,d\mu_{A_d}(x) = m_d(z). \end{align*} This proves the stated resolvent formula.[/step]

[guided]The point of the argument is to compute the resolvent in a basis where $A_d$ is diagonal. Since $A_d$ is real symmetric, the finite-dimensional spectral theorem for real symmetric matrices gives an orthonormal eigenbasis. Equivalently, there are an orthogonal matrix $Q \in \mathbb{R}^{d \times d}$ and a diagonal matrix $\Lambda \in \mathbb{R}^{d \times d}$ satisfying \begin{align*} A_d = Q\Lambda Q^\top, \qquad \Lambda_{jj} = \lambda_j(A_d) \quad \text{for } j \in \{1,\dots,d\}, \end{align*} with $Q^\top Q = QQ^\top = I_d$. Now shift by $z$. Since $z \notin \mathbb{R}$ while every eigenvalue $\lambda_j(A_d)$ is real, no denominator $\lambda_j(A_d)-z$ is zero. Therefore $\Lambda-zI_d$ is invertible, and its inverse is the diagonal matrix with diagonal entries \begin{align*} (\Lambda-zI_d)^{-1}_{jj} = \frac{1}{\lambda_j(A_d)-z}. \end{align*} The same orthogonal conjugation gives \begin{align*} A_d-zI_d = Q(\Lambda-zI_d)Q^\top, \end{align*} so multiplying directly by $Q(\Lambda-zI_d)^{-1}Q^\top$ on either side gives \begin{align*} (A_d-zI_d)^{-1} = Q(\Lambda-zI_d)^{-1}Q^\top. \end{align*} It remains to understand what happens to the trace. For finite matrices $B,C \in \mathbb{C}^{d \times d}$, the trace of $BC$ is \begin{align*} \operatorname{tr}(BC)=\sum_{i=1}^{d}\sum_{k=1}^{d} B_{ik}C_{ki}. \end{align*} Since both sums are finite, we may interchange their order and obtain \begin{align*} \sum_{i=1}^{d}\sum_{k=1}^{d} B_{ik}C_{ki}=\sum_{k=1}^{d}\sum_{i=1}^{d} C_{ki}B_{ik}=\operatorname{tr}(CB). \end{align*} This shows that we may cyclically move one factor across the trace. Applying this to the product $Q(\Lambda-zI_d)^{-1}Q^\top$ gives \begin{align*} \operatorname{tr}\left((A_d-zI_d)^{-1}\right)=\operatorname{tr}\left((\Lambda-zI_d)^{-1}Q^\top Q\right). \end{align*} Using $Q^\top Q=I_d$, we get \begin{align*} \operatorname{tr}\left((A_d-zI_d)^{-1}\right)=\operatorname{tr}\left((\Lambda-zI_d)^{-1}\right). \end{align*} Because $(\Lambda-zI_d)^{-1}$ is diagonal, its trace is the sum of its diagonal entries: \begin{align*} \operatorname{tr}\left((\Lambda-zI_d)^{-1}\right) = \sum_{j=1}^{d}\frac{1}{\lambda_j(A_d)-z}. \end{align*} Finally, the empirical spectral measure is \begin{align*} \mu_{A_d} = \frac{1}{d}\sum_{j=1}^{d}\delta_{\lambda_j(A_d)}. \end{align*} Define the function $g_z:\mathbb{R}\to \mathbb{C}$ by $g_z(x)=(x-z)^{-1}$ for every $x\in\mathbb{R}$. Therefore integrating $g_z$ against $\mu_{A_d}$ gives the finite average \begin{align*} m_d(z) = \int_{\mathbb{R}} g_z(x)\,d\mu_{A_d}(x) = \frac{1}{d}\sum_{j=1}^{d}\frac{1}{\lambda_j(A_d)-z}. \end{align*} Combining this with the trace computation yields \begin{align*} m_d(z) = \frac{1}{d}\operatorname{tr}\left((A_d-zI_d)^{-1}\right), \end{align*} which is the desired formula.[/guided]