Custom-made Gauss quadrature: an introduction for statisticians

An $n$-point Gauss quadrature rule approximates the weighted integral of a function by a weighted average of $n$ evaluations of this function and is exact for polynomials of degree at most $2n-1$. Such rules can be highly accurate with relatively few evaluations. For weight functions that are associated with classical orthogonal polynomials of a continuous variable (such as Legendre, Hermite and Laguerre), these rules are readily available. We suppose that this is not the case, so that these rules must be custom-made. The two most easily understood methods for the computation of these rules are (a) moment determinants and (b) the Stieltjes procedure. We implement them in the Julia package CustomGaussQuadrature, which uses type-generic numerical programming and adaptive high-precision arithmetic to assess the approximation error due to roundoff. We describe access from R via JuliaConnectoR.