I had success in computing the roots numerically for the Bessel polynomial $\theta_n(x) = x^ny_n(1/x)=\sum\limits_{k=0}^n\frac{(n+k)!}{(n-k)!k!}\frac{x^{n-k}}{2^k}$ by using this residue equation I found in an article, which relates the roots to each other:
$$0 = 1 + \frac{n}{z_j} + \sum\limits_{k\neq j} \frac{1}{z_k - z_j}$$
Is there a similar residue equation for the denominator $D_{pq}(x)$ of the Padé approximant for $e^{-x}?$
$$D_{pq}(x) = \sum\limits_{j=0}^q\frac{(p+q-j)!q!}{(p+q)!j!(q-j)!}x^j$$
Both this and the Bessel polynomials are hypergeometric functions that are solutions to a differential equation, but I don't know how to translate that into a residue equation.