Skip to content
Prev 279943 / 398506 Next

efficiently finding the integrals of a sequence of functions

JeffND <Zuofeng.Shang.5 <at> nd.edu> writes:
We had the same discussion last month under the heading "performance of
adaptIntegrate vs. integrate", see

   https://stat.ethz.ch/pipermail/r-help/2011-November/295260.html

It also depends on the accuracy you want to achieve. Adaptive methods will
be slower as the adaptation will be different in each single step for each
function, i.e. no vectorization here.

Non-adaptive Gaussian quadrature appears to be a good candidate. Assume you
have found grid points x_i and weights w_i for your interval [a, b], then if
F is the matrix with F_ij = f_j(x_i) amd the integrals will be computed all
at once with w %*% F .

Example: A function that returns x, x^2, ..., x^5 in columns

  f <- function(x) cbind(x, x^2, x^3, x^4, x^5)

The grid points and weights for the interval [0, 1] are:

  x <- c(0.025446, 0.129234, 0.297077, 0.500000, 0.702923, 0.870766, 0.974554)
  w <- c(0.064742, 0.139853, 0.190915, 0.208980, 0.190915, 0.139853, 0.064742)

and the integrals for these five functions are

    w %*% f(x)      # 0.5 0.3333334 0.25 0.2 0.1666667

Functions to calculate Gaussian points and weights are mentioned in the thread
above.

Hans Werner
sequence-of-functions-tp4179452p4179452.html