[Rd] minor flaw in integrate()

Duncan Murdoch murdoch at stats.uwo.ca
Tue Jul 3 03:56:23 CEST 2007


On 28/06/2007 5:05 PM, Peter Ruckdeschel wrote:
> Hi,
> 
> I noticed a minor flaw in integrate() from package stats:
> 
> Taking up arguments lower and upper from integrate(),
> 
>    if (lower ==  Inf) && (upper ==  Inf)
> 
>        or
> 
>    if (lower == -Inf) && (upper == -Inf)
> 
> integrate() calculates the value for (lower==-Inf) && (upper==Inf).
> 
> Rather, it should return 0.

Wouldn't it be better to return NA or NaN, for the same reason Inf/Inf 
doesn't return 1?

Duncan Murdoch

> 
> Quick fix:
> 
> ### old code ###
> ### [snip]
>     else {
>             if (is.na(lower) || is.na(upper))
>                 stop("a limit is missing")
>             if (is.finite(lower)) {
>                 inf <- 1
>                 bound <- lower
>             }
>             else if (is.finite(upper)) {
>                 inf <- -1
>                 bound <- upper
>             }
>             else {
>                 inf <- 2
>                 bound <- 0
>             }
>             wk <- .External("call_dqagi", ff, rho = environment(),
>                 as.double(bound), as.integer(inf), as.double(abs.tol),
>                 as.double(rel.tol), limit = limit, PACKAGE = "base")
>     }
> ### [snip]
> 
> ### new code  to replace the old one ###
> 
> ### [snip]
>     else {
>             if (is.na(lower) || is.na(upper))
>                 stop("a limit is missing")
> 
>             if (lower == upper){
> 
>                 wk <- list("value" = 0, "abs.error" = 0,
>                         "subdivisions" = subdivisions,
>                         "ierr" = 0 )
> 
>             } else {
>                     if (is.finite(lower)) {
>                         inf <- 1
>                         bound <- lower
>                     }
>                     else if (is.finite(upper)) {
>                             inf <- -1
>                         bound <- upper
>                     }
>                     else {
>                         inf <- 2
>                         bound <- 0
>                     }
>                     wk <- .External("call_dqagi", ff, rho = environment(),
>                             as.double(bound), as.integer(inf),
>                             as.double(abs.tol), as.double(rel.tol),
>                             limit = limit, PACKAGE = "base")
> 
>             }
>     }
> ### [snip]
> 
> Best, Peter
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel



More information about the R-devel mailing list