I'm reading a paper which said that the Green function for $\left(-\Delta_g\right)^m$ on $2m$-dimensional closed manifold is of the form $$\tag{1} G_y(x)=\frac{2}{\Lambda_1} \log \frac{1}{d_g(x, y)}+\gamma(x, y), $$ where $\gamma$ is smooth on $M \times M$ and $\Lambda_1=(2 m-1) ! \operatorname{vol}\left(S^{2 m}\right)$. I searched this result in many textbooks, they all only showed the existence of Green function of Laplacian and the estimate of it on closed manifold of dimension $n$, for example: There exists a constant $k$ such that: $|\mathbf{G}(P, Q)|<k(1+|\log r|)$ for $n=2$ and
$$ \begin{aligned} & |\mathbf{G}(P, Q)|<k r^{2-n} \text { for } n>2, \text{with} \ r=d(P, Q) . \end{aligned} $$ I want to ask where fomula(1) comes from and what is this $\gamma$?
I found a reference Here, which calculate the Green function on compact locally harmonic Blaschke manifolds, the proof may take a long time to read, could you recommend me some references on calculating the Green function on sphere? So that I could understand the difference of calculating the Green function on closed manifold and on domain with boundary by specific calculation.