{'fraction_non_alphanumeric': 0.04519257801518133, 'fraction_numerical': 0.014736200918376909, 'mean_word_length': 4.204853066699183, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Wavelet Compressibility of Compound Poisson Processes

In this paper, we precisely quantify the wavelet compressibility of compound Poisson processes. To that end, we expand the given random process over the Haar wavelet basis and we analyse its asymptotic approximation properties. By only considering the nonzero wavelet coefficients up to a given scale, what we call the greedy approximation, we exploit the extreme sparsity of the wavelet expansion that derives from the piecewise-constant nature of compound Poisson processes. More precisely, we provide lower and upper bounds for the mean squared error of greedy approximation of compound Poisson processes. We are then able to deduce that the greedy approximation error has a sub-exponential and super-polynomial asymptotic behavior. Finally, we provide numerical experiments to highlight the remarkable ability of wavelet-based dictionaries in achieving highly compressible approximations of compound Poisson processes.


I. INTRODUCTION A. Sparsity and the Limits of Gaussian Models
The statistical modelling of data plays a central role in numerous research domains, such as signal processing [1] and pattern recognition [2]. In that regard, Gaussian models have been the first and by far the most considered ones, thanks to their desirable mathematical properties and relatively simple characterization. For instance, the Karhunen-Loève transform (KLT) identifies the optimal basis for representing data with Gaussian priors [3] and Kalman filters are optimal denoisers of Gaussian signals [4], both in the mean-square sense. These facts, among others, have made Gaussian statistical priors very convenient in practice. They also reveal the fundamental relationship between Fourier-based signal representations and Gaussian models.
However, it has been a long standing observation that Gaussian models fail to capture several key statistical properties of most naturally-occurring signals [5], [6]. Indeed, the latter frequently have heavy-tailed marginals [7]- [10] or richer structure of dependencies than Gaussian ones [11], [12]. Realworld signals are highly structured and often admit concise representations, typically on wavelet bases that appear to be genuinely versatile [13], [14]. This has led to the current Shayan  paradigm in modern data science where sparsity plays one of the central roles in statistical learning [15], [16] and signal modelling [8], [17], [18]. Classical Gaussian priors cannot model sparsity as they tend to produce poorly compressible signals [19], [20]. Many recent efforts in signal processing have been directed towards the development of deterministic frameworks that are better tailored to the reconstruction or synthesis of sparse signals, such as traditional compressed sensing [21]- [23] and its infinite-dimensional extensions [24]- [27].

B. Wavelets and Signal Representations
The development of wavelet methods, based on the pioneering works of I. Daubechies, Y. Meyer, and S. Mallat in the late 80's [28]- [30], has shed new lights on signal representation. Repeated numerical observations confirmed that wavelet-based compression techniques such as JPEG-2000 [31] outperform classical Fourier-based standards (e.g., JPEG) for natural images. This is despite the fact that the discrete Fourier transform (DFT) and its real-valued counterpart, the discrete cosine transform (DCT) [32], are asymptotically equivalent to KLT and, hence, are optimal for representing signals with Gaussian prior [33].
Wavelets are celebrated for their excellent approximation properties for large classes of signals and functions [34]. They revived the field of functional analysis [30], [35], culminating with the Abel prize of Yves Meyer in 2017 and feeding remarkable applications to various scientific and engineering fields [36]. One of the remarkable aspects of wavelets is that they are unconditional bases for many function spaces, including Hölder, Sobolev, and Besov spaces [30], [35] which is a key property for studying the best M -term approximation in a given basis [34], [37].

C. Probabilistic Models for Sparse and Analog Signals
As we have seen, probabilistic models beyond the Gaussian paradigm are of special interest for the modelling of sparse signals. A systematic attempt in this direction has been developed in the monograph [8]. In this work, analog signals are modelled as solutions of stochastic differential equations driven by non-Gaussian Lévy white noises. The so-called sparse stochastic processes have been used to develop novel techniques for essential signal processing tasks, such as denoising [38] and estimation [39] for signals with non-Gaussian priors. These methods have also been used in biomedical image reconstruction [40], highlighting the practical aspects of this new statistical framework.
The simplest class of non-Gaussian model is the one of compound Poisson processes. The latter are random piecewise constant functions with independent and stationary increments. As such, they are part of the family of Lévy processes [41], which also includes the Brownian motion. Compound Poisson processes are fully determined by the heights and locations of their countably many jumps. Contrary to Brownian motion, they are part of the class of signals with finite rate of innovations [42]- [44], meaning that their realizations on compact intervals can be fully encoded by finitely many numbers. This makes them particularly appealing for the modelling of highly-compressible piecewise constant signals. It has also been shown recently that any Lévy process is the limit in law of compound Poisson processes whose rate of innovation tends to infinity [45]. This theoretical observation permits the development of methods for generating trajectories of Lévy processes from compound Poisson processes, as exploited in [46].

D. Gaussian versus Poisson: Two Extreme Compressibility Behaviors
The aforementioned class of Lévy processes (see Section II-A for a formal definition) allows for various compressibility behaviors: the Brownian motion is the less compressible, while the compound Poisson ones are at the other extreme. This compressibility hierarchy has been recently revealed in two different theoretical frameworks.
In the first one, the compressibility is measured via the speed of convergence of the best M -term approximation in wavelet bases. The decay rate of the best M -term error is known to be directly linked to the Besov regularity [34], [47], which has been quantified for a broad class of Lévy processes [48]- [53]. Hence, the compressibility of Lévy processes has already been characterized using this approach [54], [55] and synthesized in [56,Chapter 6]. In a nutshell, state-of-the-art results show that the best M -term quadratic approximation error of the Brownian motion behaves asymptotically like 1/M 1 , while the same quantity decays faster than any polynomial for compound Poisson processes [55,Theorems 4 and 5].
In the second framework, the compressibility of a Lévy process is quantified in the information theoretic sense through the entropy of the underlying Lévy white noise, as in [57], [58]. These two frameworks are complementary and based on totally different tools, but they are consistent and lead to the same compressibility hierarchy.

E. Contributions and Outline
This paper contributes to the analysis of the compressibility of Lévy processes, focusing on the compound Poisson and Gaussian cases. We consider the Haar wavelet approximations of these random processes and quantify the decay rate of their approximation error in the mean squared sense.
More precisely, we focus on quantities such as  I  THE DECAY RATE OF THE MSE (1) OF LINEAR, GREEDY, AND BEST  M -TERM APPROXIMATION SCHEMES FOR COMPOUND POISSON  PROCESSES AND THE BROWNIAN MOTION. Brownian Motion Compound Poisson where s is a compound Poisson process or the Brownian motion and P M : is a possibly nonlinear approximation operator based on M ≥ 1 Haar wavelet coefficients of the input function. We compare various approximation schemes, depending on which wavelet coefficients are chosen. The two best-known schemes are the linear and the best M -term approximation, albeit both suffer from practical limitations. On one hand, the linear scheme does not capture the sparsity that might be inherent in the signal of interest (see Proposition 2). On the other hand, in order to exactly implement a compression scheme based on the best M -term approximation of the random process, one needs to have access to the full infinite set of wavelet coefficients. Without additional knowledge on the wavelet expansion, the implementation may become cumbersome and not memory efficient, if not impossible. This is why alternative approximation schemes have been proposed, most notably the "tree approximation" scheme which has brought significant attention in the literature [59]- [62].
In the same spirit, we consider a very simple greedy approximation scheme, in which only the first M nonzero wavelet coefficients are preserved. This scheme is well-suited to compound Poisson processes, for which most of the wavelet coefficients are zero due to their piecewise constancy.
Our main result is to provide lower and upper-bounds for the greedy approximation error in the mean-squared sense (Theorem 1). It essentially states that the mean-square error of the Haar greedy approximation of the compound Poisson process s behaves roughly as where N is the (random) number of jumps of s (see (24) for the precise meaning of (2)). This allows us to deduce that the mean-square error decays faster than any polynomial, and slower than any exponential (Theorem 2). We also perform a similar analysis for the linear approximation of the compound Poisson process, as well as for the linear and greedy approximations of the Brownian motion. This highlights the specificity of the compound Poisson processes: the greedy approximation dramatically outperforms the linear scheme for compound Poisson processes, contrary to the Gaussian case. We summarize this situation in Table I, where the main contribution is highlighted in bold. Finally, we illustrate our theoretical findings with numerical examples in various cases. Specifically, we highlight that the approximation error obtained within our method is close to the best M -term approximation. Moreover, we highlight the role of the wavelet dictionary by comparing the linear and best M -term schemes for compound Poisson processes and the Brownian motion in a Fourier-type dictionary corresponding to the discrete cosine transform (DCT). These empirical observations raise interesting theoretical questions which we briefly expose and can be exploited in future works.

F. Outline
The paper is organized as follows: in Section II, we present the relevant mathematical concepts. We then discuss our approximation scheme and compare it with the linear and best M -term methods in Section III. We present our main theoretical results in Section IV and finally, we demonstrate our theoretical results within numerical examples in Section V.

II. MATHEMATICAL PRELIMINARIES
In this section, we recall the relevant mathematical concepts and state preparatory results that we will use throughout the paper.

A. Lévy Processes and Lévy White Noises
Brownian motions and compound Poisson processes are members of the general family of Lévy processes, which are continuous-domain random processes characterized by their independent and stationary increments [41], [63]. Lévy processes are defined 2 , is not the most common one. However, it is proven to be equivalent to more classical ones in [64]. Consequently, we shall view the wavelet coefficients as linear functionals acting on the Lévy white noises which in effect allows us to simply characterize their probability laws. as the solutions of the stochastic differential equation with the boundary condition s(0) = 0. In (3), D denotes the (weak) derivative operator and w is a Lévy white noise. We choose to only consider zero-mean white noises. However, this comes with no loss of generality as the results are readily extendable to the general case. The formal construction of the family of Lévy white noises as generalized random processes have been exposed in [65,Chapter 3]. In this framework, the Lévy white noises are defined based on their observation through smooth test functions. For each adequate test function ϕ, w, ϕ is then a zeromean random variable. The collection of random variables ( w, ϕ ) ϕ satisfy two important properties: • (Stationarity) For any test function ϕ and any shift value t 0 ∈ R, the random variables w, ϕ and w, ϕ(· − t 0 ) have the same law. • (Whiteness) For any pair of test functions (ϕ 1 , ϕ 2 ) with disjoint support, the random variables w, ϕ 1 and w, ϕ 2 are independent. The class of valid test functions for a Lévy white noise have been characterized in [66], [67]. It is sufficient for us to know 2 This definition, based on the theory of distributions that w, ϕ is well defined for any Lévy white noise w and any square-integrable and compactly supported test function ϕ [67, Proposition 5.10].
The most studied example of Lévy processes is the Brownian motion, for which w is a Gaussian white noise. In this case, for any ϕ ∈ L 2 (R), the random variable w, ϕ has a normal distribution with zero mean and variance σ 2 ϕ 2 2 , where σ 2 is the variance of the noise [65, Section 2.5].
Another prominent subfamily of Lévy processes are the compound Poisson processes. They are piecewise constant processes and their statistics is characterized by their probability law of jumps P and their Poisson parameter λ > 0 that controls the sparsity of the random process (see Figure 1). More precisely, the compound Poisson white noise w with law of jumps P and Poisson parameter λ > 0 can be written as a sum of non-uniform Dirac impulses, as where the sequence {a k } k∈Z of height of Diracs is i.i.d. with law P and the sequence {τ k } k∈Z of locations of Diracs is a stationary Poisson point-process with parameter λ > 0 (see [68] for a formal definition of point processes), the a k and the τ k being independent. The key property regarding the Dirac locations is that the number for any vector u = (u 1 , . . . , u n ) ∈ R n . It is worth noting that the probability density of x, once condition to N = n ≥ 1, does not depend on λ anymore. Throughout the paper, we shall write the ordered jump positions of compound Poisson processes with the letter x, and the unordered ones with the letter τ .
In Lemma 1, we characterize the law of the minimal distance ∆ between two consecutive jumps of a compound Poisson process. The proof of Lemma 1 can be found in Appendix A.
Then, almost surely, ∆ ≤ (b − a)/N , and for any n ≥ 1 and δ ∈ [0, (b − a)/n], we have In this paper, we shall consider compound Poisson processes and white noises that are zero-mean, finite variance (which is equivalent to say that the jumps themselves are zero-mean with finite variance), and whose probability law of jumps P has a PDF (in particular, it has no atoms, what will be used in our analysis). The prototypical example is a compound Poisson process with Gaussian jump heights.
Despite the fact that their sample paths have very distinct behaviors (see Figure 1), finite-variance compound Poisson processes have the same second-order statistics as the Brownian motion. Indeed, for any test function ϕ, the random variable w, ϕ has zero-mean and variance σ 2 ϕ 2 2 for any Lévy white noise with variance σ 2 and zero mean [8,Proposition 4.15].

B. Haar Wavelets
For a pair of functions ψ, φ ∈ L 2 (R), that are referred to as the mother and father wavelets, respectively, the wavelet family Ψ contains all (normalized) dyadic scales and integer shifts of ψ plus the integer shifts of φ. In other words, we have for all scaling factor j ≥ 0 and all shifting parameter k ∈ Z.
We consider the family of Haar wavelets whose mother and father wavelets are respectively Haar wavelets are known to form an orthonormal basis for L 2 (R) [28]. This means that any function f ∈ L 2 (R) admits the unique expansion where ·, · denotes the standard inner product in L 2 (R), defined as ϕ 1 , The simple characteristics and implementation of Haar wavelets make them favorable in practice [69], [70]. They are also compactly supported, which is of great importance in our analysis, due to the whiteness property of Lévy white noises (see above). Last but not least, the family consists of piecewise constant functions. Hence, it is natural to represent compound Poisson processes (that are themselves almost surely piecewise constant) in this basis.

C. Haar Decomposition of Lévy Processes
In the sequel, we restrict both the random processes and the wavelet transforms to [0, 1] and study the local compressibility of compound Poisson processes over this compact interval.
Due to the support localization of the Haar wavelets, we readily see that the family Ψ = {ψ j,k } j≥0,0≤k≤2 j −1 ∪ {φ} forms an orthonormal basis of L 2 ([0, 1]), hence the Lévy process s can be almost surely written as The probability law of the Haar wavelet coefficients of s has been characterized in [9], where their characteristic functions have been explicitly computed. Here, we study the law of wavelet coefficients using the properties of the underlying Lévy white noise. In order to achieve this goal, we introduce the auxiliary functions defined for t ∈ [0, 1] as for any j ≥ 0 and k = 0, . . . , 2 j − 1. We conclude this part with Proposition 1, that expresses the Haar wavelet coefficients of s using the underlying Lévy white noise and the auxiliary functions (12) and (13). The proof is available in Appendix B.
Proposition 1. Let s be a Lévy process. Then, for any j ≥ 0 and 0 ≤ k ≤ 2 j − 1, we have where w is the Lévy white noise such that Ds = w.

III. WAVELET-BASED APPROXIMATION SCHEMES
In this section, we consider three different approximation schemes for square-integrable functions over [0, 1]: the linear, best M -term, and greedy approximations. Our main goal and contribution is to precisely quantify the approximation power of the greedy scheme.

A. Wavelet-Based Approximation Schemes
In what follows, we consider the natural indexing of wavelets by defining the indexing function Ind : Ψ → N as for all j ≥ 0 and k = 0, . . . , 2 j − 1.
We denote by the linear approximation of f , that is obtained by keeping the first M wavelet coefficients (with respect to the indexing function Ind) of f in the expansion (11). • P best M (f ), the best M -term approximation of f , that is obtained by keeping the M largest wavelet coefficients of f .
The first scheme in Definition 1 is called linear due to the fact that P lin M (f ) depends linearly on f . However, the best M -term approximation is adaptive to the signal and is therefore nonlinear. One can hope that the adaptiveness of the best M -term approximation significantly improve the quality of the approximation when compared with the linear one, what appears to be the case for some classes of functions [34].
As an alternative approach, we consider a compression scheme for compound Poisson processes that can be performed in an online fashion with respect to the stream of the wavelet coefficients. The main idea is to exploit the tremendous sparsity of the expansion of compound Poisson processes over the Haar wavelet basis, that is done by retaining only the nonzero wavelet coefficients and is called the greedy approximation.
the greedy approximation of f , where only the M first nonzero wavelet coefficients are preserved (the ordering being understood with respect to the indexing function Ind in (15)).
As for the best M -term approximation, the greedy approximation of f is nonlinear with respect to f . However, it is greedy in the sense that it can be computed by simply looking at the ordered wavelets coefficients. Hence, it does not necessitate to observe the complete set of wavelet coefficients, contrary to the best M -term approximation. It therefore shares the simplicity of the linear scheme and the adaptiveness of the optimal scheme (the best M -term).
The three approximation schemes introduced in this section clearly satisfy the relations for any function f ∈ L 2 ([0, 1]) and any M ≥ 0.

B. Mean-Squared Error of the Wavelet Approximations
Let s be a Lévy process. To quantify the performance of an approximation scheme, we consider the mean-squared error (MSE), which we denote by MSE method M for the approximation scheme method ∈ {lin, greedy, best} and is defined as It is clear from (16) that C. The Linear Scheme In Proposition 2, we determine the MSE lin M of any Lévy process that has finite variance. Its proof is available in Appendix C.
Proposition 2. Let s be a Lévy process with finite variance σ 2 0 . Then, for every M ≥ 1, we have where Proposition 2 shows that the linear approximations of Lévy processes with finite variance share the same mean-square error. Let us also remark that if s is a Brownian motion, then the random variables X j,k = s, ψ j,k = w,ψ j,k are all Gaussian. Hence, P(∃j, k : X j,k = 0) ≤

D. The Greedy Approximation of Compound Poisson Processes
When the wavelet coefficients are sparse (i.e. when at each scale, only a few of them are nonzero), the linear and greedy approximation schemes are no longer identical. In Proposition 3, we study the sparsity of the wavelet coefficients of compound Poisson processes. Precisely, we first characterize when a specific wavelet coefficient vanishes, depending on the presence of jumps. Using this primary result, we provide upper and lower bounds for the minimal (random) scale at which at least M wavelet coefficients are nonzero. The proof of Proposition 3 is provided in Appendix D.
Proposition 3. Let s be a compound Poisson process whose law of jumps admits a PDF with zero-mean and finite variance.
1) For all j ≥ 0 and k = 0, . . . , 2 j − 1, denote K j,k as the random number of jumps of s in the support of ψ j,k . Then, we almost surely have s, ψ j,k = 0 ⇔ K j,k = 0.
In other words, the symmetric difference between the events s, ψ j,k = 0 and K j,k = 0 has probability zero. 2) Consider the wavelet expansion (11) of s and denote by N J , the random number of nonzero wavelet coefficients with scale no larger than J. Furthermore, condition to {N ≥ 1}, let J M be the smallest random value of J such that N J ≥ M ; that is, J M is characterized by where the random variable ∆ is defined in (6).

