I was given a function $f(x)=\mbox{Li}_{-n}(x)$, where Li is the polylogarithm of order $-n$ ($n>0\in\mathbb{N}$) and $x\in(-\infty,0)$. The function in this domain is bounded and has some extremes. I would like o obtain analytically the extreme (maximum/minimum) of this function. I have tried to use the property $$\frac{d}{dx}\mbox{Li}_{-n}(x)=\frac{1}{x}\mbox{Li}_{-n-1}(x)=0$$ but I got stacked, since I don't know how to solve the zeros of the polylogarithm.
Note: Just in case it is useful, I have checked with numerics that the maximum scales as $\sim n!$, but I would like to obtain it more analytically.
We have: $$\begin{eqnarray*} \text{Li}_{-k}(x) = \sum_{n\geq 1}n^k x^n = \left.\frac{d^k}{dt^k}\sum_{n\geq 1}e^{nt} x^n\right|_{t=0}&=&\sum_{j=0}^{k}j!{k+1\brace j+1}\left(\frac{x}{1-x}\right)^{j+1}\\[0.2cm]&=&\frac{1}{(1-x)^{k+1}}\sum_{j=0}^{k-1}\left\langle k\atop j\right\rangle x^{k-j}\end{eqnarray*} $$ as a consequence of Worpitzky's identity, too. Since $f(x)=\frac{x}{1-x}$ is a bijective function from $\mathbb{R}^-$ to $(-1,0)$, in order to compute the stationary points of $\text{Li}_{-k}(x)$ over $\mathbb{R}^-$ it is enough to compute the stationary points of: $$ g_k(x) = \sum_{j=0}^{k}j!{k+1\brace j+1}x^{j+1}=\text{Li}_{-k}\left(\frac{x}{x+1}\right) $$ over $(-1,0)$. Dobinski's formula gives: $$B_k(z)=\sum_{j=0}^{k}{k \brace j}z^j $$ where the LHS is a Bell polynomial. So we have: $$ B_{k+1}(xz) = \sum_{j=0}^{k+1}{k+1\brace j+1} z^{j+1} x^{j+1} $$ and: $$ g_k(x) = -(k+1)!x^{k+2}+\int_{0}^{+\infty}\frac{B_{k+1}(xz)}{ze^{z}}\,dz.$$ It follows that a more accurate approximation for $\sup_{u<0}\left|\text{Li}_{-k}(u)\right|$ should be given in terms of Bell numbers and complementary Bell numbers. Work in progress.