7

我遇到了一个有趣但相当烦人的问题。

我正在尝试集成一个从数据集中计算出来的函数。数据可以在这里找到:链接到 sample.txt

我首先为我的数据拟合一条线。这可以与 线性approxfun或非线性完成splinefun。在下面的示例中,我使用后者。现在,当我尝试集成拟合函数时,我遇到了错误

  • maximum number of subdivisions reached

但是当我增加细分时,我得到

  • roundoff error

从我的示例代码中的值可以看出,对于这个特定的数据集,阈值是 754->755。

我的同事将这个数据集集成到 Matlab 中没有问题。有没有办法操纵我的数据进行整合?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 either
Fx.spline<-integrate(fx.spline,min(data$wave),max(data$wave))
#Above: Number of subdivision reached
Fx.spline<-integrate(fx.spline,min(data$wave),max(data$wave),subdivisions=754)
#Above: Number of subdivision reached
Fx.spline<-integrate(fx.spline,min(data$wave),max(data$wave),subdivisions=755)
#Above: Roundoff error
4

1 回答 1

7

R 中有许多集成例程,您可以通过“RSiteSearch”或使用“sos”包找到其中的一些。

例如,包pracma有几个实现,例如

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

请注意,这些是纯 R 脚本,因此比integrate具有这种振荡功能的编译例程慢。

于 2012-05-04T19:04:16.210 回答