---
title: 'Bessel Distribution: Theory & Applications'
url: https://www.emergentmind.com/topics/bessel-distribution
type: topic
---

# Bessel Distribution: Theory & Applications

The term **Bessel distribution** denotes several distinct constructions whose common feature is the appearance of modified Bessel functions in a density, transform, or kernel formula. In probability and statistics, two prominent uses are a zero-order modified Bessel, or **Bessel-\(K\)**, law on \(\mathbb R\) obtained by compounding a zero-mean Gaussian with a \(\chi^2_1\)-distributed variance, and a univariate normalized inverse-Gaussian law on \((0,1)\) used for regression with bounded responses [2507.21812], [2003.05157]. Related work studies broader modified Bessel families on \([0,\infty)\), including McKay, \(K\)-distribution, and generalized inverse-Gaussian models, with emphasis on infinite divisibility and allied structural properties [2406.17721]. In p-adic representation theory, by contrast, a **Bessel distribution** is a bi-equivariant distribution attached to a generic representation of \(GL_n(F)\), rather than a probability law [1512.02365].

## 1. Terminological scope and mathematical setting

In the recent literature, the phrase **Bessel distribution** is not reserved for a single canonical probability law. The construction studied in "Properties and approximations of a Bessel distribution for data science applications" is a real-valued law generated from Gaussian variance mixing and expressed through the zero-order modified Bessel function \(K_0\); that paper also places the law as a special case of the McKay family of Bessel distributions and of a family of generalized Laplace distributions [2507.21812]. The construction studied in "Bessel regression model: Robustness to analyze bounded data" is instead the univariate normalized inverse-Gaussian distribution on \((0,1)\), parameterized either by \((\alpha,\beta)\) or by mean and precision \((\mu,\phi)\) [2003.05157].

The broader survey literature uses **modified Bessel distributions** for continuous laws on \([0,\infty)\) whose densities or Laplace transforms involve modified Bessel functions of the first kind \(I_\nu\) or the second kind \(K_\nu\). This includes McKay distributions, the \(K\)-distribution, generalized inverse-Gaussian distributions, and several new families obtained from Stieltjes-transform representations [2406.17721]. A separate representation-theoretic literature uses the same expression for a distribution on \(GL_n(F)\) attached to Whittaker data; in that setting, the object is a distribution in the sense of harmonic analysis, not a random variable [1512.02365].

This multiplicity of usage is important because formulas, support, and inferential roles differ sharply across the constructions. The probabilistic families are linked by Bessel-function structure, but they are not reparameterizations of one another.

## 2. Zero-order modified Bessel law on \(\mathbb R\)

For the real-valued construction in [2507.21812], let \(X\sim N(0,\sigma^2)\), and conditionally on \(X=x\), let \(Y\mid X=x\sim N(0,x^2)\). The marginal law of \(Y\) is the zero-order modified Bessel distribution with scale parameter \(\sigma\), denoted \(Y\sim K(\sigma)\). An equivalent representation is
\[
Y=\sigma Z_1 Z_2,
\]
where \(Z_1\) and \(Z_2\) are independent standard normals [2507.21812].

Its density admits the standard form
\[
f_Y(y)=\frac{1}{\pi \sigma}K_0(|y|/\sigma), \qquad y\in\mathbb R,
\]
where \(K_0(\cdot)\) is the modified Bessel function of the second kind of order zero. The law is symmetric, so \(E[Y]=0\), and its variance is
\[
\operatorname{Var}[Y]=\sigma^2.
\]
The even moments are
\[
E[Y^{2m}]=(\sigma^2)^m\frac{(2m)!}{2^m m!}, \qquad E[Y^{2m+1}]=0.
\]

A closed form is also available for the cumulative distribution function. With \(L_\nu\) denoting the modified Struve functions,
\[
F_Y(y)=\frac12+\frac{y}{2\sigma}\left[K_0(|y|/\sigma)L_{-1}(|y|/\sigma)+K_1(|y|/\sigma)L_0(|y|/\sigma)\right], \qquad y\in\mathbb R.
\]
The stated motivation for this representation is computational: it avoids numerical integration of \(K_0\) and exploits fast routines for \(K_\cdot\) and \(L_\cdot\) [2507.21812].

The construction is analytically notable because it sits at the intersection of Gaussian scale mixtures, generalized Laplace families, and Bessel-function models. The product representation \(Y=\sigma Z_1Z_2\) also makes explicit that the law is generated by multiplicative Gaussian structure rather than by additive exponential tails alone.

## 3. Approximations, computation, and hypothesis testing

