Computing fast and accurate convolutions

Abstract

The analysis of data often models random components as a sum of in- dependent random variables (RVs). These RVs are often assumed to be lattice-valued, either implied by the problem or for computational efficiency. Thus, such analysis typically requires computing, or, more commonly, ap- proximating a portion of the distribution of that sum. Computing the underlying distribution without approximations falls un- der the area of exact tests. These are becoming more popular with continuing increases in both computing power and the size of data sets. For the RVs above, exactly computing the underlying distribution is done via a convolu- tion of their probability mass functions, which reduces to convolving pairs of non-negative vectors. This is conceptually simple, but practical implementations must care- fully consider both speed and accuracy. Such implementations fall prey to the round-off error inherent to floating point arithmetic, risking large rela- tive errors in computed results. There are two main existing algorithms for computing convolutions of vectors: naive convolution (NC) has small bounds on the relative error of each element of the result but has quadratic runtime; while Fast Fourier Transform-based convolution (FFT-C) has almost linear runtime but does not control the relative error of each element, due to the accumulation of round-off error. This essay introduces two novel algorithms for these problems: aFFT-C for computing convolution of two non-negative vectors, and sisFFT for com- puting p-values of sums of independent and identically-distributed lattice- valued RVs. Through a rigorous analysis of round-off error and its accumula- tion, both aFFT-C and sisFFT provide control of the relative error similar to NC, but are typically closer in speed to FFT-C by careful use of FFT-based convolutions and by aggressively discarding irrelevant values. Both accuracy and performance are demonstrated empirically with a variety of examples.