Update I want to add a "philosophical" comment. The question you asked is a special case of the following more general question.
Suppose that $f:\mathbb{R}\to\mathbb{R}$ is a continuous function. For simplicity, let us assume it is also bounded. We can define the bounded symmetric operator $f(\sqrt{\Delta})$ where $\Delta$ is the Laplacian on an $m$-dimensional manifold $M$. Investigate the behavior of $f(\varepsilon\Delta)$ as $\varepsilon \to 0$.Your case corresponds to $f(x)=(1+x^2)^{-1}$. The heat equation problems correspond to $f(x)=e^{-x^2}$. Suppose that $f$ is a symbol of order $k$, where $k$ could be $-\infty$. For example $(1+x^2)^{-k}$ is a symbol of order $-2k$, while $e^{-x^2}$ is a symbol of order $-\infty$.
In any case, when $f$ is a symbol, then $f(\Delta)$ is a pseudodifferential operator, and as such it has a Schwartz kernel which is a distribution on $M\times M$. $\newcommand{\ve}{\varepsilon}$ Your question is about the behavior as $\ve\to 0$ of the Schwartz kernel of $f(\ve \Delta)$ along normal directions to the diagonal of $M\times M$.
If $f$ is rapidly decaying at $\infty$, say $f(x) < (1+x^2)^{-m}$, $m=\dim M$, then the Schwartz kernel of $f(\Delta)$ is given by a continuous function and one can be quite precise about the behavior of the kernel of $f(\ve \Delta)$. In fact, the faster the decay of $f$ at $\infty$, the more accurate one can be about the behavior of the Schwartz kernel of $f(\ve \Delta)$. The radial symmetry you are talking about is then a simple consequence if $f$ decays faster than $|x|^{-N}$, $N$ sufficiently large. (I believe that $N>2m$ ought to do it but I don't want to be too firm.) If $f$ has exponential decay at $\infty$ one can be remarkable accurate and recover the radial symmetry you are mentioning. Your question involves the symbol $(1+x^2)^{-1}$ that isn't decaying fast enough at $\infty$. Translation: your problem requires a bit of care.