A central theme in [2507.21812] is that the Bessel-\(K\) law can often be replaced by simpler approximations when quantiles or \(p\)-values are the main computational target. The first approximation is Laplace. If the \(\chi^2_1\) variance in the Gaussian-mixture representation is replaced by a \(\chi^2_2\) variable, the resulting law has density
\[
f_L(x;b)=\frac{1}{2b}\exp(-|x|/b), \qquad x\in\mathbb R.
\]
Matching variances gives \(b=\sigma/\sqrt{2}\), since \(\operatorname{Var}[Y]=\sigma^2\) and \(\operatorname{Var}[L]=2b^2\). Optimizing a Kolmogorov-Smirnov distance over \(b\) gives \(b\approx \sigma/1.83\), while minimizing a Wasserstein distance gives \(b\approx \sigma/1.54\). For hypothesis-testing applications, especially \(90\%\), \(95\%\), and \(99\%\) one-sided critical values, taking \(b\approx \sigma/1.5\) yields quantiles within \(10\%\) of the exact Bessel values [2507.21812].

A second approximation is the empirical power-series or Martin-Maas approximation. For \(x\ge 0\),
\[
K_0(x)\simeq \sqrt{\pi/(2x)}\,e^{-x}.
\]
Replacing \(K_0(|y|/s)\) by this leading term and renormalizing gives
\[
f_M(y)=\frac{1}{2\sqrt{\pi s}}\frac{e^{-|y|/s}}{\sqrt{|y|}}, \qquad y\in\mathbb R,
\]
with elementary-form CDF
\[
F_M(y)=\frac12\left[1+\operatorname{erf}\!\bigl(\sqrt{y/s}\bigr)\right], \qquad y\ge 0,
\]
and the symmetric counterpart for \(y<0\). The approximation is stated to be most accurate in the tails, \(|y|\gg s\), while retaining the \(1/\sqrt{|y|}\) singularity at the origin [2507.21812].

The computational motivation is explicit. Direct evaluation of \(K_0\) and Struve-\(L\) functions is described as typically tens of times slower than elementary exponentials and square roots. By contrast, the Laplace density and CDF are \(O(1)\)-cost elementary operations, and the Martin-Maas CDF requires one error-function call but avoids Struve-\(L\). The paper gives a p-value sketch for an observed \(Y_{\text{obs}}\): choose the Laplace or Martin-Maas approximation, compute \(u=|Y_{\text{obs}}|/b\) or \(u=|Y_{\text{obs}}|/s\), take \(p=\exp(-u)\) for Laplace or \(p=1-\operatorname{erf}(\sqrt{u})\) for Martin-Maas, and return one-sided \(p/2\) [2507.21812].

These approximations are motivated by applications in statistical hypothesis testing. Exact critical values for \(K(\sigma)\) require inversion of the Struve-\(L\) representation. In regression of Poisson data with systematic errors, one often convolves \(\chi^2_\nu\) with \(K(\sigma)\); replacing \(K(\sigma)\) by \(\operatorname{Laplace}(0,\sigma/1.5)\) yields closed-form likelihood-ratio distributions and rapid p-value calculation. For large-scale inference, including imaging and genomics, the reported gain is that replacing special-function calls by \(O(1)\) exponentials can reduce computation time by orders of magnitude with negligible loss of accuracy [2507.21812].

## 4. Univariate normalized inverse-Gaussian Bessel distribution and regression

In the regression literature, the **Bessel distribution** is a bounded-response law derived from inverse-Gaussian variables rather than a symmetric real-valued law. Let \(Y_1\) and \(Y_2\) be independent inverse-Gaussian random variables with common scale parameter \(1\) and shape parameters \(\alpha>0\) and \(\beta>0\). Then
\[
Z=\frac{Y_1}{Y_1+Y_2}
\]
has the univariate normalized inverse-Gaussian distribution, written \(Z\sim \mathrm{Bessel}(\alpha,\beta)\), with support \(0<Z<1\) [2003.05157].

A mean-precision reparameterization is
\[
\mu=\frac{\alpha}{\alpha+\beta}\in(0,1), \qquad \phi=\alpha+\beta>0,
\]
equivalently \(\alpha=\mu\phi\) and \(\beta=(1-\mu)\phi\), so that \(Z\sim \mathrm{Bessel}(\mu,\phi)\). In this parameterization, \(\mu\) is the mean and \(\phi\) is a precision parameter. The density is written through
\[
\zeta_\mu(z)=\sqrt{1+\frac{(z-\mu)^2}{z(1-z)}},
\]
and involves the modified Bessel function \(K_1\bigl(\phi\,\zeta_\mu(z)\bigr)\) together with the factor \([z(1-z)]^{-3/2}\) [2003.05157].

