Integrate: the integral is probably divergent

问题:

I was doing some integration into a loop using integrate and I come up with an error I can't understand neither get rid of. Here is a MWE I could extract:

u_min = 0.06911363
u_max = 1.011011 
m = 0.06990648
s = 0.001092265
integrate(f = function(v){pnorm(v, mean = m, sd = s, lower.tail =  FALSE)}, u_min, u_max)

this returns an error "the integrale is probably divergent" which is obviously false. I tried to modify the parameters a little bit and got this working for example:

u_min <- 0.07
u_max <- 1.1
m <- 0.0699
s <- 0.00109
integrate(f = function(v){pnorm(v, mean = m, sd = s, lower.tail =  FALSE)}, u_min, u_max)

I tried to have a look into the integrate function with debugbut it's a wrapper of a Ccode. Also I'm not a specialist of quadrature techniques. I saw this SO post but couldn't make anything from it.

thanks

回答1:

The default tolerance of .Machine$double.eps^0.25 (= 0.0001220703) needs to be lowered. Try, for example, this:

f <- function(v) pnorm(v, mean = m, sd = s, lower.tail =  FALSE)
integrate(f, u_min, u_max, rel.tol = 1e-15)

## 0.0009421867 with absolute error < 1.1e-17

回答2:

I'd use this work-around:

integrate(f = function(v){pnorm(v, mean = m, sd = s, lower.tail =  FALSE)}, 
      max(u_min,m-10*s),min(u_max,m+10*s))$value  + (u_min-m+10*s)*(u_min<m+10*s)

What I've done:

• pnorm with lower.tail=FALSE is basically zero when very far on the right from the mean. So there is no point in "stretching" the right limit of the integral. So, when u_max > m+10*s, you just integrate to m + 10*s. You can of course change the 10 factor to add precision;
• on the other hand, on the left pnorm is basically always 1; so you can enhance the left limit and the missing part is just u_min - m+10*s. Same logic as above.