Yoann MocquinFollowTowards Data Science--ListenShareThe Fourier transform algorithm is considered one of the greatest discoveries in all of mathematics. French mathematician Jean-Baptiste Joseph Fourier laid the foundation for harmonic analysis in his book “Théorie analytique de la chaleur” in 1822. Today, the Fourier transform and all its variants form the basis of our modern world, powering technologies like compression, communication, image processing.This wonderful framework also provides great tools for analysing time-series… and that’s why we’re here!This post is part of a series on the Fourier transform. Today we will talk about convolution and how the Fourier transform provides the fastest way you can do it.All figures and equations are made by the author.Let’s start with basic definitions. The discrete Fourier transform for a discrete time sequence x of N elements is :where k denotes the k-th frequency of the spectrum of x. Note that some author add a scaling factor of 1/N to that definition, but is not of great importance for this post — all in all it is just a matter of definition and sticking it to it.The inverse Fourier transform is then (given the definition for the forward Fourier transform):That being said, one of the most important theorem on Fourier transform is that convolution in one space is equivalent to multiplication in the other. In other words, the Fourier transform of a product is the convolution of the respective Fourier spectra, and the Fourier transform of a convolution is the product of the respective Fourier spectra.andwhere the dot represents the standard product (multiplication) and the circled star represents circular convolution.Two important notes:If you’re familiar with linear convolution, often simply referred to as ‘convolution’, you won’t be confused by circular convolution. Basically, circular convolution is just the way to convolve periodic signals. As you can guess, linear convolution only makes sense for finite length signals, that do not span from minus infinity to infinity. In our case, in the context of Fourier analysis, our signals are periodic therefore do not satisfy this condition. We can’t talk about (linear) convolution.Yet we still can intuit a linear-convolution-like operation on our periodic signals : just convolve the periodic signal on a one-period length. That’s what circular convolution does : it convolves 2 periodic signals of the same length along a one-period span.To further convince yourself of the difference, compare both formulas for discrete linear convolution and discrete circular convolution :Notice the differences :- bounds : linear convolution uses samples from minus-infinity to plus infinity — as stated previously, in this context x and y have finite energy the sum makes sense. For the circular convolution, we only need to what happened in a period span, so the sum only spans one period- circular indexes : in the circular convolution, we “wrap” the index of y using a modulo operation with a length of N. That’s just a way to ensure that y is considered periodic with a period of N : when we want to know the value of y at position k, we just use the value of y at position k%N — since y is N periodic, we get the right value. Again, this is just a mathematical way to handle periodic infinite-length sample sequences.Numpy provides great tools for finite length signals: and that’s a good news because as we just saw, our infinite-length periodic signal can be represented with just a period.Let’s create a simple class to represent such signals. We add a method to quickly plot the array, as well as an additional period before and after the “base” array, to keep in mind that we are dealing with a periodic sequence.Let’s look at two examples : first with a sampled sinus sequence, then with a linear sequence. Both are considered to be N-periodic (N=10 in this case).Let’s now implement the circular convolution equation seen above. Using indexing and the modulo operator, it is pretty straightforward:That’s great, we can now see what the circular convolution between two signals looks like. Putting everything in a single figure :Now this solution works really well, but it has a major flaw : it is slow. As you can see, we have to go through 2 nested loops to compute the result : one for each position in the result array, and one to compute the result at that position : we say that the algorithm is O(N²), as N increases the number of operations will increase as the square of N.For small arrays like those in the example, it is not a problem, but as your array will grow it’ll become a major problem.Furthermore, loops on numerical data is most of the time considered a bad practice in python. There must be a better way…And that’s where the Fourier transform and the convolution theorem come into play. Because of the way the discrete Fourier transform is implemented, in a very fast and optimized way using the Fast Fourier Transform (FFT), the operation is **very** fast (we say the FFT is O(N log N), which is way better than O(N²)).Using the convolution theorem, we can use the fact the product of the DFT of 2 sequences, when transformed back into the time-domain using the inverse DFT, we get the convolution of the input time sequences. In other words, we have :where DFT represents the discrete Fourier transform, and IDFT the inverse operation.We can then implement this algorithm to compute the convolution of x and y using numpy very easily :To finish, let’s verify that both approaches yield the same results, and compare the time it take python to compute the circular convolution using both techniques:That’s a perfect match! Both are rigorously equivalent in term of numerical values.Now for the time comparison:And the results are:That’s huge! Consider now what it can brings you when you analyze your time-series with thousands and thousands of samples!We have seen in this post that the Fourier transform is powerful tool, especially thanks to the convolution theorem that allows us to compute convolution in a very efficient manner. We’ve seen that linear and circular are not exactly the same operation but are both based on a convolution.Also, check out my other post and if you like any of them, please subscribe it helps me a lot to reach my goal of 100 subscribers:300-Times Faster Resolution of Finite-Difference Method Using NumPy | by Yoann Mocquin | Towards Data Science (medium.com)PCA/LDA/ICA : a components analysis algorithms comparison | by Yoann Mocquin | Towards Data Science (medium.com)Wrapping numpy’s array. The container approach. | by Yoann Mocquin | Towards Data Science (medium.com)Deep Dive into Seaborn: Color Palettes | by Yoann Mocquin | Analytics Vidhya | Medium----Towards Data SciencePhysics engineer by training, python by passionYoann MocquininTowards Data Science--Khuyen TraninTowards Data Science--22Miriam SantosinTowards Data Science--13Yoann MocquininTowards Data Science--2Monit Sharma--Shawhin TalebiinTowards Data Science--4Monit Sharma--1Mathcube--1Mathcube--1Monit Sharma--HelpStatusWritersBlogCareersPrivacyTermsAboutText to speechTeams