Hubbry Logo
search
logo
2210702

Digamma function

logo
Community Hub0 Subscribers
Read side by side
from Wikipedia
The digamma function ,
visualized using domain coloring
Plots of the digamma and the next three polygamma functions along the real line (they are real-valued on the real line)

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

It is the first of the polygamma functions. This function is strictly increasing and strictly concave on ,[4] and it asymptotically behaves as[5]

for complex numbers with large modulus () in the sector for any .

The digamma function is often denoted as or Ϝ[6] (the uppercase form of the archaic Greek consonant digamma meaning double-gamma).

Relation to harmonic numbers

[edit]

The gamma function obeys the equation

Taking the logarithm on both sides and using the functional equation property of the log-gamma function gives:

Differentiating both sides with respect to z gives:

Since the harmonic numbers are defined for positive integers n as

the digamma function is related to them by

where H0 = 0, and γ is the Euler–Mascheroni constant. For half-integer arguments the digamma function takes the values

Integral representations

[edit]

If the real part of z is positive then the digamma function has the following integral representation due to Gauss:[7]

Combining this expression with an integral identity for the Euler–Mascheroni constant gives:

The integral is Euler's harmonic number , so the previous formula may also be written

A consequence is the following generalization of the recurrence relation:

An integral representation due to Dirichlet is:[7]

Gauss's integral representation can be manipulated to give the start of the asymptotic expansion of .[8]

This formula is also a consequence of Binet's first integral for the gamma function. The integral may be recognized as a Laplace transform.

Binet's second integral for the gamma function gives a different formula for which also gives the first few terms of the asymptotic expansion:[9]

From the definition of and the integral representation of the gamma function, one obtains

with .[10]

Infinite product representation

[edit]

The function is an entire function,[11] and it can be represented by the infinite product

Here is the kth zero of (see below), and is the Euler–Mascheroni constant.

Note: This is also equal to due to the definition of the digamma function: .

Series representation

[edit]

Series formula

[edit]

Euler's product formula for the gamma function, combined with the functional equation and an identity for the Euler–Mascheroni constant, yields the following expression for the digamma function, valid in the complex plane outside the negative integers (Abramowitz and Stegun 6.3.16):[1]

Equivalently,

Evaluation of sums of rational functions

[edit]

The above identity can be used to evaluate sums of the form

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,

For the series to converge,

otherwise the series will be greater than the harmonic series and thus diverge. Hence

and

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

provided the series on the left converges.

Taylor series

[edit]

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

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

[edit]

The Newton series for the digamma, sometimes referred to as Stern series, derived by Moritz Abraham Stern in 1847,[12][13][14] reads

where (s
k
)
is the binomial coefficient. It may also be generalized to

where m = 2, 3, 4, ...[13]

Series with Gregory's coefficients, Cauchy numbers and Bernoulli polynomials of the second kind

[edit]

There exist various series for the digamma containing rational coefficients only for the rational arguments. In particular, the series with Gregory's coefficients Gn is

where (v)n is the rising factorial (v)n = v(v+1)(v+2) ... (v+n-1), Gn(k) are the Gregory coefficients of higher order with Gn(1) = Gn, Γ is the gamma function and ζ is the Hurwitz zeta function.[15][13] Similar series with the Cauchy numbers of the second kind Cn reads[15][13]

A series with the Bernoulli polynomials of the second kind has the following form[13]

where ψn(a) are the Bernoulli polynomials of the second kind defined by the generating equation

It may be generalized to

where the polynomials Nn,r(a) are given by the following generating equation

so that Nn,1(a) = ψn(a).[13] Similar expressions with the logarithm of the gamma function involve these formulas[13]

and

where and .

Reflection formula

[edit]

The digamma and polygamma functions satisfy reflection formulas similar to that of the gamma function:

.
.
.

Recurrence formula and characterization

[edit]

The digamma function satisfies the recurrence relation

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

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

where γ is the Euler–Mascheroni constant.

Actually, ψ is the only solution of the functional equation

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

Some finite sums involving the digamma function

[edit]

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

are due to Gauss.[16][17] More complicated formulas, such as

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

We also have [19]

Gauss's digamma theorem

[edit]

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[20] [21]

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

Multiplication theorem

[edit]

The multiplication theorem of the -function is equivalent to[22]

Asymptotic expansion

[edit]

The digamma function has the asymptotic expansion

where Bk is the kth Bernoulli number and ζ is the Riemann zeta function. The first few terms of this expansion are:

Although the infinite sum does not converge for any z, any finite partial sum becomes increasingly accurate as z increases.

The expansion can be found by applying the Euler–Maclaurin formula to the sum[23]

The expansion can also be derived from the integral representation coming from Binet's second integral formula for the gamma function. Expanding as a geometric series and substituting an integral representation of the Bernoulli numbers leads to the same asymptotic series as above. Furthermore, expanding only finitely many terms of the series gives a formula with an explicit error term:

Inequalities

[edit]

When x > 0, the function

is completely monotonic and in particular positive. This is a consequence of Bernstein's theorem on monotone functions applied to the integral representation coming from Binet's first integral for the gamma function. Additionally, by the convexity inequality , the integrand in this representation is bounded above by . Consequently

is also completely monotonic. It follows that, for all x > 0,

This recovers a theorem of Horst Alzer.[24] Alzer also proved that, for s ∈ (0, 1),

Related bounds were obtained by Elezovic, Giordano, and Pecaric, who proved that, for x > 0 ,

where is the Euler–Mascheroni constant.[25] The constants ( and ) appearing in these bounds are the best possible.[26]

The mean value theorem implies the following analog of Gautschi's inequality: If x > c, where c ≈ 1.461 is the unique positive real root of the digamma function, and if s > 0, then

Moreover, equality holds if and only if s = 1.[27]

Inspired by the harmonic mean value inequality for the classical gamma function, Horzt Alzer and Graham Jameson proved, among other things, a harmonic mean-value inequality for the digamma function:

for

Equality holds if and only if .[28]

Computation and approximation

[edit]

The asymptotic expansion gives an easy way to compute ψ(x) when the real part of x is large. To compute ψ(x) for small x, the recurrence relation

can be used to shift the value of x to a higher value. Beal[29] suggests using the above recurrence to shift x to a value greater than 6 and then applying the above expansion with terms above x14 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(x1/2) and ln x. Going down from x + 1 to x, ψ decreases by 1/x, ln(x1/2) decreases by ln(x + 1/2) / (x1/2), which is more than 1/x, and ln x decreases by ln(1 + 1/x), which is less than 1/x. From this we see that for any positive x greater than 1/2,

or, for any positive x,

The exponential exp ψ(x) is approximately x1/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

or

