R devel - Mar 2007 - integrate function (PR#9557)

Full_Name: Bert De Boeck
Version: R 2.2.0
OS: Windows
Submission from: (NULL) (157.193.193.152)


I think there is a bug when using integrate for integrating a function which is
0 in a whole sub-interval. For example:

#define uniform function
f<-function(x){ifelse(x<1,0,ifelse(x<3,1,0))}

#this is the correct integral
integrate(f,-10,10)

#here there is a problem
integrate(f,-50,50)
integrate(f,-10,50)
integrate(f,-50,10)
integrate(f,-50,50,sub=10000)

# I noticed this for a more complex function, but as you see even for a trivial

# function there is a serious problem

I think the behaviour is somehow documented, see ?integrate, in any case
I take the opportunity to report some more on integrate:


Under Win XP (AMD 64 3700+ 2Gb RAM)
> integrate(dnorm,0,2,sub=1e+09)
Error: cannot allocate vector of size 1889872 Kb
In addition: Warning messages:
1: Reached total allocation of 1536Mb: see help(memory.size)
2: Reached total allocation of 1536Mb: see help(memory.size)
> integrate(dnorm,0,2,sub=1e+10)
Error in if (limit < 1 || (abs.tol <= 0 && rel.tol < max(50 *
.Machine$double.eps,  :
        missing value where TRUE/FALSE needed
In addition: Warning message:
NAs introduced by coercion

but .....
> integrate(dnorm,0,2,sub=1e+08)
crashes the RGUI with an access violation error
> R.version
               _
platform       i386-pc-mingw32
arch           i386
os             mingw32
system         i386, mingw32
status
major          2
minor          4.1
year           2006
month          12
day            18
svn rev        40228
language       R
version.string R version 2.4.1 (2006-12-18)


*****************************************************************************************
Under LINUX SUSE 10.2 (Quad Opteron 8218 32Gb RAM)
> integrate(dnorm,0,2,sub=1e+08)
0.4772499 with absolute error < 5.3e-15
> integrate(dnorm,0,2,sub=1e+09)
 *** caught segfault ***
address (nil), cause 'memory not mapped'

Traceback:
 1: .External("call_dqags", ff, rho = environment(), as.double(lower),
    as.double(upper), as.double(abs.tol), as.double(rel.tol),
limit = limit, PACKAGE = "base")
 2: integrate(dnorm, 0, 2, sub = 1e+09)


> R.version
               _
platform       x86_64-unknown-linux-gnu
arch           x86_64
os             linux-gnu
system         x86_64, linux-gnu
status
major          2
minor          4.1
year           2006
month          12
day            18
svn rev        40228
language       R
version.string R version 2.4.1 (2006-12-18)

______________________________________________________

Simone Giannerini
Dipartimento di Scienze Statistiche "Paolo Fortunati"
Universita' di Bologna
Via delle belle arti 41 - 40126  Bologna,  ITALY
Tel: +39 051 2098262  Fax: +39 051 232153
______________________________________________________


On 3/8/07, Bert.DeBoeck at ugent.be <Bert.DeBoeck at ugent.be>
wrote:
> Full_Name: Bert De Boeck
> Version: R 2.2.0
> OS: Windows
> Submission from: (NULL) (157.193.193.152)
>
>
> I think there is a bug when using integrate for integrating a function
which is
> 0 in a whole sub-interval. For example:
>
> #define uniform function
> f<-function(x){ifelse(x<1,0,ifelse(x<3,1,0))}
>
> #this is the correct integral
> integrate(f,-10,10)
>
> #here there is a problem
> integrate(f,-50,50)
> integrate(f,-10,50)
> integrate(f,-50,10)
> integrate(f,-50,50,sub=10000)
>
> # I noticed this for a more complex function, but as you see even for a
trivial
>
> # function there is a serious problem
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>

-- 
______________________________________________________

Simone Giannerini
Dipartimento di Scienze Statistiche "Paolo Fortunati"
Universita' di Bologna
Via delle belle arti 41 - 40126  Bologna,  ITALY
Tel: +39 051 2098262  Fax: +39 051 232153

This is what numerical integration functions are like, unfortunately. 
That's why the help page says:
  "Like all numerical integration routines, these evaluate the function on
a
   finite set of points. If the function is approximately constant (in
   particular, zero) over nearly all its range it is possible that the
   result and error estimate may be seriously wrong."

If you know anything about where the function has jumps you can use this 
knowledge to pick better intervals.  Changing the subdvisions argument has 
no effect, since this is the *maximum* number, not the minimum number.

I don't think we know of numerical quadrature routines that are uniformly 
better than these and GPL-compatible.

 	-thomas


On Thu, 8 Mar 2007, Bert.DeBoeck at Ugent.be wrote:
> Full_Name: Bert De Boeck
> Version: R 2.2.0
> OS: Windows
> Submission from: (NULL) (157.193.193.152)
>
>
> I think there is a bug when using integrate for integrating a function
which is
> 0 in a whole sub-interval. For example:
>
> #define uniform function
> f<-function(x){ifelse(x<1,0,ifelse(x<3,1,0))}
>
> #this is the correct integral
> integrate(f,-10,10)
>
> #here there is a problem
> integrate(f,-50,50)
> integrate(f,-10,50)
> integrate(f,-50,10)
> integrate(f,-50,50,sub=10000)
>
> # I noticed this for a more complex function, but as you see even for a
trivial
>
> # function there is a serious problem
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
Thomas Lumley			Assoc. Professor, Biostatistics
tlumley at u.washington.edu	University of Washington, Seattle