R: Integrate: Max number of subdivisions reached, roundoff error

Question

I came across an interesting but rather annoying problem. I am trying to integrate a function which has been calculated from a dataset. The data can be found here: Link to sample.txt. I start by fitting a line to my data. this can be done linear with `approxfun` or non-linear with `splinefun`. In my example below I use the latter. Now, when I try to integrate the fitted function I run into the error - `maximum number of subdivisions reached` but when I increase the subdivision I get - `roundoff error` From the values in my sample code, you can see that for this specific dataset the threshold is 754->755. My colleague has no problem to integrate this dataset in Matlab. Is there a way to manipulate my data to integrate? Is there another method for numerical integration in R? ``` data<-read.table('sample.txt',sep=',') colnames(data)<-c('wave','trans') plot(data$wave,data$trans,type='l') trans<- -1 * log(data$trans) plot(data$wave,trans,type='l') fx.spline<-splinefun(data$wave,trans) #Try…

Accepted Answer

There are many integration routines in R, and you can find some of them by 'RSiteSearch'ing or by using the 'sos' package.

For example, package `pracma` has several implementations, for instance

```
quad(fx.spline,min(data$wave),max(data$wave))   # adaptive Simpson
# [1] 2.170449                                  # 2.5 sec
quadgk(fx.spline,min(data$wave),max(data$wave)) # adaptive Gauss-Kronrod
# [1] 2.170449                                  # 0.9 sec
quadl(fx.spline,min(data$wave),max(data$wave))  # adaptive Lobatto
# [1] 2.170449                                  # 0.8 sec
```

Please not that these are pure R scripts and therefore slower than, e.g., the compiled `integrate` routine with such an oscillating function.

R: Integrate: Max number of subdivisions reached, roundoff error

Solution

Problem