8

I want to calculate/evaluate the convolution

$$g(x)=\int_D f(x-t) \phi(t) dt,$$

where $f$ is a density and $\phi$ is a smooth function with compact support $D$. The convolution is not available in closed-form and I need to integrate it numerically. My question is: Is there an efficient way to do this? I want to implement it in R, so, I would like to see if there is a better way than using the command integrate().

Community
  • 1
Cook
  • 83
  • 1
  • 4
  • 2
    Depending on circumstances I usually either discretize to a large power of 2 bins, and use the fast fourier transform (`?fft`) or use `convolve`. The fft approach takes a bit more work to set up but is better if you need to convolve with something several times. Sometimes it takes a while to figure out the right argument settings with convolve. – Glen_b Nov 05 '14 at 13:52
  • @Glen_b Thanks. For a univariate function $f$, I think the direct integration may be faster, then. – Cook Nov 05 '14 at 13:53
  • 2
    You asked for an efficient way - fft is *really fast*; it just requires a little bit of setup (binning, padding with zeros). – Glen_b Nov 05 '14 at 13:55
  • @Glen_b Yeah, I agree the fft is really fast, but the previous step may slow down the process. I will compare both methods, anyway. Thank you. – Cook Nov 05 '14 at 13:56
  • I recall using `convolve` for this purpose several times. Simple working examples appear at http://stats.stackexchange.com/a/41263, http://stats.stackexchange.com/a/41255, and http://stats.stackexchange.com/a/49444. – whuber Nov 05 '14 at 15:47

1 Answers1

11

Have you taken a look at dedicated R packages for that? Like convolve, https://stat.ethz.ch/R-manual/R-devel/library/stats/html/convolve.html

jmnavarro
  • 446
  • 5
  • 12
  • Thanks. If I properly understand, this package is for convolving numeric sequences, rather than two functions. Am I missing something? – Cook Nov 05 '14 at 13:51
  • 2
    Well, I've not used that package, but if you're trying to do the convolution in R, then you'll basically have two numeric sequences which will be the values from your functions, wouldn't they? I don't think R lets you work with the function definition. – jmnavarro Nov 05 '14 at 13:53