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.
Hmm. Let's see if I can muddle my way through this one, too.
Grünbaum has an article stating that for polynomials with simple roots
In other words,
$$\frac{P''(x_k)}{2P'(x_k)} = \sum\limits_{j \neq k}\frac{1}{x_k-x_j}$$
at each root $x_k$.
According to this Wikipedia page on the Padé table which has $R_{m,n}(z) = f_{m,n}(z)/g_{m,n}(z)$ with $g_{m,n}(z) = {}_1F_1(-n;-m-n;-z)$for the approximation of $e^{x}$,
$$D_{pq}(x) = g_{p,q}(-x) = {}_1F_1(-q;-p-q;x)$$
and the Wikipedia page on ${}_1F_1$ has $w={}_1F_1(a;b;z)$ satisfying the differential equation
$$z \frac{d^2w}{dz^2} + (b-z)\frac{dw}{dz} - aw = 0$$
and at the roots where $w=0$,
$$z \frac{d^2w}{dz^2} + (b-z)\frac{dw}{dz} = 0$$
so that
$$ \frac{\frac{d^2w}{dz^2}}{2\frac{dw}{dz}} = \frac{z-b}{2z} $$
and with the formula in Grünbaum's article, that means
$$\frac{x_k-b}{2x_k} = \sum\limits_{j\neq k}\frac{1}{x_k-x_j}$$
or
$$0 = -\frac{1}{2} + \frac{b}{2x_k} + \sum\limits_{j\neq k}\frac{1}{x_k-x_j}$$
for each root $x_k$, where $b = -p-q$.