From the above asymptotic series for ψ, one can derive an asymptotic series for exp(−ψ(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.

This is similar to a Taylor expansion of exp(−ψ(1 / y)) at y = 0, but it does not converge.[30] (The function is not analytic at infinity.) A similar series exists for exp(ψ(x)) which starts with

If one calculates the asymptotic series for ψ(x+1/2) it turns out that there are no odd powers of x (there is no x−1 term). This leads to the following asymptotic expansion, which saves computing terms of even order.

Similar in spirit to the Lanczos approximation of the -function is Spouge's approximation.

Another alternative is to use the recurrence relation or the multiplication formula to shift the argument of into the range and to evaluate the Chebyshev series there.[31][32]

Special values

[edit]

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

Moreover, by taking the logarithmic derivative of or where is real-valued, it can easily be deduced that

Apart from Gauss's digamma theorem, no such closed formula is known for the real part in general. We have, for example, at the imaginary unit the numerical approximation OEISA248177

Roots of the digamma function

[edit]

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 x0 = 1.46163214496836234126.... All others occur single between the poles on the negative axis:

x1 = −0.50408300826445540925...
x2 = −1.57349847316239045877...
x3 = −2.61072086844414465000...
x4 = −3.63529336643690109783...

Already in 1881, Charles Hermite observed[33] that

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

and using a further term it becomes still better

which both spring off the reflection formula via

and substituting ψ(xn) by its not convergent asymptotic expansion. The correct second term of this expansion is 1/2n, where the given one works well to approximate roots with small n.

Another improvement of Hermite's formula can be given:[11]

Regarding the zeros, the following infinite sum identities were recently proved by István Mező and Michael Hoffman[11][34]

In general, the function

can be determined and it is studied in detail by the cited authors.

The following results[11]

also hold true.

Regularization

[edit]

The digamma function appears in the regularization of divergent integrals

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

In applied mathematics

[edit]

Many notable probability distributions use the gamma function in the definition of their probability density or mass functions. Then in statistics when doing maximum likelihood estimation on models involving such distributions, the digamma function naturally appears when the derivative of the log-likelihood is taken for finding the maxima.

See also

[edit]
  • Polygamma function
  • Trigamma function
  • Chebyshev expansions of the digamma function in Wimp, Jet (1961). "Polynomial approximations to integral transforms". Math. Comp. 15 (74): 174–178. doi:10.1090/S0025-5718-61-99221-3.

References

[edit]
[edit]
Revisions and contributorsEdit on WikipediaRead on Wikipedia
from Grokipedia
The digamma function, denoted ψ(z)\psi(z), is a special function in complex analysis defined as the logarithmic derivative of the gamma function Γ(z)\Gamma(z), specifically ψ(z)=Γ(z)Γ(z)\psi(z) = \frac{\Gamma'(z)}{\Gamma(z)} for zC{0,1,2,}z \in \mathbb{C} \setminus \{0, -1, -2, \dots \}.[1] This definition extends the concept of the derivative of the logarithm, making ψ(z)\psi(z) meromorphic across the complex plane with simple poles of residue 1-1 at the non-positive integers.[1] The digamma function serves as the basis for the polygamma functions of order n1n \geq 1, defined as the successive derivatives ψ(n)(z)=dndznψ(z)\psi^{(n)}(z) = \frac{d^n}{dz^n} \psi(z), where ψ(1)(z)=ψ(z)\psi^{(1)}(z) = \psi'(z) is the trigamma function. It captures essential behaviors of the gamma function's growth and oscillation.[2] Key properties of the digamma function include its reflection formula ψ(1z)ψ(z)=πcot(πz)\psi(1 - z) - \psi(z) = \pi \cot(\pi z) and the recurrence relation ψ(z+1)=ψ(z)+1z\psi(z+1) = \psi(z) + \frac{1}{z}, which facilitate computation and analytic continuation. For positive integers nn, it relates directly to the harmonic numbers Hn=k=1n1kH_n = \sum_{k=1}^n \frac{1}{k} via ψ(n+1)=Hnγ\psi(n+1) = H_n - \gamma, where γ0.57721\gamma \approx 0.57721 is the Euler-Mascheroni constant, and ψ(1)=γ\psi(1) = -\gamma.[3] Series representations, such as ψ(z)=γ+k=0(1k+11k+z)\psi(z) = -\gamma + \sum_{k=0}^\infty \left( \frac{1}{k+1} - \frac{1}{k+z} \right) for z>0\Re z > 0, and asymptotic expansions like ψ(z)lnz12zk=1B2k2kz2k\psi(z) \sim \ln z - \frac{1}{2z} - \sum_{k=1}^\infty \frac{B_{2k}}{2k z^{2k}} as z|z| \to \infty in argz<π|\arg z| < \pi, highlight its utility in approximations and limits. These features underscore its foundational role in special functions theory. The digamma function finds broad applications across mathematics and related fields, including number theory through connections to the Riemann zeta function via ψ(1)=ζ(2)=π26\psi'(1) = \zeta(2) = \frac{\pi^2}{6}, and in statistics where it appears in the expected value of the logarithm of gamma-distributed variables, such as E[lnW]=ψ(α)\mathbb{E}[\ln W] = \psi(\alpha) for WΓ(α,1)W \sim \Gamma(\alpha, 1).[2] Its derivative, the trigamma function, computes variances like Var(lnW)=ψ(α)\mathrm{Var}(\ln W) = \psi'(\alpha), aiding inference in gamma and beta distributions.[4] In physics, it supports models in statistical mechanics and stellar evolution by providing exact expressions for partition functions and energy distributions derived from hypergeometric forms.[5] Additionally, generalizations such as the k-digamma function, defined via the k-gamma function, extend its applications to studies involving complete monotonicity and infinite series analysis.[6]

Definition and Basic Properties

Definition

The digamma function, commonly denoted by ψ(z)\psi(z), is defined as the logarithmic derivative of the gamma function Γ(z)\Gamma(z):
ψ(z)=ddzlnΓ(z)=Γ(z)Γ(z), \psi(z) = \frac{d}{dz} \ln \Gamma(z) = \frac{\Gamma'(z)}{\Gamma(z)},
where zz is a complex number not equal to a non-positive integer. This definition positions the digamma function as a key analytic tool for studying properties of the gamma function through differentiation. The digamma function was first studied by James Stirling in 1730 and later contributed to by Carl Friedrich Gauss in 1813 through his investigations into series expansions related to the gamma function, marking it as the first in the family of polygamma functions.[7] The digamma function naturally arises in the differentiation of expressions involving the gamma function, such as in the analysis of factorial generalizations and their applications in probability distributions and special function theory.[8]

Relation to Gamma function

The digamma function ψ(z)\psi(z) is defined as the logarithmic derivative of the gamma function Γ(z)\Gamma(z), given by
ψ(z)=ddzlnΓ(z)=Γ(z)Γ(z). \psi(z) = \frac{d}{dz} \ln \Gamma(z) = \frac{\Gamma'(z)}{\Gamma(z)}.
This representation holds for Rez>0\operatorname{Re} z > 0, with analytic continuation to the meromorphic function on the complex plane excluding the non-positive integers. The Weierstrass canonical product form of the gamma function implies that Γ(z)\Gamma(z) has simple poles at z=0,1,2,z = 0, -1, -2, \dots, and consequently ψ(z)\psi(z) inherits simple poles at these points, each with residue 1-1. From the functional equation of the gamma function, Γ(z+1)=zΓ(z)\Gamma(z+1) = z \Gamma(z), taking the logarithmic derivative on both sides yields the basic recurrence relation
ψ(z+1)=ψ(z)+1z, \psi(z+1) = \psi(z) + \frac{1}{z},
valid for z0,1,2,z \neq 0, -1, -2, \dots. This relation highlights the digamma function's role in extending properties of the gamma function to difference equations. The digamma function serves as the zeroth-order case in the family of polygamma functions, denoted ψ(0)(z)=ψ(z)\psi^{(0)}(z) = \psi(z). Higher-order polygamma functions are defined as the successive derivatives,
ψ(m)(z)=dmdzmψ(z),m1, \psi^{(m)}(z) = \frac{d^m}{dz^m} \psi(z), \quad m \geq 1,
which are also meromorphic with poles of order m+1m+1 at the non-positive integers.[2] The digamma function is uniquely characterized as the meromorphic solution to the recurrence ψ(z+1)=ψ(z)+1/z\psi(z+1) = \psi(z) + 1/z that satisfies the asymptotic growth condition ψ(z)lnz1/(2z)\psi(z) \sim \ln z - 1/(2z) as z|z| \to \infty in argzπδ|\arg z| \leq \pi - \delta for any fixed δ>0\delta > 0. This uniqueness parallels the Bohr–Mollerup characterization of the gamma function itself via its recurrence and log-convexity.

Connection to harmonic numbers

The digamma function provides a natural extension of the harmonic numbers to non-integer arguments. For positive integers nn, the digamma function evaluates to ψ(n+1)=γ+Hn\psi(n+1) = -\gamma + H_n, where Hn=k=1n1kH_n = \sum_{k=1}^n \frac{1}{k} is the nnth harmonic number and γ\gamma is the Euler-Mascheroni constant.
\] This relation follows from the recurrence property of the digamma function and the limiting definition of $\gamma = \lim_{m \to \infty} (H_m - \ln m)$.\[
Equivalently, γ=ψ(1)\gamma = -\psi(1), a characterization that links the constant directly to the digamma function at the origin of the positive reals.
\] Leonhard Euler introduced $\gamma$ in 1734 while studying the divergence of the harmonic series, recognizing its role in bridging discrete sums and continuous logarithms.\[
This integer case generalizes to complex arguments zz with Re(z)>0\operatorname{Re}(z) > 0 through the series representation ψ(z+1)=γ+k=1(1k1k+z)\psi(z+1) = -\gamma + \sum_{k=1}^\infty \left( \frac{1}{k} - \frac{1}{k+z} \right). This expression interpolates the harmonic numbers, as the partial sum up to NN approximates HNk=1N1k+zH_N - \sum_{k=1}^N \frac{1}{k+z}, which converges to the full form as NN \to \infty. The series embodies the conceptual shift from finite harmonic sums to an infinite difference that defines the digamma for non-integers, preserving the asymptotic behavior near integers. A useful consequence arises in differences of digamma values: for Re(a)>0\operatorname{Re}(a) > 0 and Re(b)>0\operatorname{Re}(b) > 0, ψ(a)ψ(b)=k=0(1b+k1a+k)\psi(a) - \psi(b) = \sum_{k=0}^\infty \left( \frac{1}{b+k} - \frac{1}{a+k} \right). This formula, derived from the series expansion of ψ\psi, facilitates evaluations in applications such as summation identities and integral approximations where harmonic-like differences appear. It highlights the digamma's role in unifying discrete and continuous harmonic structures.

Mathematical Representations

Integral representations

One of the fundamental integral representations of the digamma function ψ(z)\psi(z) for Re(z)>0\operatorname{Re}(z) > 0 is given by
ψ(z)=0(ettezt1et)dt. \psi(z) = \int_0^\infty \left( \frac{e^{-t}}{t} - \frac{e^{-zt}}{1 - e^{-t}} \right) \, dt.
This form arises from regularizing the divergent integral 0et/tdt\int_0^\infty e^{-t}/t \, dt against the geometric series expansion of 1/(1et)1/(1 - e^{-t}), providing a principal value interpretation that converges in the specified half-plane.
An alternative representation, also valid for Re(z)>0\operatorname{Re}(z) > 0, expresses the digamma function in terms of an integral over the unit interval:
ψ(z)=γ+011tz11tdt, \psi(z) = -\gamma + \int_0^1 \frac{1 - t^{z-1}}{1 - t} \, dt,
where γ\gamma is the Euler-Mascheroni constant. This form is particularly useful for computational purposes when zz is close to positive integers, as the integrand simplifies near those points.
The digamma function can also be derived directly from the integral representation of the gamma function Γ(z)=0tz1etdt\Gamma(z) = \int_0^\infty t^{z-1} e^{-t} \, dt for Re(z)>0\operatorname{Re}(z) > 0 by differentiating under the integral sign with respect to zz, yielding
ψ(z)=Γ(z)Γ(z)=1Γ(z)0tz1etlntdt. \psi(z) = \frac{\Gamma'(z)}{\Gamma(z)} = \frac{1}{\Gamma(z)} \int_0^\infty t^{z-1} e^{-t} \ln t \, dt.
This expression highlights the connection to the logarithmic derivative and is justified by the dominated convergence theorem applied to the parameter differentiation, ensuring validity in the right half-plane.
For analytic continuation to the entire complex plane except non-positive integers, a contour integral representation using the Hankel contour HH (which starts at ++\infty, proceeds to the origin just above the positive real axis, encircles the origin counterclockwise, and returns to ++\infty just below the positive real axis) is
ψ(z)=Hπcot(πζ)z+ζdζ. \psi(z) = \int_H \frac{\pi \cot(\pi \zeta)}{z + \zeta} \, d\zeta.
[9] This form captures the poles of cot(πζ)\cot(\pi \zeta) at integers, reproducing the recurrence properties of ψ(z)\psi(z), and extends the function meromorphically.

Series representations

The digamma function admits a fundamental series representation
ψ(z+1)=γ+k=1(1k1k+z), \psi(z+1) = -\gamma + \sum_{k=1}^{\infty} \left( \frac{1}{k} - \frac{1}{k+z} \right),
valid for Re(z)>1\operatorname{Re}(z) > -1 and z0,1,2,z \neq 0, -1, -2, \dots. This form arises from the Weierstrass product for the gamma function and is useful for numerical evaluation away from poles. An equivalent expression is
ψ(z)=γ1z+k=1zk(k+z), \psi(z) = -\gamma - \frac{1}{z} + \sum_{k=1}^{\infty} \frac{z}{k(k+z)},
also holding under the same conditions. The Taylor series expansion of the digamma function about z=1z=1 is
ψ(1+z)=γ+k=2(1)kζ(k)zk1, \psi(1+z) = -\gamma + \sum_{k=2}^{\infty} (-1)^{k} \zeta(k) z^{k-1},
convergent for z<1|z| < 1, where ζ(k)\zeta(k) denotes the Riemann zeta function at positive integers k2k \geq 2. This series reflects the analytic continuation and provides coefficients directly tied to zeta values, facilitating approximations near z=1z=1. Equivalently, shifting the index yields
ψ(1+z)=γ+n=1(1)n+1ζ(n+1)zn, \psi(1+z) = -\gamma + \sum_{n=1}^{\infty} (-1)^{n+1} \zeta(n+1) z^{n},
with the same radius of convergence. The Newton series, also known as the Stern series, expresses the digamma function in the binomial basis as an interpolation form using forward differences:
ψ(s+1)=γk=1(1)kk(sk), \psi(s+1) = -\gamma - \sum_{k=1}^{\infty} \frac{(-1)^{k}}{k} \binom{s}{k},
valid for s<1|s| < 1. This representation leverages the recurrence ψ(s+2)ψ(s+1)=1/(s+1)\psi(s+2) - \psi(s+1) = 1/(s+1) and converges to the Taylor expansion in the binomial transform.[10] Another series expansion involves the Bernoulli polynomials of the second kind bn(x)b_{n}(x), defined via the generating function tet1ext=n=0bn(x)tnn!\frac{t}{\mathrm{e}^{t}-1} \mathrm{e}^{xt} = \sum_{n=0}^{\infty} b_{n}(x) \frac{t^{n}}{n!}:
ψ(x)=γ1x+k=1(1)k+1kbk1(x), \psi(x) = -\gamma - \frac{1}{x} + \sum_{k=1}^{\infty} \frac{(-1)^{k+1}}{k} b_{k-1}(x),
for x>0x > 0. The Bernoulli polynomials of the second kind connect to Gregory coefficients Gn=(1)nbn(0)G_{n} = (-1)^{n} b_{n}(0) and Cauchy numbers of the first kind, providing a basis for expressing the digamma in terms of signed stirling-like structures. These series representations enable evaluation of certain rational sums; for example, the infinite sum k=11k(k+m)=ψ(m+1)+γm\sum_{k=1}^{\infty} \frac{1}{k(k+m)} = \frac{\psi(m+1) + \gamma}{m} for positive integer mm, derived by partial fractions and the general series form.

Infinite product representation

The infinite product representation of the digamma function ψ(z)\psi(z) is derived from the Weierstrass form of the gamma function Γ(z)\Gamma(z). The reciprocal of the gamma function is expressed as
1Γ(z)=zeγzk=1(1+zk)ez/k, \frac{1}{\Gamma(z)} = z e^{\gamma z} \prod_{k=1}^\infty \left(1 + \frac{z}{k}\right) e^{-z/k},
where γ\gamma is the Euler-Mascheroni constant. Taking the natural logarithm yields
ln(1Γ(z))=lnz+γz+k=1[ln(1+zk)zk]. \ln\left(\frac{1}{\Gamma(z)}\right) = \ln z + \gamma z + \sum_{k=1}^\infty \left[ \ln\left(1 + \frac{z}{k}\right) - \frac{z}{k} \right].
Differentiating both sides with respect to zz gives
Γ(z)Γ(z)=1z+γ+k=1[1k+z1k], -\frac{\Gamma'(z)}{\Gamma(z)} = \frac{1}{z} + \gamma + \sum_{k=1}^\infty \left[ \frac{1}{k + z} - \frac{1}{k} \right],
and since ψ(z)=Γ(z)/Γ(z)\psi(z) = \Gamma'(z)/\Gamma(z), it follows that
ψ(z)=γ1z+k=1[1k1k+z]. \psi(z) = -\gamma - \frac{1}{z} + \sum_{k=1}^\infty \left[ \frac{1}{k} - \frac{1}{k + z} \right].
This series converges for all complex zz except the non-positive integers. An equivalent form is obtained by rewriting the terms in the sum:
k=1[1k1k+z]=k=1zk(k+z), \sum_{k=1}^\infty \left[ \frac{1}{k} - \frac{1}{k + z} \right] = \sum_{k=1}^\infty \frac{z}{k(k + z)},
so
ψ(z)=γ+k=0[1k+11k+z]=γ+k=0z(k+1)(k+z). \psi(z) = -\gamma + \sum_{k=0}^\infty \left[ \frac{1}{k+1} - \frac{1}{k + z} \right] = -\gamma + \sum_{k=0}^\infty \frac{z}{(k+1)(k + z)}.
This representation also holds for Rez>0\operatorname{Re} z > 0, with analytic continuation to the rest of the complex plane. The Weierstrass product ensures the meromorphicity of Γ(z)\Gamma(z) on the complex plane, with simple poles at the non-positive integers and no zeros, implying that ψ(z)\psi(z) is meromorphic with simple poles at these points, each of residue 1-1.

Functional Equations

Recurrence relation

The digamma function satisfies the functional recurrence relation
ψ(z+1)=ψ(z)+1z, \psi(z+1) = \psi(z) + \frac{1}{z},
which follows directly from differentiating the recurrence equation Γ(z+1)=zΓ(z)\Gamma(z+1) = z \Gamma(z) for the gamma function and using the definition ψ(z)=Γ(z)Γ(z)\psi(z) = \frac{\Gamma'(z)}{\Gamma(z)}. This relation holds for all complex zz not equal to a non-positive integer, where ψ(z)\psi(z) has simple poles. By iterated application of the recurrence, a generalization holds for any positive integer nn:
ψ(z+n)=ψ(z)+k=0n11z+k. \psi(z + n) = \psi(z) + \sum_{k=0}^{n-1} \frac{1}{z + k}.
This finite sum expresses the shift by nn units in terms of the original value plus harmonic-like terms, facilitating computations for arguments with positive real part.[8] The digamma function is the unique meromorphic function on the complex plane satisfying the recurrence ψ(z+1)=ψ(z)+1/z\psi(z+1) = \psi(z) + 1/z together with the asymptotic behavior
ψ(z)lnz12z \psi(z) \sim \ln z - \frac{1}{2z}
as z|z| \to \infty in the sector argz<π|\arg z| < \pi. This characterization, combining the iterative property with the leading asymptotic terms derived from Stirling's series for the gamma function, uniquely determines ψ(z)\psi(z) among meromorphic functions with simple poles at the non-positive integers. The recurrence is equivalently expressed in difference form as the first forward difference ψ(z+1)ψ(z)=1/z\psi(z+1) - \psi(z) = 1/z. Higher-order forward differences of ψ(z)\psi(z) connect to the polygamma functions of positive order, where the mm-th polygamma function ψ(m)(z)\psi^{(m)}(z) (for m1m \geq 1) satisfies ψ(m)(z+1)ψ(m)(z)=(1)m+1m!z(m+1)\psi^{(m)}(z+1) - \psi^{(m)}(z) = (-1)^{m+1} m! \, z^{-(m+1)}.

Reflection formula

The reflection formula for the digamma function relates the values at zz and 1z1 - z:
ψ(1z)ψ(z)=πcot(πz) \psi(1 - z) - \psi(z) = \pi \cot(\pi z)
for complex zz not equal to an integer. This formula is derived from the reflection formula for the gamma function,
Γ(z)Γ(1z)=πsin(πz), \Gamma(z) \Gamma(1 - z) = \frac{\pi}{\sin(\pi z)},
by taking the natural logarithm of both sides and differentiating with respect to zz:
Γ(z)Γ(z)+Γ(1z)Γ(1z)(1)=πcos(πz)sin(πz). \frac{\Gamma'(z)}{\Gamma(z)} + \frac{\Gamma'(1 - z)}{\Gamma(1 - z)} \cdot (-1) = -\frac{\pi \cos(\pi z)}{\sin(\pi z)}.
Substituting the definition ψ(w)=Γ(w)/Γ(w)\psi(w) = \Gamma'(w)/\Gamma(w) yields
ψ(z)ψ(1z)=πcot(πz), \psi(z) - \psi(1 - z) = -\pi \cot(\pi z),
which rearranges to the desired result. The reflection formula enables analytic continuation of the digamma function across the complex plane, particularly to regions near negative non-integer values, by relating values across the unit interval. It also aids in evaluating the digamma function at half-integer points; for instance, combined with the recurrence relation ψ(z+1)=ψ(z)+1/z\psi(z + 1) = \psi(z) + 1/z, it allows computation of ψ(3/2)=ψ(1/2)+2\psi(3/2) = \psi(1/2) + 2 once ψ(1/2)=γ2ln2\psi(1/2) = -\gamma - 2\ln 2 is known from other representations. Differentiating the reflection formula with respect to zz gives the corresponding relation for the trigamma function ψ1(z)=ddzψ(z)\psi_1(z) = \frac{d}{dz} \psi(z):
ψ1(1z)+ψ1(z)=π2csc2(πz). \psi_1(1 - z) + \psi_1(z) = \pi^2 \csc^2(\pi z).
This extension follows directly from the chain rule applied to the left side and the derivative of cot(πz)\cot(\pi z) on the right.[11]

Multiplication theorem

The multiplication theorem for the digamma function provides a relation between the value of ψ(mz)\psi(mz) and the average of shifted digamma values for positive integers mm. Specifically, for mz0,1,2,mz \neq 0, -1, -2, \dots,
ψ(mz)=lnm+1mk=0m1ψ(z+km). \psi(mz) = \ln m + \frac{1}{m} \sum_{k=0}^{m-1} \psi\left(z + \frac{k}{m}\right).
This formula generalizes the scaling properties of the digamma function and is valid for complex zz in the appropriate domain where the terms are defined. The theorem is derived from Gauss's multiplication formula for the gamma function, which states that for positive integers mm,
Γ(mz)=(2π)(1m)/2mmz1/2k=0m1Γ(z+km). \Gamma(mz) = (2\pi)^{(1-m)/2} m^{mz - 1/2} \prod_{k=0}^{m-1} \Gamma\left(z + \frac{k}{m}\right).
Taking the natural logarithm of both sides and differentiating with respect to zz yields ψ(mz)mΓ(mz)=ddz[ln((2π)(1m)/2mmz1/2)Γ(mz)+k=0m1lnΓ(z+km)]\psi(mz) \cdot m \Gamma(mz) = \frac{d}{dz} \left[ \ln \left( (2\pi)^{(1-m)/2} m^{mz - 1/2} \right) \cdot \Gamma(mz) + \sum_{k=0}^{m-1} \ln \Gamma\left(z + \frac{k}{m}\right) \right], which simplifies via the definition ψ(w)=Γ(w)Γ(w)\psi(w) = \frac{\Gamma'(w)}{\Gamma(w)} to the multiplication theorem after algebraic rearrangement. A special case occurs for m=2m=2, known as the duplication formula:
ψ(2z)=12ψ(z)+12ψ(z+12)+ln2. \psi(2z) = \frac{1}{2} \psi(z) + \frac{1}{2} \psi\left(z + \frac{1}{2}\right) + \ln 2.
This relation is particularly useful for computing digamma values at even multiples and follows directly from substituting m=2m=2 into the general formula. The multiplication theorem connects to higher-order generalizations, such as the multiple gamma functions introduced by Barnes, where the logarithmic derivative of the Barnes GG-function or multiple gamma yields sums involving digamma functions at fractional shifts, extending the theorem to non-integer orders.

Gauss's digamma theorem

Gauss's digamma theorem provides an explicit closed-form expression for the value of the digamma function ψ(z)\psi(z) at rational arguments z=p/qz = p/q, where pp and qq are positive integers with 1p<q1 \leq p < q and gcd(p,q)=1\gcd(p, q) = 1. The theorem states that
ψ(pq)=γlog(2q)π2cot(πpq)+2j=1q/2cos(2πpjq)logsin(πjq), \psi\left( \frac{p}{q} \right) = -\gamma - \log(2q) - \frac{\pi}{2} \cot\left( \frac{\pi p}{q} \right) + 2 \sum_{j=1}^{\lfloor q/2 \rfloor} \cos\left( \frac{2\pi p j}{q} \right) \log \sin\left( \frac{\pi j}{q} \right),
where γ\gamma is the Euler-Mascheroni constant.[12] This formula expresses ψ(p/q)\psi(p/q) in terms of elementary functions, the Euler-Mascheroni constant, and logarithms of sines, which are algebraic numbers related to the cyclotomic field Q(ζq)\mathbb{Q}(\zeta_q), where ζq=e2πi/q\zeta_q = e^{2\pi i / q}. The theorem holds for both odd and even qq, though the sum simplifies naturally for odd qq without a middle term contributing nontrivially. The formula originates from the work of Carl Friedrich Gauss, who derived it in his 1813 memoir on the hypergeometric series, where he analyzed the infinite product representations and limiting behaviors connected to the gamma function. Gauss's result was a significant advancement in understanding special values of the digamma function, building on Euler's earlier studies of the gamma function and harmonic numbers. Subsequent expositions, such as those by Jensen and Lehmer, refined and verified the expression using complex analysis and properties of the sine product.[12] A special case of the theorem arises when p=1p = 1, giving an explicit form for ψ(1/q)\psi(1/q):
ψ(1q)=γlog(2q)π2cot(πq)+2j=1q/2cos(2πjq)logsin(πjq). \psi\left( \frac{1}{q} \right) = -\gamma - \log(2q) - \frac{\pi}{2} \cot\left( \frac{\pi}{q} \right) + 2 \sum_{j=1}^{\lfloor q/2 \rfloor} \cos\left( \frac{2\pi j}{q} \right) \log \sin\left( \frac{\pi j}{q} \right).
This can alternatively be expressed using sums over quadratic residues modulo qq, leveraging the real parts of logarithms of roots of unity, though the trigonometric form is more direct for computation.[12] The derivation of Gauss's theorem relies on the multiplication theorem for the digamma function and the reflection formula. The multiplication theorem states that for positive integer mm,
ψ(mz)=logm+1mk=0m1ψ(z+km). \psi(m z) = \log m + \frac{1}{m} \sum_{k=0}^{m-1} \psi\left( z + \frac{k}{m} \right).
Setting z=1/mz = 1/m yields ψ(1)=γ=logm+1mk=0m1ψ(1+km)\psi(1) = -\gamma = \log m + \frac{1}{m} \sum_{k=0}^{m-1} \psi\left( \frac{1 + k}{m} \right), or equivalently,
k=1m1ψ(km)=mγmlogm. \sum_{k=1}^{m-1} \psi\left( \frac{k}{m} \right) = -m \gamma - m \log m.
The reflection formula ψ(1x)ψ(x)=πcot(πx)\psi(1 - x) - \psi(x) = \pi \cot(\pi x) then pairs terms: for k=1k = 1 to (m1)/2(m-1)/2 (assuming mm odd for simplicity),
ψ(km)+ψ(1km)=πcot(πkm). \psi\left( \frac{k}{m} \right) + \psi\left( 1 - \frac{k}{m} \right) = \pi \cot\left( \frac{\pi k}{m} \right).
Summing these pairs gives the cotangent contributions. To isolate individual ψ(p/m)\psi(p/m), one applies the Weierstrass infinite product for the gamma function, takes the logarithmic derivative, and uses Fourier expansion or limiting processes on the sine product sin(πz)=πzk=1(1z2/k2)\sin(\pi z) = \pi z \prod_{k=1}^\infty (1 - z^2/k^2), leading to the cosine-log-sine sum after separating real and imaginary parts in the complex logarithm.[13] This approach, originally due to Gauss, connects the rational evaluation to the periodicity and symmetry of the cotangent and sine functions.

Special Values and Sums

Special values

The digamma function admits closed-form expressions at positive integers, half-integers, and certain rational arguments, often involving the Euler-Mascheroni constant γ, natural logarithms, and trigonometric constants. These values arise from the functional equations and integral representations of the function.[14] At positive integers n ≥ 1, the digamma function is expressed using harmonic numbers as ψ(n) = -γ + H_{n-1}, where H_m = ∑_{k=1}^m 1/k is the mth harmonic number. Specific cases include ψ(1) = -γ, ψ(2) = -γ + 1, and ψ(3) = -γ + 3/2.[14] At half-integer arguments, the values are ψ(1/2) = -γ - 2 \ln 2 and ψ(3/2) = -γ + 2 - 2 \ln 2.[14] Closed-form expressions for the digamma function at other rational points follow from Gauss's digamma theorem, which provides a summation formula reducible to elementary terms for small denominators. For instance, at 1/3, ψ(1/3) = -γ - \frac{\pi}{2\sqrt{3}} - \frac{3}{2} \ln 3.[14] The digamma function has simple poles at non-positive integers z = 0, -1, -2, \dots, each with residue -1, so ψ(z) \to -\infty as z \to 0^+. Near a pole at z = -n for nonnegative integer n and small positive ε, the leading behavior is ψ(-n + ε) \approx -1/ε - γ + H_n.[14] Exact values at reciprocals of small integers 1/k (k = 1 to 6) are given in the following table, derived via the reflection formula or Gauss's digamma theorem:
kψ(1/k)
1
2-γ - 2 \ln 2
3-γ - \frac{3}{2} \ln 3 - \frac{\pi}{2 \sqrt{3}}
4-γ - 3 \ln 2 - \frac{\pi}{2}
5-γ - \frac{1}{2} \ln 5 - \frac{\pi}{2} \cot \frac{\pi}{5} + 2 \sum_{j=1}^{2} \cos \frac{2 \pi j}{5} \ln \sin \frac{\pi j}{5}
6-γ - 2 \ln 2 - \frac{1}{2} \ln 3 - \frac{\pi}{2} \cot \frac{\pi}{6} + 2 \sum_{j=1}^{2} \cos \frac{2 \pi j}{6} \ln \sin \frac{\pi j}{6}

Finite sums

The digamma function provides closed-form expressions for various finite sums related to harmonic numbers. Specifically, the (n1)(n-1)th harmonic number Hn1=k=1n11kH_{n-1} = \sum_{k=1}^{n-1} \frac{1}{k} is given by Hn1=ψ(n)+γH_{n-1} = \psi(n) + \gamma, where γ\gamma is the Euler-Mascheroni constant.[8] This relation follows from the definition ψ(n)=γ+Hn1\psi(n) = -\gamma + H_{n-1} for positive integers n1n \geq 1.[8] More generally, partial sums of the harmonic series can be expressed using differences of the digamma function. The sum k=mn1k=HnHm1=ψ(n+1)ψ(m)\sum_{k=m}^{n} \frac{1}{k} = H_n - H_{m-1} = \psi(n+1) - \psi(m), for positive integers 1mn1 \leq m \leq n.[8] This identity arises directly from the harmonic number representation and the recurrence property of the digamma function, ψ(z+1)=ψ(z)+1z\psi(z+1) = \psi(z) + \frac{1}{z}.[8] The sum of the first nn digamma values at positive integers also admits a closed form: k=1nψ(k)=nψ(n)n+1\sum_{k=1}^n \psi(k) = n \psi(n) - n + 1. This expression is derived from the relation ψ(k)=γ+Hk1\psi(k) = -\gamma + H_{k-1} and the known summation formula for harmonic numbers, k=1n1Hk=nHn1(n1)\sum_{k=1}^{n-1} H_k = n H_{n-1} - (n-1), yielding an exact result without asymptotic terms. A notable multidimensional-like finite sum arises from the multiplication theorem of the digamma function: k=0m1ψ(1+km)=mψ(m)mlnm\sum_{k=0}^{m-1} \psi\left(1 + \frac{k}{m}\right) = m \psi(m) - m \ln m, for positive integer m2m \geq 2. Excluding the k=0k=0 term, where ψ(1)=γ\psi(1) = -\gamma, gives k=1m1ψ(1+km)=mψ(m)mlnm+γ\sum_{k=1}^{m-1} \psi\left(1 + \frac{k}{m}\right) = m \psi(m) - m \ln m + \gamma. This sum evaluates rational-point arguments of the digamma function and connects to Gauss's digamma theorem for related rational evaluations.

Sums involving rational functions

The digamma function provides a closed-form expression for certain infinite sums arising from differences of harmonic-like terms, which often appear in the partial fraction decomposition of rational functions. Specifically, the sum k=1[1k+a1k+b]=ψ(b+1)ψ(a+1)\sum_{k=1}^{\infty} \left[ \frac{1}{k+a} - \frac{1}{k+b} \right] = \psi(b+1) - \psi(a+1) holds for (a)>1\Re(a) > -1, (b)>1\Re(b) > -1, and aba \neq b, where the convergence follows from the asymptotic behavior of the digamma function. This identity derives directly from the infinite series representation of the digamma function, ψ(z)=γ+k=1[1k1k+z1]\psi(z) = -\gamma + \sum_{k=1}^{\infty} \left[ \frac{1}{k} - \frac{1}{k+z-1} \right] for (z)>0\Re(z) > 0, by subtracting the series for ψ(a+1)\psi(a+1) and ψ(b+1)\psi(b+1). Rational functions of the form 1/(k(k+m))1/(k(k+m)), where mm is a positive integer, can be decomposed using partial fractions as 1k(k+m)=1m(1k1k+m)\frac{1}{k(k+m)} = \frac{1}{m} \left( \frac{1}{k} - \frac{1}{k+m} \right). The infinite sum then becomes k=11k(k+m)=1mk=1[1k1k+m]=ψ(m+1)+γm\sum_{k=1}^{\infty} \frac{1}{k(k+m)} = \frac{1}{m} \sum_{k=1}^{\infty} \left[ \frac{1}{k} - \frac{1}{k+m} \right] = \frac{\psi(m+1) + \gamma}{m}, since ψ(1)=γ\psi(1) = -\gamma. This evaluation leverages the aforementioned difference formula with a=0a = 0 and b=mb = m, highlighting the digamma function's role in summing first-order rational terms. More generally, partial fraction expansions of rational functions with simple poles lead to linear combinations of such digamma differences, enabling closed-form results for a broad class of series. For finite sums, telescoping series provide an initial approach that connects to the digamma function in the limit. For instance, k=1n1k(k+1)=k=1n(1k1k+1)=11n+1\sum_{k=1}^{n} \frac{1}{k(k+1)} = \sum_{k=1}^{n} \left( \frac{1}{k} - \frac{1}{k+1} \right) = 1 - \frac{1}{n+1}, which approaches 1 as nn \to \infty and aligns with the infinite sum formula using the recurrence relation ψ(z+1)=ψ(z)+1/z\psi(z+1) = \psi(z) + 1/z. Extending to higher-degree rationals, such as those decomposable into first-order differences, the partial sums can be expressed as ψ(n+a+1)ψ(a+1)[ψ(n+b+1)ψ(b+1)]\psi(n+a+1) - \psi(a+1) - [\psi(n+b+1) - \psi(b+1)], though the focus here remains on the infinite case where the digamma directly yields the exact value. While higher-order polygamma functions handle sums like k=1(k+z)2=ψ1(z)+ψ1(1z)\sum_{k=-\infty}^{\infty} \frac{1}{(k+z)^2} = \psi_1(z) + \psi_1(1 - z) (the trigamma function), the digamma specializes to first-order rational sums, emphasizing differences rather than derivatives. These representations are foundational in evaluating series in complex analysis and special functions, often appearing in integral transforms and residue computations.[15]

Asymptotic Behavior and Inequalities

Asymptotic expansion

The asymptotic expansion of the digamma function ψ(z)\psi(z) for large z|z| is derived from the Stirling series for the logarithm of the gamma function, LnΓ(z)\operatorname{Ln} \Gamma(z), by term-by-term differentiation, since ψ(z)=ddzLnΓ(z)\psi(z) = \frac{d}{dz} \operatorname{Ln} \Gamma(z).[16] As zz \to \infty in the sector argz<π|\arg z| < \pi, the leading terms of the expansion are given by
ψ(z)lnz12zk=1B2k2kz2k, \psi(z) \sim \ln z - \frac{1}{2z} - \sum_{k=1}^\infty \frac{B_{2k}}{2k z^{2k}},
where B2kB_{2k} denotes the 2k2k-th Bernoulli number. This series is divergent but provides a useful asymptotic approximation when truncated appropriately. The full Stirling series is
LnΓ(z)(z12)lnzz+12ln(2π)+m=1B2m2m(2m1)z2m1, \operatorname{Ln} \Gamma(z) \sim \left(z - \frac{1}{2}\right) \ln z - z + \frac{1}{2} \ln (2\pi) + \sum_{m=1}^\infty \frac{B_{2m}}{2m(2m-1) z^{2m-1}},
valid under the same sector condition. Differentiating this yields the expansion for ψ(z)\psi(z), confirming the form above, with the sum arising from the derivatives of the higher-order terms in the Stirling approximation.[16] The remainder after truncating the series at the term involving B2NB_{2N} is O(1/z2N+1)O(1/|z|^{2N+1}) as z|z| \to \infty. For positive real z>0z > 0, the remainder has the same sign as the first neglected term and is bounded by its absolute value; in the complex case, the bound involves a factor of sec2N+1(12argz)\sec^{2N+1}(\frac{1}{2} \arg z) times the first neglected term.[16]

Inequalities

The digamma function ψ(x)\psi(x) is strictly increasing on (0,)(0, \infty), since its first derivative, the trigamma function ψ(x)\psi'(x), satisfies ψ(x)>0\psi'(x) > 0 for all x>0x > 0. This monotonicity follows from the integral representation ψ(x)=0text1etdt>0\psi'(x) = \int_0^\infty \frac{t e^{-x t}}{1 - e^{-t}} \, dt > 0.[17] Additionally, ψ(x)\psi(x) is strictly concave on (0,)(0, \infty), as its second derivative, the polygamma function of order 2, is negative there. This concavity implies Jensen-type inequalities: for weights λi0\lambda_i \geq 0 with λi=1\sum \lambda_i = 1 and points xi>0x_i > 0, ψ(λixi)λiψ(xi)\psi\left( \sum \lambda_i x_i \right) \geq \sum \lambda_i \psi(x_i).[18] Several useful bounds exist for ψ(x)\psi(x). For x>1x > 1, the inequality ln(x1/2)<ψ(x)<lnx\ln(x - 1/2) < \psi(x) < \ln x holds, providing simple logarithmic enclosures. From the asymptotic expansion, truncating after the 1/x21/x^2 term yields ψ(x)=lnx12x112x2+θ(x)x2\psi(x) = \ln x - \frac{1}{2x} - \frac{1}{12 x^2} + \frac{\theta(x)}{x^2}, where θ(x)1120|\theta(x)| \leq \frac{1}{120} for x>0x > 0.[19] For differences, since the trigamma function satisfies ψ(z)>1/z2\psi'(z) > 1/z^2 for z>0z > 0, integration gives ψ(x)ψ(y)>yx1z2dz=xyxy\psi(x) - \psi(y) > \int_y^x \frac{1}{z^2} \, dz = \frac{x - y}{x y} whenever 0<y<x0 < y < x. This follows from the integral representation ψ(z)=0tezt1etdt>0teztdt=1z2\psi'(z) = \int_0^\infty \frac{t e^{-z t}}{1 - e^{-t}} \, dt > \int_0^\infty t e^{-z t} \, dt = \frac{1}{z^2}.[20]

Computation and Approximation

Numerical computation methods

Numerical computation of the digamma function ψ(z)\psi(z) typically begins by using the recurrence relation ψ(z+1)=ψ(z)+1/z\psi(z+1) = \psi(z) + 1/z to reduce the argument to the fundamental interval (0,1](0,1] or [1,2][1,2], where more efficient evaluation methods can be applied. For arguments in (0,1](0,1], the recurrence is applied forward to shift to [1,2][1,2], avoiding direct evaluation near the pole at z=0z=0. Once reduced, values are computed using series expansions, such as the Taylor series around z=1z=1, or integral representations like ψ(z)=0(ettezt1et)dt\psi(z) = \int_0^\infty \left( \frac{e^{-t}}{t} - \frac{e^{-zt}}{1-e^{-t}} \right) dt for z>0\Re z > 0, evaluated via numerical quadrature along paths of steepest descent to ensure convergence.[21] These techniques are particularly effective for moderate precision and are detailed in comprehensive surveys on special function computation.[22] For higher precision and efficiency, especially in arbitrary-precision arithmetic, Spouge's approximation for the gamma function is adapted to the digamma via logarithmic differentiation or finite differences. The core approximation for Γ(z)\Gamma(z) is Γ(z+1)2π(z+a)z+1/2e(z+a)[c0+k=1a1ck(z+a)k]\Gamma(z+1) \approx \sqrt{2\pi} (z+a)^{z+1/2} e^{-(z+a)} \left[ c_0 + \sum_{k=1}^{a-1} c_k (z+a)^{-k} \right], where aa controls the accuracy and ckc_k are precomputed coefficients, yielding relative errors bounded by a1/2(2π)(a+1/2)a^{-1/2} (2\pi)^{-(a+1/2)} for Re(z)0\operatorname{Re}(z) \geq 0. For ψ(z)\psi(z), this leads to ψ(z+1)ln(z+a)a1/2z+ak=1a1kck(z+a)(k+1)k=0a1ck(z+a)k\psi(z+1) \approx \ln(z+a) - \frac{a-1/2}{z+a} - \frac{ \sum_{k=1}^{a-1} k c_k (z+a)^{-(k+1)} }{ \sum_{k=0}^{a-1} c_k (z+a)^{-k} }, with absolute errors under D(ln2a)a1/2(2π)(a+1/2)/Re(z+a)D (\ln 2a) a^{1/2} (2\pi)^{-(a+1/2)} / \operatorname{Re}(z+a), where D1.021D \approx 1.021. This method offers O(1)O(1) complexity per evaluation after O(a)O(a) precomputation, making it suitable for large-scale computations.[23] Software libraries implement these and related algorithms for practical evaluation. In Mathematica, the function PolyGamma[0, z] (equivalent to ψ(z)\psi(z)) employs a combination of series expansions for small z|z|, asymptotic series for large z|z| in the right half-plane, and reflection formulas ψ(1z)ψ(z)=πcot(πz)\psi(1-z) - \psi(z) = \pi \cot(\pi z) for the left half-plane, with automatic handling of high precision.[24] Similarly, Python's SciPy library provides scipy.special.digamma(z), which uses the Cephes library's implementation relying on recurrences to the interval [1,2][1,2], followed by power series or continued fraction expansions for accuracy up to double precision, and asymptotic approximations for large z|z|.[25] These implementations ensure robust performance across real and complex arguments. Computing ψ(z)\psi(z) presents challenges due to its simple poles at non-positive integers z=0,1,2,z = 0, -1, -2, \dots, where residues are 1-1, requiring careful avoidance or special handling in algorithms. For complex arguments, the principal branch is defined with a branch cut along the negative real axis, leading to discontinuities that must be navigated, particularly near the poles; implementations often use analytic continuation via reflection or recurrence to stay in regions of analyticity. Numerical stability is further ensured by backward recurrence from asymptotic expansions for large positive real parts.[22]
Argument RangePrimary MethodKey Features
Small zz near (0,1]Power series (e.g., around z=1)Converges quickly; used after recurrence reduction; suitable for low to medium precision.
Large z|z|, Re(z)>0\operatorname{Re}(z) > 0Asymptotic expansion ψ(z)lnz1/(2z)k=1B2k/(2kz2k)\psi(z) \sim \ln z - 1/(2z) - \sum_{k=1}^\infty B_{2k}/(2k z^{2k})High accuracy for large arguments; combined with backward recurrence.
General complexContinued fraction or Spouge adaptationEfficient for arbitrary precision; handles wide ranges with bounded error; alternative to series for avoiding slow convergence.[23][26]

Approximation formulas

One common approximation for the digamma function arises from truncating the asymptotic expansion derived from Stirling's series for the logarithm of the gamma function. Specifically, for large positive xx, the digamma function satisfies
ψ(x+1)lnx12x112x2+1120x4, \psi(x+1) \approx \ln x - \frac{1}{2x} - \frac{1}{12x^2} + \frac{1}{120x^4} - \cdots,
where the series terms involve Bernoulli numbers and converge asymptotically as xx \to \infty.[23] This truncated form provides a practical approximation with error decreasing rapidly for x>10x > 10, achieving relative errors below 10610^{-6} in typical applications.[23] Rational approximations, particularly Padé-type approximants, offer improved accuracy for the digamma function near integer values or in regions where the asymptotic series converges slowly. A notable rational approximant for ψ(x)\psi(x) is constructed by accelerating the series using sequence transformations, yielding a form that matches the first few terms of the asymptotic expansion while bounding the transformation error explicitly.[27] For instance, such approximants achieve relative errors less than 101010^{-10} for xx in the interval [1,10][1, 10], making them suitable for numerical evaluations in that range.[28] A recent bounded-error approximation refines the Stirling-type formula by isolating the remainder term. For s>0s > 0,
ψ(s)=lns12s112s2+θ(s)s2, \psi(s) = \ln s - \frac{1}{2s} - \frac{1}{12 s^2} + \frac{\theta(s)}{s^2},
where the remainder satisfies 136<θ(s)<15\frac{1}{36} < \theta(s) < \frac{1}{5}.[19] This provides explicit error bounds without requiring further series truncation, with the maximum relative error bounded by approximately 0.028/s20.028/s^2 for s>1s > 1.[29] Integral-based quadrature approximations, such as those adapting Windschitl-type expansions for the gamma function to its logarithmic derivative, further extend these methods for moderate ss. These yield relative errors under 10810^{-8} for 1<s<1001 < s < 100 by integrating the asymptotic form over finite intervals.

Advanced Topics

Roots of the digamma function

The zeros of the digamma function ψ(z)\psi(z) lie entirely on the real axis in the complex plane and are all simple. There are no non-real zeros, and the function has exactly one positive zero with all others negative. This property follows from the fact that the digamma function is real-valued on the real line (away from its poles) and its behavior in the complex plane, governed by the reflection formula ψ(1z)=ψ(z)+πcot(πz)\psi(1 - z) = \psi(z) + \pi \cot(\pi z), ensures no off-axis zeros exist.[30][31] The positive zero occurs at z1.461632z \approx 1.461632. The negative zeros are located one in each open interval (n1,n)(-n-1, -n) for nonnegative integers n=0,1,2,n = 0, 1, 2, \dots, accumulating toward -\infty along the real axis. For large nn, the nnth negative zero αn\alpha_n satisfies the asymptotic αnn+1πarctan(πlogn)\alpha_n \approx -n + \frac{1}{\pi} \arctan\left(\frac{\pi}{\log n}\right).[30][32] The first few zeros can be approximated numerically as follows:
Zero indexApproximate location
Positive1.461632
n=0n=0-0.504083 (in (1,0)(-1,0))
n=1n=1-1.573498 (in (2,1)(-2,-1))
n=2n=2-2.610721 (in (3,2)(-3,-2))
n=3n=3-3.635293 (in (4,3)(-4,-3))
These values are computed via series expansions or numerical root-finding near the poles, with higher precision available through dedicated algorithms. The zeros are all irrational, with at most one possible rational exception (which does not occur).[33][34] The zeros of ψ(z)\psi(z) correspond to the points where the derivative of the gamma function vanishes, Γ(z)=0\Gamma'(z) = 0, since ψ(z)=Γ(z)/Γ(z)\psi(z) = \Gamma'(z)/\Gamma(z) and Γ(z)\Gamma(z) has no zeros. In terms of the Weierstrass infinite product representation of Γ(z)\Gamma(z), these occur where the logarithmic derivative of the product is zero, marking saddle points (critical points) of Γ(z)\Gamma(z) in the complex plane. Although the real zeros of ψ(z)\psi(z) exhibit spacing that accumulates asymptotically like the negative integers, their distribution has been analogized in some studies to the zeros of the Riemann zeta function due to similar accumulation patterns, but no direct connection to the Riemann hypothesis exists.[8][30]

Regularization

The digamma function provides the analytic continuation of the generalized harmonic numbers to complex arguments, defined for positive integer nn as Hn=k=1n1kH_n = \sum_{k=1}^n \frac{1}{k} and extended via Hz=ψ(z+1)+γH_z = \psi(z+1) + \gamma, where γ\gamma is the Euler-Mascheroni constant.[35] This extension assigns finite values to formally divergent harmonic-like sums through the meromorphic properties of ψ(z)\psi(z), enabling the regularization of expressions that diverge logarithmically.[15] A key application in regularization arises from the series representation of the digamma function, which interprets it as the Hadamard finite part of a divergent sum:
ψ(z)=γ+n=0(1n+11n+z), \psi(z) = -\gamma + \sum_{n=0}^\infty \left( \frac{1}{n+1} - \frac{1}{n+z} \right),
valid for (z)>0\Re(z) > 0. Here, the infinite sum n=1(1n1n+z1)\sum_{n=1}^\infty \left( \frac{1}{n} - \frac{1}{n+z-1} \right) diverges, but its finite part equals ψ(z)+γ\psi(z) + \gamma, extracting the regularized value by discarding the divergent terms associated with the Euler constant. This approach is particularly useful for regularizing differences of harmonic series in analytic continuations.[8] The digamma function also connects to zeta regularization through its Taylor series expansion around z=0z=0:
ψ(1+z)=γ+k=2(1)kζ(k)zk1,z<1, \psi(1 + z) = -\gamma + \sum_{k=2}^\infty (-1)^k \zeta(k) z^{k-1}, \quad |z| < 1,
where ζ(k)\zeta(k) is the Riemann zeta function evaluated at positive integers. This relation links the digamma's behavior near unity to the values of ζ(s)\zeta(s) for (s)>1\Re(s) > 1, facilitating the regularization of products or series involving the gamma function in contexts where zeta functions handle divergences at s=1s=1. For instance, in the analytic continuation of integrals like 0xz1ψ(x+1)dx=πsin(πz)ζ(1z)\int_0^\infty x^{z-1} \psi(x+1) \, dx = -\frac{\pi}{\sin(\pi z)} \zeta(1-z), the digamma regularizes the contribution from the zeta pole.[36] In quantum field theory, the digamma function appears in the regularization of spectral zeta functions used to compute Casimir energies. Spectral zeta functions ζ(s)=λns\zeta(s) = \sum \lambda_n^{-s} over eigenvalues λn\lambda_n of operators encounter divergences that are tamed via analytic continuation, often involving the gamma function whose logarithmic derivative is the digamma; for example, expressions for the effective action include terms like ψ(s)\psi(s) in the Mellin transform representations.[37] This role underscores the digamma's utility in extracting finite vacuum energies from formally infinite sums in curved spacetimes or bounded geometries.[38]

Applications in applied mathematics

In statistics, the digamma function plays a key role in modeling uncertainty and entropy for distributions involving compositional data, such as the Dirichlet and multinomial distributions. For instance, the expected logarithm of probabilities under a Dirichlet prior, which is essential for computing expected log-likelihoods in Bayesian models of categorical data, is given by the difference of digamma functions: E[logθi]=ψ(αi)ψ(jαj)\mathbb{E}[\log \theta_i] = \psi(\alpha_i) - \psi(\sum_j \alpha_j), where θi\theta_i are the components and α\alpha the concentration parameters.[39] This expression facilitates entropy calculations in topic models and sequential sampling, where the digamma function quantifies variations in information content, such as one-step entropy changes during species abundance estimation.[40] In Bayesian updating for multinomial parameters, ratios of digamma functions approximate posterior means of log-ratios, enabling efficient inference in hierarchical models without full integration.[41] In physics, particularly random matrix theory, the digamma function arises in the analysis of spectral densities and stability properties of ensembles. It contributes to the computation of Lyapunov exponents for products of random matrices, which describe the exponential growth rates of perturbations and underlie universal behaviors in local spectral statistics, such as level spacing distributions in non-Hermitian systems.[42] For black hole entropy, asymptotic approximations via Stirling's formula for large factorials in microstate counting often invoke the digamma function through its relation to the logarithmic derivative of the gamma function, providing corrections in high-dimensional partition sums that align with thermodynamic limits.[43] In signal processing, the digamma function supports estimation in time series models with compositional or positive-valued observations, such as Dirichlet autoregressive (AR) processes for multivariate signals. The expected values in these models, used for centering innovations in ARMA frameworks, are expressed in closed form using digamma functions, aiding phase and amplitude recovery in narrowband signals with gamma-modulated noise. This is particularly useful in autoregressive models where gamma-distributed variances model heteroscedasticity, with digamma terms appearing in the deviance for generalized linear AR structures.[44] Recent developments from 2020 to 2025 highlight the digamma function's role in machine learning, especially variational inference for approximate Bayesian computation. In optimizing the evidence lower bound (ELBO) for models with Beta or Dirichlet components, such as latent variable mixtures, digamma functions compute expectations of log-densities, enabling scalable inference in deep clustering and uncertainty quantification tasks.[45] For example, extended variational methods for Dirichlet process mixtures of Beta-Laplace distributions use digamma-based bounds to handle non-conjugate priors, improving posterior approximations in neural network evidential learning.[46] In hydrology, gamma distribution models for basin response and unit hydrograph derivation incorporate digamma functions to assess convergence rates in rainfall-runoff simulations, optimizing parameters for flood prediction via genetic algorithms on gamma-parameterized convergence.[47] Beyond these areas, the digamma function aids combinatorial analyses of partition functions by providing modular relations and asymptotic expansions for generalized forms, facilitating counts of integer partitions with constraints like distinct parts.[48] In numerical solutions to partial differential equations (PDEs), digamma functions serve as basis functions in mesh-free methods like the complex variable boundary element method (CVBEM), offering compact support and improved accuracy for boundary value problems without traditional meshing.[49]

References

User Avatar
No comments yet.