r - R中正支持函数的卷积

Question

我想要在 [0,Inf) 上定义的两个函数的卷积， 说

f=function(x)
    (1+0.5*cos(2*pi*x))*(x>=0)

和

g=function(x)
    exp(-2*x)*(x>0)

使用 R I 的积分功能可以做到这一点，

cfg=function(x)
    integrate(function(y) f(y)*g(x-y),0,x)$value

通过在网上搜索，似乎有更有效(更准确)的方法来做到这一点(比如使用 fft() 或 convolve() )。有这种经历的人可以解释一下吗？
谢谢!

Answer 1

convolve或 fft解决方案是获得离散结果，而不是您在 cfg 中定义的函数.他们可以为您提供 cfg 的数字解决方案在一些常规的离散输入上。
fft用于周期函数(仅)，因此不会有帮助。然而，convolve有一种称为“open”的操作模式，它模拟由 cfg 正在执行的操作。 .

请注意，使用 type="open" ，您必须反转第二个序列(请参阅 ?convolve ，“详细信息”)。您还必须只使用结果的前半部分。这是c(2,3,5)的卷积结果的图例示例与 c(7,11,13)由 convolve(c(2,3,5), rev(c(7,11,13)), type='open') 执行:

             2  3  5       2  3  5   2  3  5    2  3  5      2  3  5
      13 11  7         13 11  7     13 11  7      13 11  7        13  11  7
 Sum:       14              43         94           94            65

请注意，前三个元素的评估类似于您的集成结果。最后三个将用于反向卷积。

这是与您的功能的比较。您的函数，矢量化，绘制

y <- seq(0,10,by=.01)
plot(y, Vectorize(cfg)(y), type='l')

以及convolve的申请用以下代码绘制。注意y中每单位间隔有100个点所以除以 100是合适的。

plot(y, convolve(f(y), rev(g(y)), type='open')[1:1001]/100, type='l')

这些不太同意，但卷积要快得多:

max(abs(Vectorize(cfg)(y) - convolve(f(y), rev(g(y)), type='open')[1:1001]/100))
## [1] 0.007474999

benchmark(Vectorize(cfg)(y), convolve(f(y), rev(g(y)), type='open')[1:1001]/100, columns=c('test', 'elapsed', 'relative'))
##                                                   test elapsed relative
## 2 convolve(f(y), rev(g(y)), type = "open")[1:1001]/100   0.056        1
## 1                                    Vectorize(cfg)(y)   5.824      104