A series about $n!$ and Riemann zeta function
Compute $$ \sum_{n=1}^{\infty}{\left( \frac{n^n}{n!e^n}-\frac{1}{\sqrt{2\pi n}} \right)}. $$ By the software Mathematica, I find $$ \sum_{n=1}^{\infty}{\left( \frac{n^n}{n!e^n}-\frac{1}{\sqrt{2\pi n}} \right)}=-\frac{2}{3}-\frac{\zeta \left( 1/2 \right)}{\sqrt{2\pi}}. $$
Well, $-\frac{1}{\sqrt{2\pi}}\zeta\left(\tfrac{1}{2}\right)$ is the $\zeta$-regularization of the divergent series $\sum_{n\geq 1}\frac{1}{\sqrt{2\pi n}}$, hence the problem boils down to finding the $\zeta$-regularization of the divergent series $\sum_{n\geq 1}\frac{n^n}{n!e^n}$. As pointed out in the comments,
$$ W(x) = \sum_{n\geq 1}\frac{n^{n-1}(-1)^{n-1}}{n!}x^n $$
holds for any $x\in\left(-\frac{1}{e},\frac{1}{e}\right)$ by Lagrange inversion theorem, hence
$$ ze^{-z} W'(-ze^{-z})=\sum_{n\geq 1}\frac{n^n}{n!e^{nz}}z^{n}=\frac{z}{1-z} =\sum_{n\geq 1}z^n\tag{1}$$
holds for any $z\in(-W(e^{-1}),1)$. Pretty strange identity, I can give you that.
Similarly, over the same interval
$$ -W(-z e^{-z})=\sum_{n\geq 1}\frac{n^{n-1}}{n!e^{nz}}z^n = z \tag{2}$$
$$ 1=\sum_{n\geq 1}\frac{n^{n}}{n!e^{nz}}z^{n-1}-\sum_{n\geq 1}\frac{n^{n}}{n!e^{nz}}z^{n}=\frac{1}{1-z}-\frac{z}{1-z}.\tag{3}$$
Since $\zeta(0)=-\frac{1}{2}$, it should not be difficult to prove from $(1)$ and $(2)$ that the $\zeta$-regularization of $\sum_{n\geq 1}\frac{n^n}{n!e^n}$ equals $-\frac{2}{3}$ as wanted, for instance by computing $\sum_{n\geq 1}\frac{n^{n-1-k}}{n!e^n}$ for any $k\in\mathbb{N}$:
$$ \sum_{n\geq 1}\frac{n^{n-2}}{n!e^n}=\int_{-1/e}^{1}\frac{W(x)}{x}\,dx = \frac{1}{2},\qquad \sum_{n\geq 1}\frac{n^{n-3}}{n!e^n}=-\int_{-1/e}^{1}\frac{W(x)}{x}(1+\log(-x))\,dx=\frac{5}{12} $$
$$ \sum_{n\geq 1}\frac{n^{n-4}}{n!e^n}=\frac{7}{18},\qquad \sum_{n\geq 1}\frac{n^{n-4}}{n!e^n}=\frac{1631}{4320},$$
$$ \sum_{n\geq 1}\frac{n^{n-1-k}}{n!e^n}= \frac{1}{\Gamma(k)}\int_{0}^{1}(1-x)(x-1-\log x)^{k-1}\,dx.\tag{4} $$
Indeed the substitution $x=e^{-s}$ in $(4)$ and the integral representation for the $\zeta$ function complete the proof.
Taking $$F\left(x\right)=\sum_{n\geq1}\frac{n^{n-1}}{n!e^{n}}x^{n}-\frac{1}{\sqrt{2\pi}}\sum_{n\geq1}\frac{x^{n}}{n^{3/2}}=-W\left(-\frac{x}{e}\right)-\frac{\mathrm{Li}_{3/2}\left(x\right)}{\sqrt{2\pi}},\,\left|x\right|<1$$ where $W\left(x\right)$ is the Lambert $W$ function and $\mathrm{Li}_{3/2}\left(x\right)$ is the Polylogarithm function, we obtain, differentiating both sides,that $$\sum_{n\geq1}\left(\frac{n^{n}}{n!e^{n}}-\frac{1}{\sqrt{2\pi n}}\right)x^{n-1}=-\frac{W\left(-\frac{x}{e}\right)}{x\left(W\left(-\frac{x}{e}\right)+1\right)}-\frac{\mathrm{Li}_{1/2}\left(x\right)}{x\sqrt{2\pi}}$$ so $$\sum_{n\geq1}\left(\frac{n^{n}}{n!e^{n}}-\frac{1}{\sqrt{2\pi n}}\right)=\lim_{x\rightarrow1^{-}}\left(-\frac{W\left(-\frac{x}{e}\right)}{x\left(W\left(-\frac{x}{e}\right)+1\right)}-\frac{\mathrm{Li}_{1/2}\left(x\right)}{x\sqrt{2\pi}}\right).$$ Now, we know that $$\mathrm{Li}_{v}\left(z\right)=\left(\Gamma\left(1-v\right)\left(1-z\right)^{v-1}+\zeta\left(v\right)\right)\left(1+O\left(\left|1-z\right|\right)\right),v\neq1,\,z\rightarrow1$$ and now we claim $$-\frac{W\left(-\frac{x}{e}\right)}{x\left(W\left(-\frac{x}{e}\right)+1\right)}\sim\frac{1}{\sqrt{2\left(1-x\right)}}-\frac{2}{3}$$ as $x\rightarrow1^{-}$. This is true because, since $$W\left(z\right)\sim-1+\sqrt{2ze+2}-\frac{2}{3}e\left(z+\frac{1}{e}\right)$$ as $z\rightarrow-1/e$, we have $$-\frac{W\left(-\frac{x}{e}\right)}{x\left(W\left(-\frac{x}{e}\right)+1\right)}\sim\frac{1-\sqrt{2\left(1-x\right)}+\frac{2}{3}\left(1-x\right)}{x\sqrt{2\left(1-x\right)}-\frac{2}{3}\left(1-x\right)x}$$ $$=\frac{1}{x}\left(-1+\frac{1}{\sqrt{2\left(1-x\right)}}\left(\frac{1}{1-\sqrt{2-2x}/3}\right)\right)=\frac{1}{x}\left(-1+\frac{1}{\sqrt{2\left(1-x\right)}}\sum_{k\geq0}\left(\frac{\sqrt{2-2x}}{3}\right)^{k}\right)$$ $$=\frac{1}{x}\left(-\frac{2}{3}+\frac{1}{\sqrt{2\left(1-x\right)}}+O\left(\sqrt{1-x}\right)\right)$$ then the claim.
This is a general answer to the followup question by Jack D'Aurizio.
Proposition
Let $n\in\mathbb{N}$. We have the asymptotic expansion $$ n!\sim \sqrt{2\pi n} \frac{n^n}{e^n} \left[ 1+ \frac1{12n} +\frac1{288n^2}-\frac{139}{51840n^3}-\frac{157}{2488320n^4}+\cdots \right] $$
Note that this is not a convergent series, but an asymptotic expansion. The error in the truncated series is asymptotically equal to the first omitted term. Regard the series on the right as an element of the ring of power series over rational numbers $\mathbb{Q}[[T]]$. $$ S(T)=1+ \frac1{12}T+\frac1{288}T^2-\frac{139}{51840}T^3-\frac{157}{2488320}T^4 + \cdots. $$ Consider the multiplicative inverse of $S(T)$ in $\mathbb{Q}[[T]]$. $$ S^{-1}(T)=1-\frac1{12}T+ g_2 T^2 + g_3 T^3 + g_4 T^4 + \cdots. $$ Let $Y_s(T)=\sum_{n=0}^{\infty} h_n(s) T^n \in\mathbb{Q}[s][[T]]$ be defined by $$ \left(\frac12 T^2\right)^{s-1}\sum_{n=0}^{\infty} h_n(s) T^n = \left[ \frac12 T^2 + \frac13 T^3 + \frac14 T^4+\cdots \right]^{s-1}. $$ Then we have
Theorem
$$\sum_{n=1}^{\infty} n^p\left[ \frac{n^n}{n!e^n}- \frac1{\sqrt{2\pi n}} \sum_{k=0}^p \frac{g_k}{n^k}\right]=(-2)^p p!h_{2p+1}(-p) - \frac1{\sqrt{2\pi}}\sum_{k=0}^p g_k \zeta\left(k+\frac12-p\right).$$
With $p=0$, it is the original series $$ \sum_{n=1}^{\infty} \left[\frac{n^n}{n!e^n} - \frac1{\sqrt{2\pi n}}\right]=-\frac23 - \frac{\zeta\left(\frac12\right)}{\sqrt{2\pi}} $$
With $p=1$, it gives the value of $$ \sum_{n=1}^{\infty} \left[ \frac{n^{n+1}}{n!e^n} - \sqrt{\frac{n}{2\pi}} + \frac{1}{12\sqrt{2\pi n}}\right] = -\frac 4{135} - \frac{\zeta\left(-\frac12\right)}{\sqrt{2\pi}} + \frac{\zeta\left(\frac12\right)}{12\sqrt{2\pi}}. $$