How can I use Mathematica to get $f(y)$, the PDF of a random variable given its Stieltjes transform $g(s)$?


$$g(s)=E\left[1/(s-X)\right]$$

I need to visualize probability density of a random variable with the following Stieltjes transform:

```
stieltjes=(Sqrt[2] (\[Pi] - 2 ArcTan[Sqrt[s]/Sqrt[2]]))/(2 Sqrt[s] + Sqrt[2] s ArcTan[Sqrt[2]/Sqrt[s]])
```

-------
Trying on some simpler distributions first:

Using Beta$(1/2,1)$ we have the following pdf $f(y)$ and $g(s)$

```
pdf = 1/2 y^(-1/2); (* BetaDistribution[1/2,1] *)
stieltjes = Integrate[pdf 1/(s - y), {y, 0, 1}]
```
$$f(y)=\frac{\frac{1}{\sqrt{y}}}{2} \leftrightarrow g(s)=\frac{\coth ^{-1}\left(\sqrt{s}\right)}{\sqrt{s}}$$

Also, using $1/x^2$ where $x$ is Zipf distributed as in the earlier [question](https://mathematica.stackexchange.com/a/286827/217) we have the following
```
distY = TransformedDistribution[1/z^2, 
   z \[Distributed] ZipfDistribution[1]];
Expectation[1/(s - i), i \[Distributed] distY]
```

$$f(y)=\left(1, \frac{6}{\pi ^2}\right),\left(\frac{1}{4}, \frac{3}{2 \pi ^2}\right), \left(\frac{1}{9} , \frac{2}{3 \pi ^2}\right),\ldots  \leftrightarrow g(s)=-\frac{3 \left(\pi  \cot \left(\frac{\pi }{\sqrt{s}}\right)-\sqrt{s}\right)}{\pi ^2 \sqrt{s}}$$

You can get Stieltjes Transforms by chaining two Laplace transforms together as seen below. Wondering if there's a way to use this to get symbolic or numeric inversion working for examples above.

```
dist = BetaDistribution[1/2, 1];
mgf = MomentGeneratingFunction[BetaDistribution[1/2, 1], t];
doubleLaplace = LaplaceTransform[mgf, t, s] // FullSimplify;
stieltjes = Expectation[1/(s - y), y \[Distributed] dist];
doubleLaplace == stieltjes (*True if s>1*)
```


(note you often need Feynman trick to do inverse Laplace transforms in Mathematica, documented [here](https://www.wolframcloud.com/obj/yaroslavvb/nn-linear/laplace-transform-tricks.nb))