{'fraction_non_alphanumeric': 0.020365398467683846, 'fraction_numerical': 0.003601597799751162, 'mean_word_length': 4.8401529636711285, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Interval propagation through the discrete Fourier transform

We present an algorithm for the forward propagation of intervals through the discrete Fourier transform. The algorithm yields best-possible bounds when computing the amplitude of the Fourier transform for real and complex valued sequences. We show that computing the exact bounds of the amplitude can be achieved with an exhaustive examination of all possible corners of the interval domain. However, because the number of corners increases exponentially with the number of intervals, such method is infeasible for large interval signals. We provide an algorithm that does not need such an exhaustive search, and show that the best possible bounds can be obtained propagating complex pairs only from the convex hull of endpoints at each term of the Fourier series. Because the convex hull is always tightly inscribed in the respective rigorous bounding box resulting from interval arithmetic, we conclude that the obtained bounds are guaranteed to enclose the true values.


Introduction
The analysis of engineering structures subject to harmonic stochastic excitation is often conveniently done in the Fourier domain [1]. The response of the structure is obtained for each harmonic without the need for step-wise numerical integration. This kind of analysis however, seems to be incompatible with the fact that real data time-series records are often affected by uncertainty. Some authors have successfully propagated interval uncertainty in dynamic structural analysis directly in the time domain [2]. In other instances the literature on the topic seem to contemplate only purely probabilistic methods to address the problem of characterization of the uncertainty of the power spectrum from real data, as shown in [3] and literature therein. It is currently very difficult to convert an imprecise time signal to its Fourier domain. An algorithm for the interval Fourier transform is much needed to bridge the divide between imprecise signals and the well established harmonic structural analysis, and to enlarge the spectrum of engineering applications that explicitly account for uncertainty.
The discrete Fourier transform (DFT) is used in a variety of different applications in science and engineering, such as in spectral analysis, random vibration, differential equations, data compression, signal processing, image processing, probabilistic programming and many more [1,4,5]. There are various algorithms available for transforming a signal with the DFT, of which the best known is probably the fast Fourier transfom (FFT), presented by Cooley & Tukey [6]. Due to increasing computational power, simulations and equivalent calculations can be carried out ever faster, which is particularly important for the numerical analysis. An overview of a variety of algorithms used can be found in abundance in [7,8].
Real data records are required for a large number of problems in engineering. A major challenge with data acquisition that should not be neglected is that sensors only work accurately within certain tolerances and are therefore subject to uncertainties. These uncertainties arise, for instance, due to damaged sensors, device failures and measurement errors, or if the sensors are not precisely calibrated. In addition, the data could be captured incorrectly due to threshold limitations. For more accurate simulation results when utilising real data records, uncertainties must be taken into account [3].
In this work, three different methods are investigated to provide bounds on the amplitude of the DFT: (1) a brute-force method, (2) a method based on the convex hull, (3) a method based on complex interval arithmetic and rigorous bounding [9][10][11][12]. These methods provide a direct relationship between interval signals and the corresponding amplitudes of the DFT. Thus, intervals can replace the exact values when investigating structures while the signal's uncertainty can be captured by these intervals.
When x is a time signal, each complex number z k ∈ C represents an harmonic of the signal with angular frequency 2πkn/N . The DFT function is sometimes denoted by F, while its interval extension by The DFT in practice is used to convert a time signal into its Fourier sequence ( Figure 1) to study the harmonic properties of the signal where each element corresponds to a frequency. The DFT is often used to reduce the dimensionality of the signal; the inverse Fourier transform is used to reconstruct the signal after compression or simply to convert the signal into the time domain.
The objective of this work is to provide best-possible bounds on the amplitude of F, a.k.a. amplitude spectrum, The interval signal of Figure 2 will be considered throughout the manuscript. A convenient way to evaluate the interval extension [F] is to split non-interval and interval components, as shown in (3). For example, when the signal is acquired using the same sensor (same precision) the interval uncertainty has the same value throughout. In such case the [F] becomes where ∆ = [−1, 1] is a unitary interval, and ξ ∈ R is the precision of the interval signal, expressed in the same units as the signal. This form can be used to assess the effect of the cumulation of uncertainty on the amplitude as more and more uncertainty components get added. Also the non-interval term can be moved to the left-hand side of the equal sign of (3) to yield a centred form of the interval Fourier transform.

Interval arithmetic
The interval extensions defined in (1), (2), and (3) can be evaluated with interval arithmetic. The generalized rules of arithmetic for addition, multiplication and power elevation are needed for the purpose. The addition + between any two intervals x, y with x, y ∈ R is The multiplication between an interval and a number a ∈ R, which symbol × is often omitted is The power elevation with exponent two is The square root of a positive real interval is With these rules the amplitude of the DFT interval extension can be evaluated using (2) without the need of a particular algorithm. The results provided by the interval arithmetic are guaranteed to enclose the exact bounds. However artefactual uncertainty due to the dependence problem makes the bounds obtained with interval arithmetic wider than they ought to be.
The interval extension of the DFT (1) can be computed using the complex addition operator. Consider two complex intervals z k = [z k , z k ], k = 1, 2, their addition is given by With the complex addition operator the program used to compute F can be recycled to obtain its interval extension. The input signal to the program will be an ordered sequence of intervals.  The only caveat in Listing 1, is that the length of the signal must be a power of two in order for the signal to be fully transformed to the Fourier domain. Slightly modified programs can be obtained with little extra effort to map signals that are not a power of two, see e.g. [1].

Notes on the dependence problem using interval arithmetic
The bounds on the amplitude of the interval extension F evaluated with interval arithmetic are rigorous in the sense that all the possible signals x ∈ x are mapped inside these bounds. However when the amplitude is computed, the bounds obtained with (2) will be suboptimal because of the dependence problem. The computation of (2) with interval arithmetic does not track the dependence between the real and the imaginary component, thus the bounds are outer-estimated. The bounds are rigorous because all of the possible dependencies between real and imaginary components are considered by the interval arithmetic bounds.

The brute-force method
The brute-force method is too inefficient to be used in practical applications, and it is herein used mainly for illustration purposes. The method tracks the dependency between real and imaginary components by constructing a binary tree of all the propagated endpoints at each iteration step. Clearly the number of endpoints increases exponentially in base two and with exponent given by the iteration number. The total number of iterations is the length of the signal, so the longer the signal the more inefficient is the brute-force method. Thus this method can only be used for research purposes either on a very short signal or on a signal with very few intervals. On a personal computer the recommended size of the problem is eight, more iterations may slow down the process significantly.

The selective method
We observe that only the endpoints on the convex hull on the current set of endpoints carry the information about the bounds. Numerical experiments performed using the brute-force method have confirmed this intuition, however the investigation is limited by the number of iterations allowed by the brute-force method. With interval arithmetic we can keep propagating the endpoints on the convex hull as well as the enclosing box for as many iterations as needed and for however large signals. Moreover, because interval arithmetic is the most efficient method its box can always be computed alongside the selective method.
The numerical experiments done running the selective method alongside interval arithmetic have confirmed and verified the intuition that the selective method yields rigorous best possible bounds in each tested case. Figure 4, shows on four different frequencies, the three methods altogether. The convex hull is the set of points resting on the set boundary, which in Figure 4 is depicted by a red closed line. The convex hull is obtained at each iteration by propagating only the endpoints on the boundary of the previous iteration step.
Interval arithmetic gives us the verified guarantee that the selective method is rigorous, i.e. does not underestimate the uncertainty. The fact that the selective method makes use of a subset of the endpoints ensures that its bounds are also best possible. In fact, it is always possible to index the convex-hull endpoints resting on the edge of the rigorous box and propagate these points back to identify the combination of endpoints yielding the bounds. However because there is no guarantee that theses endpoints will be unique, the algorithm cannot be used directly to perform the inverse transform.
The code produces the set of points in the convex hull as shown in Listing 4 and can be obtained as a derivation from the brute-force method. The function Selective in the code snippet Listing 4 will be invoked for each frequency k of the Fourier sequence.

The bounding-box method
There are no repeated variables in the interval extension of (1), thus the final enclosing box is the smallest box possible containing the set of points generated by the brute force method. In Listing 5 we show the procedure to obtain the enclosing interval box at each iteration of the Fourier transform. Note that if the enclosing box was not the smallest or best possible, there would have been a gap between the box and the convex hull of points, effectively making it impossible for us to verify the results produced by the selective method. Luckily in our numerical experiments there has been no instance of such a case. The rigorous bounding method can be performed alongside the selective method with very little additional computations, to provide a prompt verification certificate of the bounds outputted by the selective method. The Fourier amplitude spectrum can be obtained just using the information carried by the corners of the enclosing box. Because the set of extreme endpoints rest on the convex hull, the amplitude obtained with the enclosing box will carry inflated uncertainty. The bounds on the amplitude will be reflective of the fact that the zone outside the convex hull between the box will be propagated too, effectively making the bounds on the amplitude much puffier. The extra puffiness produced by this method however is compensated by the efficiency of this method, which practically carries no additional computational cost when compared to the precise DFT. In this section we show how to obtain the bounds on the amplitude for each frequency. A schematic representation of how the three methods interact with each other is provided in Figure 5. The cost of the selective method is O(N 2 log N ) because a convex hull is computed at each iteration. The divide and conquer algorithm for computing the convex hull has O(N log N ) cost, hence the total cost of the selective method.
The convex hull and the interval box at the last iteration for n = N − 1, can both be used to compute the amplitude at each frequency. We have shown that the endpoints on the convex hull yield the best possible bounds, thus the bounds on the amplitude are obtained by taking the minimum and maximum amplitude corresponding to this set of endpoints. The points corresponding to the extrema are shown in Figure 6. Because the amplitude of the Fourier transform is effectively the Euclidean norm on the R 2 isomorphism of the complex plane, minimum and maximum correspond to the points that are nearest and furthest from the origin of the complex plane. However, this is true as long as the convex hull does not enclose the origin. When the convex hull encloses the origin the minimum will be attained at the origin itself. The algorithm for the amplitude bounds will therefore assume that the minimum is attained at the origin of the complex plane in such situation, as shown in Figure 7.
A simple implementation of the algorithm for the amplitude bounds with the selective method is shown in Listing 6. The function AmplitudeBounds needs an external function called origin in complex hull that returns True when the convex hull encloses the origin and a False otherwise. The procedure for the amplitude spectrum bounds with the interval box is equivalent except that is much easier to check if the origin is contained in the box. Iterating over all the frequencies of the Fourier sequence as shown in Listing 7 leads to the amplitude bounds shown in Figure 8

Conclusions
We have presented three methods for computing the bounds of the Fourier amplitude spectrum when an interval signal is provided. We have shown how the problem can be tackled with a brute-force method, but only on relatively small interval signals (< 20 intervals), with limited applications in practice. Despite its inefficiency the brute-force method has provided the inspiration for the second more efficient selective method, which yields the best possible bounds on the amplitude. In this work, we have questioned the rigour and reliability of the selective algorithm, which has led us to our third method, which uses complex interval arithmetic to propagate the bounds through the Fourier series. The third method called bounding box, has provided a viable efficient tool to verify the rigorousness of the selective method even for large interval signals. The bounding-box method has also provided a robust alternative to obtain the amplitude bounds when efficiency is a priority at the expense of losing specificity in the bounds of the amplitude.

Replicability
The code, the example and the algorithms presented in this paper can be accessed in a single instance via Github at: https://github.com/marcodeangelis/Fourier-transform.