IV. COMPRESSIBILITY OF COMPOUND POISSON PROCESSES
In this section, we present our main result on characterizing the asymptotic behavior of the greedy approximation of compound Poisson processes.
Theorem 1. Let s be a compound Poisson process with Poisson parameter λ > 0 whose law of jumps admits a PDF with zero-mean and finite variance. Then for every M ∈ N, we have that where N is a Poisson random variable with parameter λ, and C 1 , C 2 > 0 are some constants.
The proof can be found in Appendix E. Here, we give a sketch of the proof. For an arbitrary fixed integer n ≥ 1, we work conditionally to N = n. From the definition of J M (see Proposition 3), one has that N J M −1 ≤ M − 1. Hence the M th nonzero wavelet coefficient is reached at scale J M , and therefore almost surely. From Proposition 2, we know the exact behavior of the linear approximation error. On the other hand, we have lower and upper-bounds for the quantity J M , thanks to Proposition 3. The rest of the proof leverages these two preliminary results in order to derive the announced bounds.
Theorem 1 provides lower and upper bounds for the greedy approximation error of any finite-variance compound Poisson process. In Theorem 2, we use these bounds to deduce subexponential super-polynomial behaviors for the greedy approximation error of compound-Poisson processes.
Theorem 2. Let s be a compound Poisson process whose law of jumps admits a PDF with zero-mean and finite variance. Then the greedy approximation error MSE greedy M of s follows a sub-exponential and super-polynomial asymptotic behavior. Precisely, for any k ∈ N, we have that and for any α > 0, The proof of Theorem 2 can be found in Appendix F. An enlightening consequence of the super-polynomial behavior of the greedy approximation error is that it demonstrates that our provided lower-and upper-bounds are asymptotically comparable. Specifically from the upper-bound provided in Theorem 1, we deduce that for any > 0. Moreover, Theorem 2 implies that the quantity M E[2 − M N ] tends to 0 as M → +∞ and is therefore bounded from above. Using a similar argumentation for the lowerbound of Theorem 1, we obtain the following corollary.
Corollary 2. For any > 0, there are positive constants C 1, , C 2, > 0 such that for all values of M ≥ 1.
Our theoretical analysis validates the two following observations in a rigorous manner: • A piecewise constant function with a fixed number of jumps n ≥ 1 is such that its greedy approximation in the Haar basis roughly behaves like O(2 −M/n ), which is exponential and therefore decays to 0 faster than any polynomial. Note that the exponential decay is faster for smaller values of n. • The number of jumps N of a compound Poisson is random. It is almost surely finite but can be arbitrarily large. The concrete effect is that the mean-square error of the greedy approximation roughly behaves like O(E[2 −M/N ]). The subexponential behavior of the MSE is then a consequence, as we have shown. It is worth noting that the characterization provided by Theorem 2 is not deducible from earlier works that was based on the machinery of Besov regularity, such as [55]. Previous works focus on the almost sure behavior of the approximation error, while we focus on the mean-square approximation on this paper. These are two different regimes and to the best of our knowledge, the possible link between the two has not been investigated.
By contrast, we obtain some information regarding the asymptotic behavior of best M-term approximation error of compound Poisson processes from Theorem 2. Indeed, by combining (18) and (26), one observes that for any k ∈ N. Using the fact that MSE is non-negative simply implies that lim The super-polynomial decay of the greedy approximation error shows that this method, despite being very simple and easily implemented, reaches excellent approximation performances.
In the next section, we will empirically show that the greedy scheme performs similarly to the practically uncomputable best M -term approximation scheme.
V. NUMERICAL ILLUSTRATIONS In this section, we provide a numerical demonstration of the main results of this paper. First, it is illustrative and reflects the potential practical impact of our theoretical claims in a complementary and empirical manner. Second, it shows that the results obtained for the greedy approximation method are similar to what would be obtained for the best Mterm approximation (Section V-B). Finally, it emphasizes that wavelets are able to exploit the inherent sparsity of non-Gaussian signals, which is not the case of traditional Fourierbased approximation schemes (Section V-C). These empirical observations then give rise to open theoretical questions that might be of interest to the community.
To simulate each approximation scheme, we first generate a signal that consists of 2 10 = 1024 equispaced samples of a given random process over [0, 1]. We then compute its (discrete) Haar wavelet coefficients of scale up to J max = 10. Finally, we create the approximated signal according to the given approximation scheme 3 . We repeat each experiment 1000 times and we report the average to reduce the effect of the underlying randomness (Monte Carlo method). The averaged values are then good approximations of the quantities of interest, that is, the MSEs given by (17) for different approximation schemes.

