Holonomic Gradient Descent and its Application to Fisher-Bingham Integral

We give a new algorithm to find local maximum and minimum of a holonomic function and apply it for the Fisher-Bingham integral on the sphere $S^n$, which is used in the directional statistics. The method utilizes the theory and algorithms of holonomic systems.

💡 Research Summary

The paper introduces a novel algorithm called Holonomic Gradient Descent (HGD) for locating local maxima and minima of holonomic functions, and demonstrates its effectiveness on the Fisher‑Bingham integral defined on the sphere (S^n). A holonomic function is one that satisfies a finite set of linear differential operators; equivalently, it belongs to a D‑module that can be described by a Gröbner basis. Because such functions obey a Pfaffian system—a first‑order linear system of partial differential equations—their gradients and higher‑order derivatives can be expressed exactly in terms of a small number of auxiliary variables.

HGD exploits this property by converting the gradient and Hessian of the target function into an ordinary differential equation (ODE) system. Starting from an initial parameter vector (\theta_0) and the corresponding function value and derivative vector, the algorithm integrates the Pfaffian ODE while simultaneously updating (\theta) using a standard gradient‑descent step. The ODE integration is performed with a high‑order numerical solver (e.g., Runge–Kutta 4), and because the underlying system is linear and holonomic, the solution exists uniquely and remains numerically stable even in relatively high dimensions. No separate finite‑difference approximations are required, which eliminates the usual accumulation of truncation errors in gradient estimation.

The authors apply HGD to the normalizing constant of the Fisher‑Bingham distribution, \