It appears to be wellknown that $\tanh(x)\le \mathrm{erf}(x)$ on $[0,\infty)$. It's offhandedly mentioned here, for example. Where can I find a formal proof? On the one hand, it's hard to imagine that a "classic" like this wouldn't have been proven already. On the other hand, the Taylor expansions are somewhat involved (tanh involves Bernoulli numbers) and unfortunately, the inequality does not hold termwise in the expansions  so it's certainly far from obvious.

Let $f(x)={\mathrm{erf}}(x)\tanh(x)$. It can be easily seen from Taylor series at $0$ and from asymptotics at $\infty$ that $f(x)>0$ for small $x$ and for large $x$. Let us prove that $f(x)>0$ by contradiction. Suppose that $f(x)$ is negative for some $x$, then $f'$ must have at least $3$ positive zeros, by Rolle's theorem. This means that the equation $$g(x):=e^{x^2}(e^{2x}+2+e^{2x})=2\sqrt{\pi}$$ has at least $3$ positive solutions. But this is not the case because the LHS is monotone. Indeed, differentiating $g$, dividing by $e^{x^2}$ and replacing $2x$ with $y$ we obtain $$g'(x)=\sinh(y)y\cosh(y)y<0,$$ because $\sinh(y) < y \cosh(y)$ as you can see from their Taylor series. 


First, $$\begin{align} 1\mathrm{erf}(x) &= \frac{2}{\sqrt{\pi}}\int_x^\infty e^{t^2}dt, \cr 1\tanh(x) &= \int_x^\infty \mathrm{sech}^2 t\;dt . \end{align}$$ Subtract: $$ \mathrm{erf}(x)\mathrm{tanh}(x) = \int_x^\infty \left(\mathrm{sech}^2 t  \frac{2}{\sqrt{\pi}}e^{t^2}\right)dt $$ So it suffices to show that this integrand is positive. It is positive for $t>1$ (proof needed), so we establish $\mathrm{erf}(x) > \mathrm{tanh}(x)$ for $x > 1$. 