A. Greedy Approximation Error
In the first experiment, we compute the MSE of greedy approximation for Brownian motions and compound Poisson processes with different values of λ = 10, 50, 100, 500 and with Gaussian jumps, as a function of the number M of coefficients that are preserved. We recall that, N being the random number of jumps of the compound Poisson process s over [0, 1], λ = E[N ] is the averaged number of jumps. To have a fair comparison, we unify the variance of the random processes in all cases to be σ 2 0 = 1 (which corresponds to a law of jumps with variance σ 2 0 /λ = 1/λ for compound Poisson processes).
The results are depicted in Figure 2, where in each case we plot the MSE in log scale, that is log 2 (M ) → 10 log 10 (MSE). From Proposition 2 and Corollary 1, we expect that the MSE of Brownian motion follows a global linear decay in the log scale, while decaying sub-linearly locally. Indeed, for M = 2 J , J ∈ N, we deduce from (20) that where α = 10 log 10 (σ 2 0 /6) and β = 10 log 10 (2) which shows a linear decay with respect to J = log 2 (M ). However, in the regime when J = log 2 (M ) is fixed, that is when 2 J ≤ M < 2 J+1 , we obtain from (19) that 10 log 10 (MSE greedy which shows that the error decays sub-linearly in this regime. These theoretical claims can be observed in Figure 2, as well. In addition, from Theorem 2, we know that the MSE of compound Poisson processes in the log scale should asymptotically decay faster than any straight line. This is also observable in Figure 2, indicating the dramatic difference between the compressiblity of compound Poisson processes and Brownian motions, as expected. We moreover remark in Figure 2 that the small-scale behavior (log 2 (M ) = J ≤ 3) does not distinguish between different values of λ, but also between compound Poisson processes and the Brownian motion. Again, this empirical fact has a theoretical counterpart: it is linked with the fact that the statistics of finite variance compound Poisson processes are barely distinguishable from the ones of the Brownian motion at coarse scales. This has been formalized in [71] which states, when particularized to our case, that compound Poisson processes with finite variance converge to the Brownian motion when zoomed out and correctly renormalized. Our numerical experiments are illustrative to this point, and will be confirmed in Sections V-B and V-C.
Finally, we observe in Figure 2 that as λ → +∞, the greedy approximation of compound Poisson processes converges pointwise to the one of Brownian motion. This empirical observation poses an interesting theoretical question which is also consistent with [45,Theorem 5], which states-when specialized to our problem-that the compound Poisson process with constant variance σ 2 0 and Gaussian jumps converges in law to the Brownian motion when λ → ∞. As we have seen in the introduction, it is particularly satisfactory to characterize the compressibility of Lévy processes via their best M -term approximation error in a given basis. Although our greedy approximation error only provides an upper-bound for the best M -term approximation error, we demonstrate numerically in Figure 3 that the two approximation schemes are comparable in the sense of MSE. This is also an important observation, as it reveals that the extremely simple greedy approximation performs almost as good as the best M -term approximation, the latter being a theoretical bound for M-term approximation schemes.

C. Haar vs. Fourier
We now investigate the effect of the dictionary in which we perform the approximation scheme. We consider the Haar transform and discrete cosine transform (DCT) for approximating the Brownian motion and compound Poisson processes with Gaussian jumps. The results are depicted in Figure 4, where we plot the best M -term approximation error of each setup in the log scale.
We observe that the DCT works slightly better than Haar for the Brownian motion. This is not surprising: The DCT is known to be asymptotically equivalent to the Karhunen-Loève transform (KLT), which is optimal for Gaussian stationary processes [33]. It is worth noting that this is also valid for the Brownian motion, which is not stationary but still admits stationary increments.
However, there is a dramatic difference between Haar and DCT for compound Poisson processes. We see in Figure 4 that, contrary to the Haar dictionary, the DCT is unable to take advantage of the effective sparsity of compound Poisson processes. This is of course not a surprise and is folklore knowledge, but it has not yet been justified theoretically for the best of our knowledge. This is nevertheless consistent with recent theoretical and empirical results demonstrating that wavelet methods outperform classical Fourier-based methods for the analysis of sparse stochastic processes [8], [55].

VI. CONCLUSION
The theoretical and empirical findings of this paper are reminiscent to the so-called "Mallat's heuristic" [37], which states that "Wavelets are the best bases for representing objects composed of singularities, when there may be an arbitrary number of singularities, which may be located in all possible spatial position." and which remarkably describes the compound Poisson model.
To do so, we provided a theoretical analysis to characterize the compressibility of compound Poisson processes. To that end, we introduced a simple approximation greedy scheme performed over the Haar wavelet basis. We then provided comparable lower and upper-bounds for the mean-squared approximation error. This enabled us to deduce the subexponential super-polynomial asymptotic behavior for the error. Future research direction is to investigate the compressibility of compound Poisson processes in other dictionaries (e.g. DCT or other wavelet families), to investigate the effect of the Poisson parameter λ in this analysis, specifically when λ → +∞, and finally, to theoretically compare the best and greedy approximation schemes.

ACKNOWLEDGMENT
The authors are extremely grateful to Prof. Michael Unser, who strongly inspired this work on its early stage. They also warmly thank Laurène Donati for her kind help during the writing process.

APPENDIX A PROOF OF LEMMA 1
Proof. We first remark that the inequality ∆ ≤ (b − a)/N is obviously true when N = 0, since ∆ = (b−a) in this case. As for N ≥ 1, we have by definition of ∆ that ∆ ≤ x i − x i−1 , for all i = 1, . . . , N , with the convention that x 0 = a. By summing up these equality for for all values of i, we obtain that This yields that ∆ ≤ (b − a)/N . For the second part, we define the random vector d = (d 1 , . . . , d n ) ∈ [0, 1] n as By rewriting (31) in the vectorial form, we obtain that where e 1 = (1, 0, . . . , 0) ∈ R n , x = (x 1 , . . . , x n ) and H ∈ R n×n is the lower-bidiagonal matrix Now, due to (5) and the change of variables (32), the PDF of d is where v 1 = |v 1 | + · · · + |v n | = v 1 + · · · + v n for v ∈ [0, 1] n . In addition, from the definition of ∆, the probability of {∆ ≥ x} for any x ∈ [0, (b − a)/n] can be computed as where the latter is obtained via the change of variable u i = v i − x b−a for i = 1, . . . , n. We remark that if u i ≥ 0 for i = 1, . . . , n and u 1 ≤ 1 − nx/(b − a), then we would have a) for any i = 1, . . . , n. In other words, the upper-limit of the integral in (35) is redundant and can be replaced with +∞. Doing so, we obtain that where (i) is due to the symmetry of the integrand with respect to the sign of u and where Leb denotes the Lebesgue measure. Finally, we use a known result stating that the volume of the 1 unit ball in R n is 2 n /n! [72]. This yields to
With a similar idea, we remark that Dφ = δ − φ. Combining with s, δ = s(0) = 0, we have that APPENDIX C PROOF OF PROPOSITION 2 Proof. One observes from Definition 1 that This together with (11) yields that Haar wavelets that are supported in [0, 1], form an orthonormal basis for L 2 ([0, 1]). Using this, we express the approximation error based on the wavelet coefficients, as By taking expectation over both sides and by using Proposition 1, we have that Finally, we replace ψ j,k 2 L2 = 2 −2j 12 (obtained via a direct computation; see Figure 5 for visualisation) for all j ≥ 0 and k = 0, . . . , 2 j − 1 in the summation above to deduce that

