# Digamma function

For Barnes' gamma function of two variables, see double gamma function.
The color representation of the Digamma function, ${\displaystyle \psi (z)}$, in a rectangular region of the complex plane

In mathematics, the digamma function is defined as the logarithmic derivative of the gamma function:[1][2]

${\displaystyle \psi (x)={\frac {d}{dx}}\ln {\Big (}\Gamma (x){\Big )}={\frac {\Gamma '(x)}{\Gamma (x)}}.}$

It is the first of the polygamma functions.

The digamma function is often denoted as ψ0(x), ψ0(x) or ${\displaystyle \digamma }$ (after the archaic Greek letter Ϝ digamma).

## Relation to harmonic numbers

The gamma function obeys the equation

${\displaystyle \Gamma (z+1)=z\Gamma (z).}$

Taking the derivative with respect to z gives:

${\displaystyle \Gamma '(z+1)=z\Gamma '(z)+\Gamma (z)}$

Dividing by Γ(z+1) or the equivalent zΓ(z) gives:

${\displaystyle {\frac {\Gamma '(z+1)}{\Gamma (z+1)}}={\frac {\Gamma '(z)}{\Gamma (z)}}+{\frac {1}{z}}}$

or:

${\displaystyle \psi (z+1)=\psi (z)+{\frac {1}{z}}}$

Since the harmonic numbers are defined as

${\displaystyle H_{n}=\sum _{k=1}^{n}{\frac {1}{k}}}$

the digamma function is related to it by:

${\displaystyle \psi (n)=H_{n-1}-\gamma }$

where Hn is the n-th harmonic number, and γ is the Euler-Mascheroni constant. For half-integer values, it may be expressed as

${\displaystyle \psi \left(n+{\frac {1}{2}}\right)=-\gamma -2\ln(2)+\sum _{k=1}^{n}{\frac {2}{2k-1}}}$

## Integral representations

If the real part of x is positive then the digamma function has the following integral representation

${\displaystyle \psi (x)=\int \limits _{0}^{\infty }\left({\frac {e^{-t}}{t}}-{\frac {e^{-xt}}{1-e^{-t}}}\right)dt}$.

This may be written as

${\displaystyle \psi (s+1)=-\gamma +\int \limits _{0}^{1}\left({\frac {1-x^{s}}{1-x}}\right)dx}$

which follows from Euler's integral formula for the harmonic numbers.

## Series formula

Digamma can be computed in the complex plane outside negative integers (Abramowitz and Stegun 6.3.16),[1] using

${\displaystyle \psi (z+1)=-\gamma +\sum _{n=1}^{\infty }{\frac {z}{n(n+z)}}\qquad z\neq -1,-2,-3,\ldots }$

or

${\displaystyle \psi (z)=-\gamma +\sum _{n=0}^{\infty }{\frac {z-1}{(n+1)(n+z)}}=-\gamma +\sum _{n=0}^{\infty }\left({\frac {1}{n+1}}-{\frac {1}{n+z}}\right)\qquad z\neq 0,-1,-2,-3,\ldots }$

This can be utilized to evaluate infinite sums of rational functions, i.e.,

${\displaystyle \sum _{n=0}^{\infty }u_{n}=\sum _{n=0}^{\infty }{\frac {p(n)}{q(n)}},}$

where p(n) and q(n) are polynomials of n.

Performing partial fraction on un in the complex field, in the case when all roots of q(n) are simple roots,

${\displaystyle u_{n}={\frac {p(n)}{q(n)}}=\sum _{k=1}^{m}{\frac {a_{k}}{n+b_{k}}}.}$

For the series to converge,

${\displaystyle \lim _{n\to \infty }nu_{n}=0,}$

or otherwise the series will be greater than the harmonic series and thus diverges.

Hence

${\displaystyle \sum _{k=1}^{m}a_{k}=0,}$

and

{\displaystyle {\begin{aligned}\sum _{n=0}^{\infty }u_{n}&=\sum _{n=0}^{\infty }\sum _{k=1}^{m}{\frac {a_{k}}{n+b_{k}}}\\&=\sum _{n=0}^{\infty }\sum _{k=1}^{m}a_{k}\left({\frac {1}{n+b_{k}}}-{\frac {1}{n+1}}\right)\\&=\sum _{k=1}^{m}\left(a_{k}\sum _{n=0}^{\infty }\left({\frac {1}{n+b_{k}}}-{\frac {1}{n+1}}\right)\right)\\&=-\sum _{k=1}^{m}a_{k}\left(\psi (b_{k})+\gamma \right)\\&=-\sum _{k=1}^{m}a_{k}\psi (b_{k}).\end{aligned}}}

With the series expansion of higher rank polygamma function a generalized formula can be given as

${\displaystyle \sum _{n=0}^{\infty }u_{n}=\sum _{n=0}^{\infty }\sum _{k=1}^{m}{\frac {a_{k}}{(n+b_{k})^{r_{k}}}}=\sum _{k=1}^{m}{\frac {(-1)^{r_{k}}}{(r_{k}-1)!}}a_{k}\psi ^{(r_{k}-1)}(b_{k}),}$

provided the series on the left converges.

## Taylor series

The digamma has a rational zeta series, given by the Taylor series at z = 1. This is

${\displaystyle \psi (z+1)=-\gamma -\sum _{k=1}^{\infty }\zeta (k+1)\;(-z)^{k}}$,

which converges for |z| < 1. Here, ζ(n) is the Riemann zeta function. This series is easily derived from the corresponding Taylor's series for the Hurwitz zeta function.

## Newton series

The Newton series for the digamma follows from Euler's integral formula:

${\displaystyle \psi (s+1)=-\gamma -\sum _{k=1}^{\infty }{\frac {(-1)^{k}}{k}}{s \choose k}}$

where ${\displaystyle \textstyle {s \choose k}}$ is the binomial coefficient.

## Reflection formula

The digamma function satisfies a reflection formula similar to that of the Gamma function,

${\displaystyle \psi (1-x)-\psi (x)=\pi \ \cot(\pi x)}$

## Recurrence formula and characterization

The digamma function satisfies the recurrence relation

${\displaystyle \psi (x+1)=\psi (x)+{\frac {1}{x}}.}$

Thus, it can be said to "telescope" 1/x, for one has

${\displaystyle \Delta [\psi ](x)={\frac {1}{x}}}$

where Δ is the forward difference operator. This satisfies the recurrence relation of a partial sum of the harmonic series, thus implying the formula

${\displaystyle \psi (n)=H_{n-1}-\gamma }$

where γ is the Euler-Mascheroni constant.

More generally, one has

${\displaystyle \psi (x+1)=-\gamma +\sum _{k=1}^{\infty }\left({\frac {1}{k}}-{\frac {1}{x+k}}\right).}$

Actually, ψ is the only solution of the functional equation

${\displaystyle F(x+1)=F(x)+{\frac {1}{x}}}$

that is monotone on R+ and satisfies F(1) = −γ. This fact follows immediately from the uniqueness of the Γ function given its recurrence equation and convexity-restriction. This implies the useful difference equation:

${\displaystyle \psi (x+N)-\psi (x)=\sum _{k=0}^{N-1}{\frac {1}{x+k}}}$

## Some finite sums involving the digamma function

There are numerous finite summation formulas for the digamma function. Basic summation formulas, such as

${\displaystyle \sum _{r=1}^{m}\psi \left({\frac {r}{m}}\right)=-m(\gamma +\ln(m)),}$
${\displaystyle \sum _{r=1}^{m}\psi \left({\frac {r}{m}}\right)\cdot \exp {\dfrac {2\pi rki}{m}}=m\ln \left(1-\exp {\dfrac {2\pi ki}{m}}\right)\ ,\qquad \quad k\in \mathbb {Z} \ ,\qquad m\in \mathbb {N} \ ,\qquad k\neq m.}$
${\displaystyle \sum _{r=1}^{m-1}\psi \left({\frac {r}{m}}\right)\cdot \cos {\dfrac {2\pi rk}{m}}\ =\ m\ln \left(2\sin \left({\frac {k\pi }{m}}\right)\right)+\gamma \ ,\qquad \qquad \qquad k=1,2,\ldots ,m-1}$
${\displaystyle \sum _{r=1}^{m-1}\psi \left({\frac {r}{m}}\right)\cdot \sin {\dfrac {2\pi rk}{m}}={\frac {\pi }{2}}(2k-m)\,,\qquad \qquad \qquad k=1,2,\ldots ,m-1}$

are due to Gauss.[3][4] More complicated formulas, such as

${\displaystyle \sum _{r=0}^{m-1}\psi \left({\frac {2r+1}{2m}}\right)\cdot \cos {\dfrac {(2r+1)k\pi }{m}}=m\ln \left(\tan {\frac {\,\pi k\,}{2m}}\right)\,,\qquad \qquad \qquad k=1,2,\ldots ,m-1}$
${\displaystyle \sum _{r=0}^{m-1}\psi \left({\frac {2r+1}{2m}}\right)\cdot \sin {\dfrac {(2r+1)k\pi }{m}}=-{\frac {\pi m}{2}}\,,\qquad \qquad \qquad k=1,2,\ldots ,m-1}$
${\displaystyle \sum _{r=1}^{m-1}\psi \left({\frac {r}{m}}\right)\cdot \cot {\frac {\pi r}{m}}=-{\frac {\pi (m-1)(m-2)}{6}}}$
${\displaystyle \sum _{r=1}^{m-1}\psi \left({\frac {r}{m}}\right)\cdot {\frac {r}{m}}=-{\frac {\gamma }{2}}(m-1)-{\frac {m}{2}}\ln m-{\frac {\pi }{2}}\sum _{r=1}^{m-1}{\dfrac {r}{m}}\cdot \cot {\dfrac {\pi r}{m}}}$
${\displaystyle \sum _{r=1}^{m-1}\psi \left({\frac {r}{m}}\right)\cdot \cos {\dfrac {(2l+1)\pi r}{m}}=-{\frac {\pi }{m}}\sum _{r=1}^{m-1}{\frac {r\cdot \sin {\dfrac {2\pi r}{m}}}{\,\cos {\dfrac {2\pi r}{m}}-\cos {\dfrac {(2l+1)\pi }{m}}\,}},\qquad \quad l\in \mathbb {Z} }$
${\displaystyle \sum _{r=1}^{m-1}\psi \left({\frac {r}{m}}\right)\cdot \sin {\dfrac {(2l+1)\pi r}{m}}=-(\gamma +\ln 2m)\cot {\frac {(2l+1)\pi }{2m}}+\sin {\dfrac {(2l+1)\pi }{m}}\sum _{r=1}^{m-1}{\frac {\ln \sin {\dfrac {\pi r}{m}}}{\,\cos {\dfrac {2\pi r}{m}}-\cos {\dfrac {(2l+1)\pi }{m}}\,}},\qquad \quad l\in \mathbb {Z} }$
${\displaystyle \sum _{r=1}^{m-1}\psi ^{2}\!\left({\frac {r}{m}}\right)=(m-1)\gamma ^{2}+m(2\gamma +\ln 4m)\ln {m}-m(m-1)\ln ^{2}2+{\frac {\pi ^{2}(m^{2}-3m+2)}{12}}+m\sum _{l=1}^{m-1}\ln ^{2}\sin {\frac {\pi l}{m}}}$

are due to works of certain modern authors (see e.g. Appendix B in[5]).

## Gauss's digamma theorem

For positive integers r and m (r < m), the digamma function may be expressed in terms of Euler's constant and a finite number of elementary functions

${\displaystyle \psi \left({\frac {r}{m}}\right)=-\gamma -\ln(2m)-{\frac {\pi }{2}}\cot \left({\frac {r\pi }{m}}\right)+2\sum _{n=1}^{\lfloor {\frac {m-1}{2}}\rfloor }\cos \left({\frac {2\pi nr}{m}}\right)\ln \sin \left({\frac {\pi n}{m}}\right)}$

which holds, because of its recurrence equation, for all rational arguments.

## Computation and approximation

According to the Euler–Maclaurin formula applied to[6]

${\displaystyle \sum _{n=1}^{x}{\frac {1}{n}}}$

the digamma function for x, also a real number, can be approximated by

${\displaystyle \psi (x)=\ln(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}}}+O\left({\frac {1}{x^{16}}}\right)}$

which is the beginning of the asymptotical expansion of ψ(x). The full asymptotic series of this expansions is

${\displaystyle \psi (x)=\ln(x)-{\frac {1}{2x}}+\sum _{n=1}^{\infty }{\frac {\zeta (1-2n)}{x^{2n}}}=\ln(x)-{\frac {1}{2x}}-\sum _{n=1}^{\infty }{\frac {B_{2n}}{2n\,x^{2n}}}}$

where ${\displaystyle B_{k}}$ is the k-th Bernoulli number and ζ is the Riemann zeta function. Although the infinite sum converges for no x, this expansion becomes more accurate for larger values of x and any finite partial sum cut off from the full series. To compute ψ(x) for small x, the recurrence relation

${\displaystyle \psi (x+1)={\frac {1}{x}}+\psi (x)}$

can be used to shift the value of x to a higher value. Beal[7] suggests using the above recurrence to shift x to a value greater than 6 and then applying the above expansion with terms above ${\displaystyle x^{14}}$ cut off, which yields "more than enough precision" (at least 12 digits except near the zeroes).

As x goes to infinity, ψ(x) gets arbitrarily close to both ln(x−1/2) and ln(x). Going down from x+1 to x, ψ decrease by 1/x, ln(x−1/2) decreases by ln((x+1/2)/(x−1/2)), which is more than 1/x, and ln(x) decrease by ln((x+1)/x), which is less than 1/x. From this we see that for any positive x greater than 1/2,

${\displaystyle \psi (x)\in (\ln(x-1/2),\ln x)}$

or, for any positive x,

${\displaystyle \exp(\psi (x))\in (x-1/2,x).}$

The exponential ${\displaystyle \exp(\psi (x))}$ is approximately x−1/2 for large x, but gets closer to x at small x, approaching 0 at x = 0.

For x < 1, we can calculate limits based on the fact that between 1 and 2, ψ(x) ∈ [−γ, 1−γ], so

${\displaystyle \psi (x)\in (-1/x-\gamma ,1-1/x-\gamma ),x\in (0,1)}$

or

${\displaystyle \exp(\psi (x))\in (\exp(-1/x-\gamma ),e\exp(-1/x-\gamma ).}$

From the above asymptotic series for ψ, one can derive an asymptotic series for ${\displaystyle 1/\exp(\psi (x))}$. The series matches the overall behaviour well, that is, it behaves asymptotically as it should for large arguments, and has a zero of unbounded multiplicity at the origin too.

${\displaystyle {\frac {1}{\exp(\psi (x))}}={\frac {1}{x}}+{\frac {1}{2\cdot x^{2}}}+{\frac {5}{4\cdot 3!\cdot x^{3}}}+{\frac {3}{2\cdot 4!\cdot x^{4}}}+{\frac {47}{48\cdot 5!\cdot x^{5}}}-{\frac {5}{16\cdot 6!\cdot x^{6}}}+\cdots }$

This can be considered a Taylor expansion of ${\displaystyle \exp(-\psi (1/y))}$ at y = 0, but it does not converge.[8]

Another expansion is more precise for large arguments and saves computing terms of even order.

${\displaystyle \exp(\psi (x+{\tfrac {1}{2}}))=x+{\frac {1}{4!\cdot x}}-{\frac {37}{8\cdot 6!\cdot x^{3}}}+{\frac {10313}{72\cdot 8!\cdot x^{5}}}-{\frac {5509121}{384\cdot 10!\cdot x^{7}}}+O\left({\frac {1}{x^{9}}}\right)\qquad {\mbox{for }}x>1}$

## Special values

The digamma function has values in closed form for rational numbers, as a result of Gauss's digamma theorem. Some are listed below:

{\displaystyle {\begin{aligned}\psi (1)&=-\gamma \\\psi \left({\tfrac {1}{2}}\right)&=-2\ln {2}-\gamma \\\psi \left({\tfrac {1}{3}}\right)&=-{\tfrac {\pi }{2{\sqrt {3}}}}-{\tfrac {3}{2}}\ln {3}-\gamma \\\psi \left({\tfrac {1}{4}}\right)&=-{\tfrac {\pi }{2}}-3\ln {2}-\gamma \\\psi \left({\tfrac {1}{6}}\right)&=-{\tfrac {\pi }{2}}{\sqrt {3}}-2\ln {2}-{\tfrac {3}{2}}\ln(3)-\gamma \\\psi \left({\tfrac {1}{8}}\right)&=-{\tfrac {\pi }{2}}-4\ln {2}-{\frac {1}{\sqrt {2}}}\left\{\pi +\ln \left(2+{\sqrt {2}}\right)-\ln \left(2-{\sqrt {2}}\right)\right\}-\gamma .\end{aligned}}}

Moreover, by the series representation, it can easily be deduced that at the imaginary unit,

{\displaystyle {\begin{aligned}\Re \left(\psi (i)\right)&=-\gamma -\sum _{n=0}^{\infty }{\frac {n-1}{n^{3}+n^{2}+n+1}},\\\Im \left(\psi (i)\right)&=\sum _{n=0}^{\infty }{\frac {1}{n^{2}+1}}={\frac {1}{2}}+{\frac {\pi }{2}}\coth(\pi ).\end{aligned}}}

## Roots of the digamma function

The roots of the digamma function are the saddle points of the complex-valued gamma function. Thus they lie all on the real axis. The only one on the positive real axis is the unique minimum of the real-valued gamma function on R+ at ${\displaystyle x_{0}=1.461632144968\ldots }$. All others occur single between the poles on the negative axis:

{\displaystyle {\begin{aligned}x_{1}&=-0.504083008...,\\x_{2}&=-1.573498473...,\\x_{3}&=-2.610720868...,\\x_{4}&=-3.635293366...,\\&\qquad \cdots \end{aligned}}}

Already in 1881, Hermite observed[citation needed] that

${\displaystyle x_{n}=-n+{\frac {1}{\ln n}}+o\left({\frac {1}{\ln ^{2}n}}\right)}$

holds asymptotically. A better approximation of the location of the roots is given by

${\displaystyle x_{n}\approx -n+{\frac {1}{\pi }}\arctan \left({\frac {\pi }{\ln n}}\right)\qquad n\geq 2}$

and using a further term it becomes still better

${\displaystyle x_{n}\approx -n+{\frac {1}{\pi }}\arctan \left({\frac {\pi }{\ln n+{\frac {1}{8n}}}}\right)\qquad n\geq 1}$

which both spring off the reflection formula via

${\displaystyle 0=\psi (1-x_{n})=\psi (x_{n})+{\frac {\pi }{\tan(\pi x_{n})}}}$

and substituting ${\displaystyle \psi (x_{n})}$ by its not convergent asymptotic expansion. The correct 2nd term of this expansion is of course ${\displaystyle {\tfrac {1}{2n}}}$, where the given one works good to approximate roots with small index n.

Regarding the zeros, the following infinite sum identities were recently proved by Mező[9]

{\displaystyle {\begin{aligned}\sum _{n=0}^{\infty }{\frac {1}{x_{n}^{2}}}&=\gamma ^{2}+{\frac {\pi ^{2}}{2}},\\\sum _{n=0}^{\infty }{\frac {1}{x_{n}^{4}}}&=\gamma ^{4}+{\frac {\pi ^{4}}{9}}+{\frac {2}{3}}\gamma ^{2}\pi ^{2}+4\gamma \zeta (3).\end{aligned}}}

Here ${\displaystyle \gamma }$ is the Euler–Mascheroni constant.

## Regularization

The Digamma function appears in the regularization of divergent integrals

${\displaystyle \int _{0}^{\infty }{\frac {dx}{x+a}},}$

this integral can be approximated by a divergent general Harmonic series, but the following value can be attached to the series

${\displaystyle \sum _{n=0}^{\infty }{\frac {1}{n+a}}=-\psi (a).}$

8. ^ If it converged to a function f(y) then ln(f(y)/y) would have the same Maclaurin series as ${\displaystyle \ln(1/y)-\phi (1/y).}$ But this does not converge because the series given earlier for φ(x) does not converge.