I am trying to approximate the digamma function in order to graph it in latex. I've already found an approximation of the gamma function from the Tex SE: $$ \Gamma(z)= (2.506628274631\sqrt{(1/z)} + 0.20888568(1/z)^{(1.5)} + 0.00870357(1/z)^{(2.5)} - \frac{(174.2106599*(1/z)^{(3.5)})}{25920} - \frac{(715.6423511(1/z)^{(4.5)})}{1244160)}\exp((-ln(1/z)-1)z) $$ and for the digamma, $\psi$, ive been using an approximation of the derivative: $$ \psi(x)=\frac{\ln \Gamma(x+0.0001)-\ln \Gamma(x-0.0001)}{0.0002} $$ But im wondering if there is a better way
4
$\begingroup$
$\endgroup$
4
$\begingroup$
$\endgroup$
Yes, there is a formula $$ψ(x) \approx \log(x) + \frac{1}{2x} - \frac{1}{12x^2} +\frac{1}{120x^4} -\frac{1}{252x^6} +\frac{1}{240x^8} -\frac{5}{660x^{10}}+\frac{691}{32760x^{12}} -\frac{1}{12x^{14}}$$ This is especially accurate for larger values of $x$. If $x$ is small, you can shift $x$ to a higher value using the relation $$ψ(x+1) = \frac{1}{x} + ψ(x)$$
You can find my source here, under the "Computation and Approximation" section.
Edit: Accuracy of Approximation
This series is quite accurate for $x$ in the interval $(1,\infty)$. Given the OP's bounds of $x = (0,20)$, we can find the difference between the expected and calculated results (note that $S(x)$ is the series in question)
$$ψ(1) - S(1) = 0.0674022$$
$$ψ(20) - S(20) = 4.44089×10^{-16}$$
From these results we find that we get exactly what we expected: the series is decent in this interval, but is MUCH better for larger $x$. However, we get a problem for $x = \frac{1}{2}$:
$$ψ(\frac{1}{2}) - S(\frac{1}{2}) = 1285.81$$
Obviously we are way off here. However, some context is useful. Our series blows up for $x<1$ towards $-\infty$, and so we can't expect anything useful here. Nevertheless, we can use the shift formula above to mitigate this issue, rewriting it as
$$ψ(x) = ψ(x+1) - \frac{1}{x}$$
Applying this just once for $x=\frac{1}{2}$, we get a difference of $1.9816604...$, an approximation $\approx 650$ times better. Thus, we see that we can use the series provided to get a decent approximation of the graph of $ψ(x)$ for all $x>0$.
answered Sep 19 '15 at 4:16
-
-
$\begingroup$ @dimebucker91 Sure thing! Note that this is an asymptotic series, hence why it works better for larger $x$ (and gets better and better as $x \to \infty$). I'm sure you can find series that are derived from other methods that will work for smaller values of $x$ with a little searching (if that is what you need)... it all depends on the domain you are graphing in $\endgroup$ – Brevan Ellefsen Sep 19 '15 at 4:24
-
1$\begingroup$ I'll check the bounds on the series and give you a worst case scenario on the difference between this series and the actual values... give me a minute. $\endgroup$ – Brevan Ellefsen Sep 19 '15 at 4:26
-
1$\begingroup$ @dimebucker91 I did the calculations and provided my results in an edit to my post (tripling my post length). I hope this is helpful to you! $\endgroup$ – Brevan Ellefsen Sep 19 '15 at 4:50
-
1$\begingroup$ @BrevanEllefsen, I think you are missing a "+" in your recurrence relation. $\endgroup$ – Luca Citi Jun 16 '18 at 20:58
4
$\begingroup$
$\endgroup$
Since you need an implementation for graphing, a simple approximation with relative low accuracy and optimized for absolute error seems suitable. A short search with computer program shows that $$\psi(x) \approx ln (x + a) - \frac{1}{b x}, \space a = 0.4849142940227510, \space b=1.0271785180163817$$ approximates with an absolute error of less than $0.00123$ on $[\frac{1}{2}, 2^{1024})$ when evaluated in IEEE-754 double precision. For arguments $x \lt \frac{1}{2}$, the function can be computed using the reflection formula $\psi(x) = \psi(1-x) + \pi \space cot \space(\pi (1-x))$.