APPENDIX D PROOF OF PROPOSITION 3
Proof. Item 1) Assume that K j,k = 0. This means that s is constant over the support of ψ j,k , taking the fixed (random) value s 0 . By recalling that for all j ≥ 0 and k = 0, . . . , 2 j − 1, we deduce that the corresponding wavelet coefficient is s, ψ j,k = s 0 R ψ j,k (x)dx = 0. For the converse, we show that, condition to the event {K j,k = K} for an arbitrary (but fixed) integer K ≥ 1, we have that P ( s, ψ j,k = 0|K j,k = K) = 0.
Consider the jumps that are inside the support of ψ j,k and denote their (unordered) locations and heights by {τ 1 , . . . ,τ K } and {ã 1 , . . . ,ã K }, respectively. Due to (14), we have that We recall that the jump locationsτ i are i.i.d. with a uniform law. Moreover, the jump heightsã i are independent ofτ i s and are themselves i.i.d. copies of a random variable that admits PDF. This implies that the random variables Z i =ã iψj,k (τ i ) for i = 1, . . . , K are also i.i.d. and their law has a PDF too, which we denote by p Z . Finally, the random variable s, ψ j,k also has PDF (that is the K times convolution of p Z with itself) and thus, is nonzero with probability one (no atoms).
Item 2) Recall that N is the total number of jumps of s over [0, 1]. Due to (22) and the fact that the wavelets ψ j,k for k = 0, . . . , 2 j −1 have disjoint support, at each scale j ≥ 1, at most N wavelet coefficients are nonzero. On the other hand, the support of any wavelet function of scale j is of size 2 −j . Hence, due to the definition of ∆, the number of jumps in the support of ψ j,k is either one or is upper-bounded by the length of the interval divided by the minimum distance (= 2 −j ∆ −1 ). In other words, the support of each wavelet of scale j contains at most max(1, 2 −j ∆ −1 ) jumps.
Denote by n j , the number of nonzero wavelet coefficients in the jth scale. Using the previous observation, we deduce for all j ≥ 1 that As for j = 0 (mother and father wavelets), we deduce similar to Item 1) that condition to N ≥ 1, we have n 0 = 2. By defining J lim = log 2 (∆ −1 ) , one readily verifies that for j ≤ J lim , we have min(1, 2 j ∆) = 2 j ∆. By contrast, min(1, 2 j ∆) = 1 for j ≥ J lim . Using these simple observations and by summing up lower-bounds of (40) for j = 1, . . . , J (together with n 0 = 2), we obtain, since J j=0 n j = N J , that To simplify the first lower-bound, we use the inequality 2 x ≥ x for x = J + 1 + log 2 ∆, which results to As for the the second lower-bound, we use It is now readily to verify that the two lower-bounds in (41) and (42) are indeed equal and hence, we have that Finally, using ∆N ≤ 1, we conclude that We follow the same principle to obtain an upper-bound for N J as well. By summing up upper-bounds of (40) for j = 1, . . . , J, together with n 0 = 2, we obtain that from which we deduce the upper-bound We complete the proof of (23) by combining (46) and (47), knowing that J M ∈ N.