The first two moments are
\[
E(Z)=\mu,
\]
and
\[
\operatorname{Var}(Z)=\mu(1-\mu)\,
\frac{1-\phi+\phi^2 e^\phi Ei(\phi)}{2},
\]
where
\[
Ei(\phi)=\int_1^\infty u^{-1}e^{-\phi u}\,du
\]
is the exponential integral. This variance function underlies the paper’s comparison with beta regression [2003.05157].

For regression, the proposed specification is
\[
\operatorname{logit}(\mu_i)=\mathbf x_i^\top\kappa, \qquad \ln(\phi_i)=\mathbf v_i^\top\lambda.
\]
The estimation strategy is EM. With the augmentation \(W_i=Y_{1i}+Y_{2i}\), one has
\[
W_i\mid (Z_i=z_i)\sim \mathrm{GIG}\bigl(a=1,\; b=\phi_i^2\zeta_{\mu_i}(z_i)^2,\; s=-1\bigr),
\]
which yields conditional expectations
\[
\psi_i=E[W_i^{-1}\mid Z_i=z_i], \qquad \chi_i=E[W_i^{-2}\mid Z_i=z_i].
\]
The E-step computes \(\psi_i^{(r)}\) and \(\chi_i^{(r)}\) at the current iterate. The M-step maximizes the resulting \(Q\)-function numerically in \(\theta=(\kappa,\lambda)\); closed-form updates are not available, so BFGS or Newton-Raphson is used with analytic gradients. Iteration continues until
\[
\|\theta^{(r+1)}-\theta^{(r)}\|/\|\theta^{(r)}\|<\varepsilon
\]
[2003.05157].

Under standard regularity, the EM estimator is the MLE and is \(\sqrt n\)-consistent and asymptotically normal. Observed information is obtained through Louis’ formula,
\[
I_n(\theta)=E\!\left[-\frac{\partial^2\ell_c}{\partial\theta\,\partial\theta^\top}\Bigm|\mathbf Z\right]
-E\!\left[\frac{\partial\ell_c}{\partial\theta}\frac{\partial\ell_c}{\partial\theta}^\top\Bigm|\mathbf Z\right],
\]
which is then evaluated at \(\widehat\theta\) to obtain standard errors and asymptotic confidence intervals [2003.05157].

The same paper introduces the DBB criterion, a discrimination procedure between Bessel and beta regressions based on the variance functions
\[
g_{\mathrm{beta}}(\phi)=\frac{1}{1+\phi}, \qquad
g_{\mathrm{bessel}}(\phi)=\frac{1-\phi+\phi^2 e^\phi Ei(\phi)}{2}.
\]
Simulation results reported there indicate that when data are generated from a Bessel regression, the EM estimators show small bias and correct coverage even for moderate \(n\). Under contamination of beta-regression data by a small fraction of outliers, Bessel regression shows smaller bias in the \(\kappa\)-parameters than beta regression, especially for larger \(n\) and larger contamination fractions. In three empirical illustrations, the DBB test often selects the Bessel model; in those cases, the Bessel regression shows better residual-diagnostic envelopes, smaller cross-validated residual sums of squares, and smaller first-and-second-moment distance on held-out data [2003.05157].

## 5. Broader modified Bessel families on \([0,\infty)\)

The survey "Infinitely divisible modified Bessel distributions" places Bessel-type laws in a larger probabilistic class: continuous univariate distributions on \([0,\infty)\) whose densities or Laplace transforms involve \(I_\nu\) or \(K_\nu\) [2406.17721]. The families discussed include McKay type I laws, generalized McKay laws, a squared-\(I_\mu\) law, the \(K\)-distribution, and the generalized inverse-Gaussian distribution.

One classical example is the \(K\)-distribution. If
\[
X\sim \Gamma(a,1), \qquad Y\sim \Gamma(\beta,1), \qquad Z=XY,
\]
then the density of \(Z\) is
\[
w_{a,\beta}(x)=\frac{2}{\Gamma(a)\Gamma(\beta)}\,x^{\tfrac{a+\beta}{2}-1}\,
K_{a-\beta}\!\bigl(2\sqrt{x}\bigr), \qquad x>0.
\]
In radar and wireless theory this is called the \(K\)-distribution or gamma-gamma law. Another fundamental example is the generalized inverse-Gaussian law
\[
f_{\mu,a,b}(x)=\frac{(a/b)^{\mu/2}}{2K_\mu(\sqrt{ab})}\,
x^{\mu-1}\exp\!\left(-\frac{a\,x+b/x}{2}\right),
\qquad x>0,\quad a,b>0,\ \mu\in\mathbb R.
\]

