【问题标题】:The "integrate" function in R gives erroneous resultsR 中的“积分”函数给出了错误的结果
【发布时间】:2016-01-23 18:15:18
【问题描述】:

我想在 x 上集成函数 M2_11(如下),用于固定 theta = c(2,0.8)c = 1.1a=c(1,1)A = matrix(c(1/0.8,0.03,0.03,2/0.8),nrow=2,ncol=2)

M2_11 = function(x, theta, c, a, A){
return((score1(x,theta)-a[1])^2* (weight(x, theta, c, a, A))^2 * f(x, theta))
}

R 的积分函数给出以下结果

theta = c(2,0.8)
c = 1.1
a=c(1,1)
A = matrix(c(1/0.8,0.03,0.03,2/0.8),nrow=2,ncol=2)
integrate(M2_11, lower = 1e-100, upper = 10 ,subdivisions = 10000, theta,c,a,A)

0.0006459957 绝对误差

以另一种方式进行集成得到相同的结果

fM2_11 = function(x){M2_11(x,theta,c,a,A)}
integrate(fM2_11, lower = 1e-100, upper = 10,subdivisions = 10000)

0.0006459957 绝对误差

然而,积分函数给出的结果显然是错误的:

x = seq(1e-100,10,by=0.001)
integrand = sapply(x,fM2_11)

曲线下面积明显大于0.00066

我还使用循环检查结果

loop_result = rep(NA,length(x))
for (i in 1:length(x)){
  loop_result[i] = M2_11(x[i],theta,c,a,A)  
}
table(integrand==loop_result)

是的

10001

发生了什么事?

【问题讨论】:

  • 为什么是integrand = sapply(x,fM2_11) 而不仅仅是fM2_11(x)?我没有深入研究代码,但我猜fM2_11 可能没有矢量化,而integrate 要求它是。阅读?integrate
  • 谢谢尼古拉!我欠你一个!

标签: r integrate


【解决方案1】:

非常感谢尼古拉。问题解决了。

integrate(Vectorize(fM2_11), lower = 1e-100, upper = 10 ,subdivisions = 10000) 

0.1588699,绝对误差

sum(integrand)*0.001 

0.1588705

没想到答案这么简单!

【讨论】:

  • 仅供参考,Vectorize() 几乎只是mapply() 的包装。在单个参数的情况下,单个结果函数将等效于sapply()
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 2016-06-07
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
相关资源
最近更新 更多