APPENDIX E PROOF OF THEOREM 1
Proof of Theorem 1. Let σ 2 0 < +∞ be the variance of the process s. We divide the proof and show each side of the inequality (24) separately.
Upper-bound: First, we show that for any n ≥ 1, we have Let us then work conditionally to N = n. From Proposition 3, we have (condition to N = n) that J M ≥ M −2 n . Thus, by combining with (25), we obtain that Taking expectation from both sides yields On the other hand, condition to N = n, we have the equality in law is the sequence of corresponding heights. Therefore, we have that where the random variables Z i = a iψj,k (τ i ) are i.i.d. copies of a zero-mean random variable. We recall that the law of jumps a i has zero-mean and variance σ 2 0 /λ. Hence, the second-order moment of Z i can be computed as where we used the independence of a i and τ i in (i) and the uniform law of τ i in (ii) and finally, we replaced ψ j,k 2 L2 = 2 −2j 12 in the last equality. Now, due to the independence of the Z i , we deduce that By substituting (50) in (49), we obtain that By taking the expectation, we obtain that where in the last inequality, we have used 2 − M n ≤ 1 and . Using (52) with t = ln2 (such that e t = 2) yields Hence, which is the announced upper-bound with the constant C 2 = 2σ 2 0 3λ (1 + e 2λ ). Lower-bound: Similar to the upper-bound, we show that for any n ≥ 1, we have the inequality which immediately implies the announced lower-bound. We treat the case N = 1 separately. Condition to N = 1, both wavelet coefficients of order zero (associated to mother and father wavelets) are nonzero. Moreover, for any j ≥ 1, there is exactly one wavelet coefficient of scale j that is nonzero. This implies that J M = M − 2 and in addition, we have that Similar to the proof of Proposition 2 and together with (50), we deduce that By defining δ = (2n 2 − 2n + 2) −1 > 0 (the precise value will be used later) and J = M −1 n +log 2 (δ −1 ) +1, we observe that Similar to the upper-bound, we consider the jumps of s in [0, 1] and we denote their (unordered) locations and heights by τ 1 , . . . , τ n and a 1 , . . . , a n , respectively. With regard to the convention τ 0 = 0, we consider the random variable∆ = min 0≤i<j<n−1 |τ i − τ j | and consequently, the event E = {∆ ≥ δ} ∩ {0 ≤ τ 1 , . . . , τ n−1 ≤ 1/2 − δ}.
This implies that condition to N = n, we have On the other hand, where we used the independence (condition to N = n) of jumps τ i and heights a i of s in (i) and we used the independence of a i from N and ∆ as well the fact that the law of a i has zero mean in (ii). By substituting E[a 2 i ] = where the latter is deduced from the independence of E and {1/2 ≤ τ n ≤ 1} (condition to N = n). By using Lemma 1 with a = 0 and b = 1/2 − δ (we remind that δ = (2n 2 − 2n + 2) −1 ), we can compute the conditional probability of the event E as P(E|N = n) = 1 − (n − 1)

APPENDIX F PROOF OF THEOREM 2
Proof. It is sufficient to prove that the quantity E[2 − M N ] has sub-exponential and super-polynomial asymptotic behavior.
Super-polynomiality: First note that there exists an integer number N 0 ∈ N such that for every n ≥ N 0 , we have P(N = n) ≤ 2 −n . We then consider the following decomposition for any M ≥ N 0 + 1 We separately show that each term of the previous decomposition decays faster than the inverse of any polynomial as M → ∞. For the first term, simply due to P(N = n) ≤ 1, we have that Now by fixing n 0 to be a sufficiently large integer so that log 2 eα− 1 n0 > 0, we deduce that the right hand side explodes.