Hello,

sorry for the cryptic Subject, but I didn't know how else to title this.

I don't have much experience with numerical integration, so it may very

well be that I'm doing something stupid, please advise :-)

I need to compute some numerical integrals of some functions that are

very "localized" in their domain (trying to compute some posterior

probability density). A Gaussian (normal) function is a good example.

I find that if the region in which the function is significantly non

zero is away from the center of the integration interval, the

integration routines "step over" that domain interval and do not "see"

the function having values different from zero.

For example, integrating a Gaussian centered around zero, works as

expected:

integrate.quad(stats.norm(loc=0.0).pdf, -np.inf, np.inf)

(0.9999999999999998, 1.0178191320905743e-08)

but if I shift the Gaussian such that the function value at 0 is

numerically zero, the integration fails:

integrate.quad(stats.norm(loc=100.0).pdf, -np.inf, np.inf)

(0.0, 0.0)

same if I have an integration interval not centered:

integrate.quad(stats.norm(loc=0.0).pdf, -100, 1000)

(4.471689909323402e-30, 8.454234364421014e-30)

vs:

integrate.quad(stats.norm(loc=0.0).pdf, -100, 100)

(1.0000000000000002, 1.0346447361664605e-12)

I definitely see why those examples are problematic, but I don't know

what I can do about it other than evaluating the functions on a grid to

restrict somehow the integration limits.

Does anyone have any suggestion?

Thanks! Cheers,

Daniele

_______________________________________________

SciPy-User mailing list

[hidden email]
https://mail.python.org/mailman/listinfo/scipy-user