The structural emphasis of [2406.17721] is on four classes: infinite divisibility (ID), self-decomposability (SD), generalized gamma convolutions (GGC), and hyperbolically completely monotone (HCM) densities. For the McKay-type family built from \(x^\mu e^{-bx}I_\mu(ax)\), the paper shows ID, SD, and GGC. For the \(K\)-distribution, it states that the density is HCM and that the law is in fact ID, SD, GGC, and HCM. For GIG, it shows membership in GGC and therefore also SD and ID [2406.17721].

A technical contribution of that work is a collection of Stieltjes-transform representations for products, quotients, and reciprocals involving modified Bessel functions, including formulas for \(I_\mu(bz)K_\nu(az)\), \(K_\mu(az)K_\nu(bz)\), reciprocal products, and ratios such as \(K_\nu(az)/K_\nu(bz)\). These representations are used to construct new infinitely divisible laws with Laplace transforms of the form
\[
L_X(t)=F(b+t)/F(b),
\]
where \(F\) is built from \(I_\nu\), \(K_\nu\), or their products [2406.17721].

The survey also records open problems, including determining the full parameter domain for which generalized McKay laws are ID, SD, GGC, or HCM; deciding whether the \(I^2\)-type law is HCM or GGC beyond the range \(\mu\in(-\tfrac12,0]\); characterizing noncentral \(\chi^2(\nu,\lambda)\) densities in HCM or GGC terms; extending results for Tricomi-\(\Psi\) quotients; and generalizing ratio-of-gammas arguments to matrix-variate and multivariate settings [2406.17721].

## 6. Bessel distribution in p-adic representation theory

In p-adic harmonic analysis, **Bessel distribution** has a different meaning. Let \(F\) be a non-archimedean local field, \(G_n=GL_n(F)\), \(N_n\subset G_n\) the standard upper-triangular unipotent subgroup, and \((\pi,V)\) an irreducible, admissible, generic representation with contragredient \((\tilde\pi,\tilde V)\). If
\[
\ell\in \operatorname{Hom}_{N_n}(\pi,\psi), \qquad
\ell'\in \operatorname{Hom}_{N_n}(\tilde\pi,\psi^{-1})
\]
are nonzero Whittaker functionals, then for \(f\in C_c^\infty(G_n)\),
\[
v_{f,\ell'}:=\int_{G_n} f(g)\,\tilde\pi(g^{-1})\,\ell'\,dg \in V,
\qquad
B_{\ell,\ell'}(f):=\ell(v_{f,\ell'})
\]
defines a bi-\(N_n,\psi-\psi\)-equivariant distribution on \(G_n\), called the Bessel distribution [1512.02365].

On the big Bruhat cell
\[
Y:=N_n w_n A_n N_n,
\]
this distribution is represented by a locally constant function \(j_0(g)\), in the sense that
\[
B_{\ell,\ell'}(f)=\int_Y f(g)\,j_0(g)\,dg
\]
for test functions supported in \(Y\). A second construction starts from the Whittaker model \(W(\pi,\psi)\): for \(W\in W(\pi,\psi)\) and \(g\in Y\),
\[
I(W;g):=\int_{N_n} W(gu)\,\psi(-u)\,du
\]
converges stably, defines a Whittaker functional, and therefore equals
\[
I(W;g)=j_\pi(g)\,W(1)
\]
for a scalar \(j_\pi(g)\), called the Whittaker-integral Bessel function [1512.02365].

The main theorem of [1512.02365] is a weak kernel formula for Bessel functions attached to irreducible generic representations of p-adic \(GL(n)\). For supercuspidal \(\pi\), the formula expresses \(W(bw_n)\) for \(b\in A_n\) as an iterated Fourier-Mellin integral. Chai then proves that the two Bessel functions coincide on \(Y\):
\[
j_0(g)=j_\pi(g).
\]
The paper states applications to local Bessel identities in the Waldspurger correspondence, to the comparison between Jacquet-integral and distribution-theoretic definitions in local Gross-Prasad settings, and to the study of local coefficients and Kirillov models [1512.02365].

This representation-theoretic usage shares the Bessel nomenclature because Bessel functions and Whittaker models occupy the analytic core of the construction, but the object is a distribution on a reductive \(p\)-adic group rather than a probability distribution. The overlap with the probabilistic literature is therefore terminological and special-function-theoretic, not measure-theoretic on a common sample space.

Source: https://www.emergentmind.com/topics/bessel-distribution