I want to approximate a function that I do not have an actual analytical expression for. I know that I want to compute this integral: integral a * b * c dx. Pretend that I get the a, b, and c are from observed data. How can I evaluate this integral? Can scipy do this? Is scipy.integrate.simps the right approach?
import numpy as np
from scipy.integrate import simps
a = np.random.random(10)
b = np.random.uniform(0, 10, 10)
c = np.random.normal(2, .8, 10)
x = np.linspace(0, 1, 10)
dx = x[1] - x[0]
print 'Is the integral of a * b * dx is ', simps(a * b, c, dx), ', ', simps(b * a, c, dx), ',', simps(a, b * c, dx), ', ', simps(a, c * b, dx), ', or something else?'
With your setup, the correct way to integrate is either
simps(a*b*c, x) # function values, argument values
or
simps(a*b*c, dx=dx) # function values, uniform spacing between x-values
Both yield the same result. Yes, simps is a very good choice for integrating sampled data. Most of the time it is more accurate than trapz.
If the data comes from a smooth function (even though you don't know the function), and you can somehow make the number of points to be 1 more than a power of 2, then Romberg integration will be even better. I compared trapz vs simps vs quad in this post.