Wavelet Bases Made of Piecewise Polynomial Functions : Theory and Applications *

We present wavelet bases made of piecewise (low degree) polynomial functions with an (arbitrary) assigned number of vanishing moments. We study some of the properties of these wavelet bases; in particular we consider their use in the approximation of functions and in numerical quadrature. We focus on two applications: integral kernel sparsification and digital image compression and reconstruction. In these application areas the use of these wavelet bases gives very satisfactory results.


Introduction
In the last few decades wavelets and wavelets techniques have generated much interest, both in mathematical analysis as well as in signal processing and in many other application fields.In mathematical analysis wavelet bases, whose elements have good localization properties both in the spatial and in the frequency domains, are very useful since, for example, they consent to approximate functions using translates and dilates of one or of several given functions.In signal processing, wavelets were initially used in the context of subband coding and of quadrature mirror filters, but later they have been used in a variety of applications, including computer vision, image processing and image compression.The link between the mathematical analysis and signal processing approaches to the study of wavelets was given by Coif-man, Mallat and Meyer (see [1][2][3][4][5][6]) with the introduction of multiresolution analysis and of the fast wavelet transform, and by Daubechies (see [7]) with the introduction of orthonormal bases of compactly supported wavelets.
Let A be an open subset of a real Euclidean space and let   2 L A be the Hilbert space of the square integrable (with respect to Lebesgue measure) real functions defined on A .In this paper when A is a suffici-ently simple set (i.e. a parallelepiped in the examples considered), starting from the notion of multiresolution analysis, we construct wavelet bases of   2 L A with an (arbitrary) assigned number of vanishing moments.The main feature of these wavelet bases is that they are made of piecewise polynomial functions of compact support; moreover the polynomials used are of low degree and generate bases with an arbitrarily high assigned number of vanishing moments.This fact makes possible to perform very efficiently some of the most common computations, such as, for example, differentiation and integration.However the lack of regularity of the piecewise polynomial functions used can create undesirable effects in the points where the discontinuities occur when, for example, continuous functions are approximated with discontinuous functions.Note that the wavelet bases studied here, in general, make use of more than one wavelet mother function.Thanks to these properties these wavelet bases in several applications can outperform in actual computations the classical wavelet bases and, for example, in this paper we show that they have very good approximation and compression properties.The numerical results of Section 4 corroborate these statements both from the qualitative and the quantitative point of view.
The wavelet bases introduced generalize the classical Haar's basis, that has only one vanishing moment and is made of piecewise constant functions (see, for example, [8]), and are a simple variation of the multi-wavelets bases introduced by Alpert in [9,10].The results reported here extend those reported in [11,12] and aim to show not only the theoretical relevance of these wavelet bases (also shown, for example, in [9,10,13,14]) but also their effective applicability in real problems.In particular in this paper we study some properties of the wavelet bases considered and the advantages of using some of them in simple circumstances.In fact the orthogonality of the wavelets to the polynomials up to a given degree (i.e. the vanishing moments property) plays a crucial role in producing "sparsity" in the representation using these wavelet bases of functions, integral kernels, images and so on.These wavelet bases as the wavelet bases used previously have good localization properties both in the spa-tial and in the frequency domains and can be used fruit-fully in several classical problems of functional analysis.In particular we focus on the representation of a function in the wavelet bases and we present some ad hoc quadrature formulae that can be used to compute efficiently the coefficients of the wavelet expansion of a function.
We consider also the use of these wavelet bases in some applications, initially we focus on integral kernel sparsification.This is a relevant task, see for example [10,15], since it makes possible, among other things, the approximation of some integral operators with sparse matrices allowing the approximation and the solution of the corresponding integral equations in very high dimensional subspaces at an affordable computational cost.In [11,12,16], for example, we exploit this property to solve some time dependent acoustic obstacle scattering problems involving realistic objects hit by waves of small wavelength when compared to the dimension of the objects.Let us note that these scattering problems are translated mathematically in problems for the wave equation and that they are numerically challenging, moreover thanks to the use of the wavelet bases, when boundary integral methods or some of their variations are used, they can be approximated by sparse systems of linear equations in very high dimensional spaces (i.e.linear systems with millions of unknowns and equations).Later on we focus on another important application of wavelets: digital image compression and reconstruction.In this framework, the basic idea is to distinguish between relevant and less relevant parts of the image details disregarding, if necessary, the second ones.In particular we proceed as follows: a digital image is represented as a sequence of wavelet coefficients (wavelet transform of the original image), a simple truncation procedure puts to zero some of the calculated wavelet coefficients (i.e.those that are smaller than a given threshold in absolute value) and keeps the remaining ones unaltered (com-pression).The truncation procedure is performed in such a way that the reconstructed image (i.e. the image obtained acting with the inverse wavelet transform on the truncated sequence of wavelet coefficients) is of quality comparable with the quality of the original image, but the amount of data needed to store the compressed image is much smaller than the amount of data needed to store the original image.We present some interesting numerical results in wavelet-based image compression and reconstruction.Moreover we define a compression coefficient to evaluate the performance of the compression procedure and we study the behaviour of the compression coefficients on some test problems, in particular we show that the compression coefficients increase when the number of vanishing moments of the wavelet basis used increases.This property can be exploited in several practical applications.
The paper is organized as follows.In Section 2 using a multiresolution approach we present the wavelet bases.In Section 3 some mathematical properties of the wavelet bases introduced are discussed and some applications of these bases to function approximation are shown.Furthermore some quadrature formulae that exploit the bases properties are presented.In Section 4 applications of the wavelet bases introduced to kernel sparsification and image compression are shown.In particular in Subsection 4.1 we study some interesting properties of the bases considered and we present some results about integral kernel sparsification.In Subsection 4.2 we focus on applications of the wavelet bases to image coding and compression showing some interesting numerical results.Finally in the Appendix we give the wavelet mother functions necessary to construct the wavelet bases employed in the numerical experience presented in Section 4. To build these mother functions we have used the Symbolic Math Toolbox of MATLAB.The website http://www.econ.univpm.it/recchioni/scattering/w17contains auxiliary material and animations that help the understanding of the results presented in this paper and makes available to the interested users the software programs used to generate the wavelet mother functions of the wavelet bases used to produce the numerical results presented.A more general reference to the work of the authors and of their coworkers in acoustic and electromagnetic scattering where the wavelet bases introduced have been widely used is the website http://www.econ.univpm.it/recchioni/scattering.

Let us use the multiresolution analysis introduced by
Mallat [1,2], Meyer [3][4][5] and Coifman and Meyer [6] to construct wavelet bases.Let us begin introducing some notation.Let R be the set of real numbers, given a positive integer s let s R be the s-dimensional real Euclidean space, and let   be a generic vector where the superscript T means transposed.Let ( , )    and  denote the Euclidean scalar product and the corresponding Euclidean vector norm respectively.
For simplicity we restrict our analysis to the interval (0,1).More precisely, we choose be the Hilbert space of square integrable (with respect to Lebesgue measure) real functions defined on the interval (0,1).We want to construct orthonormal wavelet bases of via a multiresolution method similar to the methods used in [9,[1][2][3][4][5][6].Note that using the ideas suggested in [13,14] it is possible to generalize the work presented here to rather general domains A in dimension greater or equal to one.
Given an integer 1, M  we consider the following decomposition of where  denotes the direct sum of two orthogonal closed subspaces . In other words, the vector space generated by the union of ,1 L and we have: We take to be the space of the polynomials of degree smaller than M defined on (0,1) and we consider a basis of made of M polynomials orthogonal in the interval (0,1) with respect to the having 2 L -norm equal to one.For example we can choose as basis of the first M Legendre orthonormal polynomials defined on (0,1) and we refer to them as ( ), q L x (0,1), we use the ideas behind the multiresolution analysis.Let us begin defining the so called "wavelet mother functions".Let 2 N  be an integer and let   2, N  we define the following piecewise polynomial functions on (0,1): , where are polynomials of degree smaller than M to be determined.The functions , , , defined in (3) will be used as "wavelet mother functions".In fact through them we generate the elements of a wavelet family via the multiresolution analysis.
For simplicity let us choose = / , We note that results analogous to the ones obtained here with this choice can be derived for the general choice of , at the price of having more involved formulae.
Let us define the functions: whose supports are the intervals , = 1, 2, , , = 0,1, , = 0,1, , 1 , where are the functions defined above when we choose We want to choose J , M , N , and the coefficients of the polynomials that constitute the functions , , , in order to guarantee that the set . This can be done when J, M, N satisfy some constraints specified later imposing the following conditions to the wavelet mother functions (3): i) the functions , , , have the first M moments vanishing, that is: ii) the functions , , , normal functions, that is: , = 1, 2, , .
Note that depending on the choice of the integers J , M, N it may not be possible to satisfy the relations ( 6), (7) with functions of the form (3).
We note that in general the number of the unknown coefficients , , , , , , is bigger than the number of distinct equations contained in ( 6) and (7).More precisely only when = =1 J M and = 2 N the number of equations is equal to the number of unknowns and the unknown coefficients are determined up to a "sign permutation".That is we can change sign to the resulting wavelet mother functions.In this case the set of functions defined in (5), ( 6), (7)


is the well-known Haar's basis (see [8]).In all the remaining cases, when the relations (3), ( 6), (7) are compatible, the functions that satisfy (6), (7) generate through the multiresolution analysis an orthonormal set of . When the integers J, M, N satisfy some relations this orthonormal set is complete, that is is an orthonormal basis of , and can be regarded as a generalization of the Haar's basis.We must choose some criterion to determine the coefficients of the polynomials contained in (3) that remain undetermined after imposing (6), (7).It will be seen later that the criterion used to choose the undetermined coefficients influences greatly the "sparsifying" properties of the resulting wavelet basis, that is influences greatly the practical value of the resulting wavelet basis.On grounds of experience we restrict our analysis to two possible criteria: 1) impose some regularity properties to the wavelet mother functions (3), 2) require that the wavelet mother functions (3) have extra vanishing moments after those required in (6).
Note that the analysis that follows considers only some simple examples of use of these criteria and can be easily extended to more general situations to take care of several other meaningful criteria besides 1) and 2) such as, for example, a combination of these two criteria or ad hoc criteria dictated by special features of the concrete problems considered.
We begin adopting criterion 1).We choose = J   and in order to determine the coefficients left undetermined by ( 6), (7) we impose the following regularity conditions to the piecewise polynomial functions , , where the sets of indices are chosen such that ( 6), ( 7), (8) (and (3)) are compatible.
= , 1, l M M  in this case, with abuse of notation, we assume that no extra conditions on , , M N j N   are added to (6), (7).In other words when for some for , , is empty.We note that only some wavelet mother functions (not all of them) can have extra vanishing moments (i.e.only for some j we can have 6), ( 7), (9) (and (3) are incompatible.In fact when constitute a vector space of dimension NM .Imposing one extra vanishing moment to the functions , , , that is requiring ( 6), ( 7) and: is equivalent to choose 1 NM  linearly independent vectors in a vector space of dimension NM and this is impossible.However, for example, it is possible to impose one extra vanishing moment to  6), ( 7), (9) (and ( 3)) correspond in general to situations where the number of conditions is smaller than the number of coefficients of the piecewise polynomial functions, that is even after adding ( 9) to ( 6), (7) (and ( 3)) some of the coefficients of the wavelet mother functions may remain undetermined.
In this case the undetermined coefficients can be chosen arbitrarily or, for example, they can be chosen using criterion 1) or some other criterion.We denote with satisfying conditions ( 6), ( 7), (9) (and ( 3)) with a non trivial choice of 0, the functions defined in (4) when we replace , , we replace , , , , We will see later that the relations ( 3), ( 6), (7) are compatible and the corresponding set ( 5) is orthonormal but is not complete, moreover if the relations ( 3) ,( 6), ( 7) are incompaible.
Note that the fact that the wavelet mother functions (3) are piecewise polynomials from one hand makes easy and efficient their use in the most common computations (i.e. for example differentiation, integration,...), but on the other hand may create undesired effects in the discontinuity points of the wavelets when regular functions are approximated with discontinuous functions.
The choice of the criteria 1) and 2) is motivated by the following reasons.When criterion 1) is adopted we try to approximate regular functions using a basis made of "as much as possible" regular wavelets.This is done in order to minimize the undesired effects coming from the fact that regular functions are approximated with non regular wavelets choosing M large and N small.The goal that we pursue when we adopt criterion 2) is the construction of wavelet bases made of piecewise polynomial wavelets made of polynomials of low degree with "as much as possible" vanishing moments so that, as will be seen better in Section 4, the "sparsifying" properties of the resulting wavelet bases are improved in comparison with those of the wavelet bases that do not satisfy criterion 2).This is done choosing M small and N large so that a great number of extra vanishing moments can be imposed to many wavelet mother functions made of polynomials of low degree.As a matter of fact, choosing M and N appropriately, it is possible to construct wavelet bases satisfying to some extent the two criteria 1) and 2) simultaneously.However this is beyond the scope of this paper.We note that increasing N the number of the wavelet mother functions and the number of their discontinuity points increase.
In the Appendix we give the wavelet mother functions necessary to construct the wavelet bases used in the numerical experience presented in Section 4.

Some Mathematical Properties of the Wavelet Bases
Let us prove that the set of functions considered previously are orthonormal bases of Then we have: Proof.Properties ( 12), ( 13) can be easily derived from (11).The proof of ( 14) follows from the density of the piecewise constant functions in (see [17] Theorem 3.13 p. 84).This concludes the proof.

Note that for
2 N  be three integers such that the conditions ( 6), ( 7) can be satisfied with functions of the form (3) and let be the functions defined in (4), we have: Proof.Property (15) follows from definition (4) and Equation ( 6).The proof of ( 16 16) is obvious, when the supports are contained one into the other condition (16) follows from (15).Finally when = , m m = ,    = j j condition ( 16) follows from Equation (7).This concludes the proof.
Note that if we choose , , , , 15) can be improved, that is we can add to it the following condition: where 0 j   are the non negative integers (not all zero) that have been used in conditions ( 6), ( 7), ( 9) to determine the wavelet mother functions.Theorem 3.3.
the conditions ( 6), ( 7) can be satisfied with functions of the form (3) and the set it is easy to see that conditions (3), ( 6), ( 7) are compatible and the set is an orthonormal set of functions such that for = 0,1, , m  the subspace m T defined in (11) is contained in the subspace generated by . This concludes the proof.
The following corollary is a particular case of Theorem 3.3: (0,1) To unify the notation let us rename the functions belonging to the basis so that the basis , 0,1 defined in ( 5) can be rewritten as: = 0 w < 0, = 0,1, w 0, = 0,1, ( , where  ) ( denotes the maximum between 0 and  .
For later convenience in the study of integral operators and images let us remark that wavelet bases in  can be easily constructed as tensor product of wavelet bases of where , M , N N  are the quantities specified previously and we have chosen . The previous construction based on the tensor product can be easily extended to the case when A is a parallelepiped in dimension 3 s  .Note that with straightforward generalizations of the material presented here, it is easy to construct wavelet bases for   2 L A when A is a sufficiently simple subset of a real Euclidean space (see [12,16] to find several choices of A useful in some scattering problems).The analysis that follows of the wavelet bases when   = 0,1 A can be extended with no substantial changes to the other choices of A considered here.
The 2 L -approximation of a function with a wavelet expansion, is based on the calculation of the inner pro-ducts of the function to be approximated with the wavelets to find the coefficients that represent the function in the chosen wavelet basis.That is the function is approximated by a weighted sum of the wavelet basis functions.Let us note that in contrast with the Fourier basis that is localized only in frequency, the wavelet basis is localized both in frequency and in space, and in the most common circumstances only a few coefficients of the wavelet expansion must be calculated to obtain a satisfactory approximation.
In order to calculate the wavelet coefficients of a generic function, we proceed as follows.Given , be the following set of indices: The wavelet expansion of a function (19) can be written as follows: where for are given by:  18) and ( 4), we can rewrite (23) as follows: .
Note that the integrals in (26) depend on the function f and on the wavelet mother functions 1, , , given in (3).The change of variable = , , where that is given by: where , are the characteristic Lagrange polynomials defined as: and characterized by the property where in (27) we obtain the following approximation: , .
This technique leads to interpolatory integration rules with weights defined in terms of the wavelet mother functions.We note that in general, the weights of these quadrature rules can have alternating signs, this impacts negatively on the stability of the computations.Nevertheless, choosing n small it is possible to minimize the number of the quadrature weights large in absolute value and not all of the same sign and it is possible in everyday computing experience to avoid the Runge phenomenon.Having in mind definition (3) and choosing a small number of quadrature nodes, the integrals (33) that define , , , , are very easy to calculate since the integrands are low degree piecewise polynomial functions.
It is easy to see that (34) is valid also for , , , In conclusion, once the wavelet basis and the nodes , (with n small) have been chosen, the integrals that define , , , , that is the integrals (33) and (35), can be calculated explicitly and tabulated so that the coefficients , ,  defined in ( 23) can be approximated with the quantities defined in (34) and therefore the wavelet expansion ( 22) can be approximated very efficiently.
Let us define the quadrature error made approximating with (34) the wavelet coefficients given in (26), that is: , .
For the quadrature error we prove the following lemma: Lemma 3.5.Let and let t be a point belonging to the domain of the function is the space of the real continuous functions defined on (0,1) 1 n  -times continuously differentiable.Then there exists a continuous function such that the quadrature error is given by: where Proof.It is sufficient to note that from ( 29) and (31) we have: The thesis follows from standard numerical analysis results, see for example [18] p. 49.This concludes the proof.
Note that a result similar to (37) is valid also for It is worthwhile to note that usually the convergence properties of general quadrature rules, such as the one presented here, depend on the smoothness of the integrand (i.e.discontinuities of the integrand or of any of its derivatives may disturb the convergence properties of the quadrature rules), when non smooth functions, such as, for example, piecewise smooth functions are involved in integrals like (23), it is convenient to split the integration interval (0,1) into subintervals corresponding to the smooth parts of the integrands and to apply a low order quadrature rule on each subinterval.Ad hoc quadrature rules for the computation of wavelet coefficients have been developed by several authors, see for example [19][20][21][22][23].It could be interesting to extend the work of these authors to the wavelet bases presented here.
The explicit computation of the integrals , , , , has been done easily with the help of the Symbolic Math Toolbox of MATLAB.More in detail, in the numerical experiments of the next section we use coefficients , , , , n  calculated with composite quadrature formulae choosing in each interval = 0 n and the node 0 t given by the middle point of each subinterval.This choice corresponds to a one-point quadrature formula in each subinterval for the evaluation of the wavelet coefficients, is motivated by the fact that in Subsection 4.2 we manipulate digital images and, due to the usual pixel structure of the digital images, a digital image can be regarded as a piecewise constant real function defined on a rectangular region of the two-dimensional Euclidean space 2 R .So that if we consider bidimensional composite quadrature formulae with as many intervals as the pixels of the image in each dimension, with the intervals coinciding with the pixel edges, and in each interval we calculate (33) and (35) with the choice = 0, n at the price of very simple calculations exact wavelet coefficients can be obtained.Moreover these quadrature formulae turn out to be useful for the rapid evaluation of the wavelet coefficients in the approximation of integral kernels.
It is worthwhile to note that composite quadrature formulae with > 0 n can be obtained with no substantial modifications of the procedure described above.

Wavelet Bases, Decay Estimates and Kernel Sparsification
We present some interesting properties of the wavelet bases introduced in Section 2. In particular we show how the representation in these wavelet bases of certain classes of linear operators acting on may be viewed as a first step for their "compression", that is as a step to approximate them with sparse matrices.After being compressed these operators can be applied to arbitrary functions belonging to in a "fast" manner and linear equations involving these operators can be approximated in high dimensional spaces with sparse systems of linear equations and solved at an affordable computational cost.In particular we show how the orthogonality of the wavelet functions to the polynomials up to a given degree (the vanishing moments property) plays a crucial role in producing these sparse approximations.
then there exists a positive constant M D such that we have: = 0,1, , = 0,1, , = 0,1, , Proof.The proof is analogous to the proof of Proposition 4.1 in [19].In fact from (39) it follows that there exists a positive constant M C such that: , m m be as above and   * * , x y be the center of mass of the set For the wavelet basis function having extra vanishing moments, the previous theorem can be improved.That is be as above, we have: = 0 where the constants , are those appearing in (9) and are such that the Equations ( 6), ( 7), (9) (and ( 3)  be a real function such that: be the set of functions defined in Section 2 and  , , ., , , be the following quantities: , , , ,  , , , ,  0 1 , , , , 0 = , , = 0,1, , = 0,1, , = 0,1, , 1, , , Proof.The proof follows the lines of the proof of Theorem 4.1.Going into details, condition (43) implies that there exists a positive constant * M E such that: it is convenient to proceed as done previously, when ( 1) with base point The estimates (41) and ( 45) are the basic properties that together with a simple truncation procedure make the wavelet bases introduced previously useful to approximate with sparse matrices linear operators.Let  be an integral operator acting on where the kernel K is a real function of two variables.
The compression algorithms we are interested in are based on the following observation.Generalizing what done in Section 2 for a function , when  is represented on a wavelet basis, i.e. when the kernel K is expanded as a function of two variables on the wavelet basis (20), the calculation of the entries of the (infinite) matrix that represents the operator in the wavelet basis involves the evaluation of the integrals (40) or (44).If the wavelet basis considered has several vanishing moments and we look at truncated wavelet expansions, that is to an expansion on a truncated wavelet basis, when the kernel K is sufficiently regular, thanks to the estimates (41) or (45), the fraction of the entries whose absolute value is smaller than a chosen threshold approaches one when the truncated expansion approaches the true expansion.The entries whose absolute value is smaller than a suitable threshold can be set to zero committing a "small" error.This property makes possible the approximation of the integral operator (47) with sparse matrices and makes possible to solve numerically integral equations very efficiently.
Note that when two different wavelet bases with the same number of vanishing moments are used to represent the operator  the estimates (41) and (45) show that the wavelets with the smaller between the constants M D or , M M j j F  will have expansion coefficients that decay faster to zero.
Similar arguments can be made in the discrete case, where we consider a piecewise constant function defined on a bounded rectangular domain of the form: , =1,2 , , , 0,1 0,1 , where r is a positive integer and > 0 h is such that = 1 h r  and where , Some applications of these estimates to specific examples are shown below.We use the Wavelet Bases 1, 2, 3, 4 of the Appendix to test the performances of the integral operators compression algorithm described previously.We apply the algorithm to a number of operators.The direct and inverse wavelet transforms used in the numerical experience associated to these wavelet bases have been implemented in FORTRAN 90.Briefly we remind that the Wavelet Bases 1 and 3 of the Appendix are done with piecewise polynomial functions made of polynomials of degree zero, while the Wavelet Bases 2 and 4 of the Appendix are done with piecewise polynomial functions made of polynomials of degree one.Moreover the Wavelet Bases 1 and 2 are obtained using the regularity criterion (i.e.criterion 1)), while the Wavelet Bases 3 and 4 are obtained using the "extra vanishing moments" criterion (i.e.criterion 2)).
The choice of reporting here and in Subsection 4.2 the numerical results obtained with the Wavelet Bases 1, 2, 3, 4 of the Appendix is motivated by the fact that these are the simplest bases among the bases introduced in Section 2 that it is possible to construct and manipulate.Moreover the Bases 1 and 3 made of piecewise constant functions are particularly well suited to deal with piecewise constant functions.This type of functions is very important since it can be identified with digital signals and images.
In particular we represent a linear operator in a wavelet basis, compress it by setting entries of its representation on the wavelet basis whose absolute value is below a chosen threshold to zero, and reconstruct the "compressed" operator using the inverse wavelet transform.The comparison between the original and the reconstructed "compressed" operator is made evaluating the relative 2  L -difference between the original and the reconstructed compressed kernels and is very satisfactory.Moreover comparing our results with those obtained by Beylkin, Coifman and Rokhlin in [19] and by Keinert in [24], we observe that the wavelet bases studied here show similar or better compression properties of the wavelet bases used in [19,24] even when we use a smaller number of vanishing moments than that used in [19] and in [24].The results of two experiments are presented in this section and are summarized in Tables 1-3 and illustrated in Figures 1-4.In these tables dim indicates the dimension of the vector space generated by the truncated wavelet basis and comp C  is the compre-ssion coefficient obtained with our algorithm when the truncation threshold > 0  is used.The compression coefficient comp C  is defined as the ratio between 2 dim and the number of matrix elements whose absolute value exceeds a threshold > 0  .Example 1.In this example we consider the kernel: and we expand it in the Wavelet Bases 1, 2, 3, 4 of the Appendix.We set to zero the entries of the resulting matrix whose absolute value is smaller than a threshold  and we obtain the results shown in Tables 1 and 2. In particular in Table 1 we show the compression coefficients comp C  for three different values of the threshold  and in Table 2 we report the relative error in Note that in Tables 1, 2 and in the tables that follow, the symbol \ indicates that using all the elements of the specified basis up to a certain "resolution level" m it is impossible to reach the dimension dim specified in the corresponding second column of the tables.
The structures of the non-zero entries after the truncation procedure in the matrices obtained using the Wavelet Bases 2 and 4 of the Appendix when    Example 2. In this example we compress the following piecewise constant function: , , , = 1, 2, , , , 0,1 0,1 .
where , i j A is the "pixel" of index , , i j defined previously.Remind that = 1 h dim  .We use the Wavelet Bases 1, 2, 3, 4 of the Appendix and three different values of the threshold  .Table 3 and Figures 3 and 4 describe the results of these experiments.The entries above the threshold From Tables 1-3 and Figures 1-4 the following observations can be made.As far as the compression property is concerned we have a really impressive improvement going from the Wavelet Basis 1 to the Wavelet Basis 3 of the Appendix, that is there is a real advantage in using the "extra vanishing moments" criterion.This fact is more evident in Example 1 where a continuous kernel is considered (in this case the compression coefficient is approximately multiplied by two).There is not the same improvement going from the Wavelet Basis 2 to the Wavelet Basis 4 of the Appendix, however using these two bases we obtain much higher compression coeffi-cients than those obtained with the Wavelet Bases 1 and 3 of the Appendix.Moreover the 2 L -relative errors obtained comparing the original and the reconstructed "compressed" operators are always at least one order of magnitude smaller when the Wavelet Bases 2, 4 of the Appendix are used instead than the Wavelet Bases 1, 3 of the Appendix.

Image Compression and Wavelets
Let us present some results about digital image compression and reconstruction.
With the growth of technology in the last decades and the entrance into the so-called "Digital-Age", a vast amount of digital information has become available for storage and/or exchange, and the efficient treatment of this enormous amount of data often present difficulties.Wavelet analysis is one way to deal with this problem.It produces several important benefits, particularly in image compression.Compression is a way of encoding data more concisely or efficiently.It relies on two main strategies: getting rid of redundant information ("redundancy reduction") and getting rid of irrelevant information ("irrelevancy reduction").Redundancy reduction concentrates on more efficient ways of encoding the image.It looks for patterns and repetitions that can be expressed more efficiently.Irrelevancy reduction aims to remove or alter information without compromising the quality of the image.These two strategies conduct to two kinds of compression schemes: lossless and lossy ones.Lossless compression is generally based on redundancy reduction and the key point is that no information is irreversibly lost in the process.Once decompressed, a lossless compressed image will always appear exactly the same as the original, uncompressed, image.Lossy compression is based on both irrelevancy and redundancy reductions.It transforms and simplifies the image in a way that a much greater compression than the com-pression obtained with lossless approaches can be achieved, but this process is by definition irreversible, that is it permanently loses information.The lossy compressions are suitable for situations where size is more crucial than quality, such as downloading via Internet.The JPEG compression is the best known lossy compression methodology, and the JPEG compression is based on the use of wavelets [25].
From the theoretical point of view wavelet compression is capable of both lossless and lossy compression.In fact the wavelet transform furnishes simplified versions of the image together with all the information necessary to reconstruct the complete original image.All the information can be kept and encoded as a lossless compressed image.Alternatively, only the most significant details can be added back into the image producing a simplified version of the image.As you might expect, this second application has a much larger area of interest.
We consider here some applications of the wavelet bases constructed in the previous sections in image compression.In particular we show some reconstructions of images and we focus on a lossy compression procedure.We limit our attention to grayscales images.Despite the appearance, compressing grayscale images is more difficult than compressing colour images.In fact human visual perception can distinguish more easily brightness (luminance) than colour (chrominance).
Going into details the key steps of our image compression and reconstruction scheme are the following: 1) digitize the original image, 2) apply the wavelet transform to represent the image through a set of wavelet coefficients, 3) only for lossy compression: manipulate (i.e.set to zero) some of the wavelet coefficients, 4) reconstruct the image with the inverse wavelet transform.
The first step, that is the digitization of the image, consists in transforming the image in a matrix of numbers.Since we consider black and white images, the digitized image can be characterized by its intensity levels, or scales of gray which range from 0 (black) to 255 (white), and its resolution, that is the number of pixels per square inch.The second step is to decompose the digitized image into a sequence of wavelet coefficients.The lossy compression step is based on the fact that many of the wavelet coefficients are small in absolute value, so that they can be set to zero with little damage to the image.This procedure is called thresholding.In practice the most simple thresholding procedure consists in choosing a fixed tolerance and in doing the following truncation procedure: the wavelet coefficients whose absolute value falls below the tolerance are put to zero.The goal is to introduce as many zeros as possible without losing a great amount of details.The crucial issue consists in the choice of the threshold.There is not a straightforward way to choose it, although the larger the threshold is chosen, the larger is the error introduced into the process.In the lossy compression scheme only the wavelet coefficients that are non zero after the thresholding procedure are used to build the so called compressed image that may be stored and transferred electronically using much less storage space than the space needed to store the original image.Obviously this type of compression is a lossy compression since it introduces error into the process.If all the wavelet coefficients are used in the inverse wavelet transform (or equivalently, if we take the threshold equal to zero) and an exact arithmetic is used the original image can be reconstructed exactly.
Mimiking what done in Subsection 4.1, we compress each test image by taking its wavelet matrix representation and setting to zero the entries of matrix wavelet representation whose absolute value is smaller than a fixed threshold.After this truncation procedure we perform an inverse wavelet transform on the resulting "compressed" matrix representation and we compare the resulting image with the original image both from the qualitative and the quantitative point of view.As done in Subsection 4.1 we compute the resulting compression coefficient comp C  as a function of the truncation threshold  and we show how the vanishing moments property plays a crucial role in producing compressed images.
Let us note that sometimes in the scientific literature the compression coefficient is not defined as done in Subsection 4.1 but it is defined dividing the original number of bytes needed to represent the image by the number of bytes needed to store the compressed image.However, for example, Wikipedia defines the compression ratio as the size of the compressed image compared to that of the uncompressed image leaving undetermined how to measure the size of an image.In [26] DeVore, Jawerth and Lucier have shown that between these two definitions of compression coefficient (i.e.number of non zero coefficients and number of bytes needed to represent these coefficients) there is a very high "correlation" (i.e.0.998).The compression coefficient (ratio) is one of the common performance indices of image compression level.
Example 3. In this example we consider the image of      Example 4. We consider the image of Figure 9.In this image are evident three different types of objects: a shape (the black star), a number ('2025') and an inscription ('black').This image has 278  268 pixels.Table 5 shows the compression coefficients      compression.In fact in this last case the use of the Wavelet Bases 1 and 3 (made of piecewise polynomial functions made of polynomials of degree zero and obtained with the two different criteria proposed in Section 2) and the use of Wavelet Bases 2 and 4 (made of piecewise polynomial functions made of polynomials of degree one and obtained with the two different criteria proposed in Section 2) furnish similar compression coefficients.This is not really surprising since the images have a natural discontinuity structure given by the pixels.Nevertheless the compression coefficients obtained with the Wavelet Basis 3 are always higher than those obtained using the Wavelet Basis 1 and the Wavelet Bases 2 and 4 give compression coefficients higher than those obtained with the Wavelet Bases 1 and 3. Furthermore the use of the Wavelet Bases 2 and 4 seems to have a regularizing effect on the images, that is the images reconstructed with these two bases appear to have a smaller number of edges and less contrast than those obtained using the Wavelet Bases 1 and 3.
As far as the reconstruction of the images is concerned, we can say that very satisfactory reconstructions are obtained when the dimension of the vector space generated by the truncated wavelet basis raised to the power two is about the same than the number of pixels of the image considered.Furthermore the relative 2  L -errors made substituting the original image with the reconstructed image have approximately the same order of magnitude independently from the basis used.
Finally the following conclusions can be made.In order to manipulate correctly operators and images it is sufficient to construct the wavelet bases with piecewise polynomial functions made of polynomials of very low degree.Really for the images seems to be adequate choosing piecewise polynomial functions made of polynomials of degree zero.As already observed this might be due to the way we calculate the wavelet coefficients and to the fact that the operators and the images are represented by piecewise constant functions.Moreover it seems to be very promising the idea of increasing the number of vanishing moments keeping low the degree of the polynomials used.Actually the bases that have a large number of extra vanishing moments, that is those constructed with the second criterion proposed in Section 2, show better compression and reconstruction properties, and in general work better than the wavelet bases constructed with the first criterion proposed in Section 2.

Appendix: Symbolic Calculus and Some Wavelet Bases
The wavelet mother functions used in the numerical experiments of Section 4 have been obtained using the Symbolic Math Toolbox of MATLAB to implement the procedure described in Section 2. In particular some of the input parameters used to determine the wavelet mother functions are: the minimum number M of vanishing moments that the wavelet basis considered must have, the number N of subintervals of the interval (0,1) employed and the points = , where the subdivision of (0,1) in subintervals takes place.The corresponding symbolic non linear system arising from Equations ( 3), ( 6), ( 7) and ( 8), or from Equations ( 3), ( 6), ( 7) and ( 9), having as unknowns the coefficients of the wavelet mother functions, has been solved with the command solve of MATLAB.Let us note that all the wavelet bases we present are uniquely determined up to a "sign permutation".
Below we show the mother functions of the wavelet bases used in the numerical experience presented in Section 4. We begin showing some functions  , , solution of ( 3), ( 6), ( 7) and ( 8), that is some wavelet basis functions obtained using as extra condition criterion 1), the regularity criterion.that are instead solution of (3), ( 6), ( 7) and (9), that is they are obtained using as extra condition criterion 2), the "extra vanishing moments" criterion.We underline that taking the number N of subin-tervals of the interval (0,1) used big enough, wavelet mother functions of piecewise polynomial functions made of polynomials of fixed degree with an arbitrary high number of vanishing moments can be constructed.
) when = , m m =    and j j  follows from the fact that the supports of sets or sets contained one into the other.When = m m and     the supports are disjoint and when , m m  let us suppose for example > , m m the supports are either disjoint sets or sets contained one into the other depending on the values of the indices  and   .In fact when the supports are disjoint condition (

2 L
-inner product of the function f either with a wavelet function, or with a Legendre polynomial depending on the values of the indices.In order to calculate efficiently the integrals(23) we construct some simple ad hoc interpolatory quadrature formulae that take into account the definition of the basis functions.In particular let y  and base point * = y y when < m m or the Taylor polynomial of degree 1 M  of   ( ) = ( , ), 0,1 , x K x y x Because of Equation(17) and assumption (46) mimicking the proof of Theorem 4.1, we can obtain the estimate (45) for  , , , , , * =x x obtaining with a similar procedure the estimate (45) for  , , , , , concludes the proof.

2 L
-norm made substituting the original operator with the reconstructed "compressed" operator after the trun- in black in Figures 1 and 2 respectively.In particular in Figures 1 and 2 the matrices shown have = 512, dim that is they are matrices of 512 rows and columns.

Figure 1 .
Figure 1.Example 1: Sparsity structure of the matrix obtained with the Wavelet Basis 2 of the Appendix for dim = 512 The entries above the threshold -6 = 10 ε in absolute value are shown in black.

Figure 2 .
Figure 2. Example 1: Sparsity structure of the matrix obtained with the Wavelet Basis 4 of the Appendix for 512 dim = .The entries above the threshold -6 = 10 ε in absolute value are shown in black.
solute value of the matrices relative to the matrix (50) obtained using the Wavelet Bases 1 and 3 of the Appendix when = 256 dim are shown in black in Figures 3 and 4 respectively.

Figure 3 . 7 = 10 ε
Figure 3. Example 2: Sparsity structure of the matrix obtained with the Wavelet Basis 1 of the Appendix for 256 dim = .The entries above the threshold -7 = 10 ε in absolute value are shown in black.

Figure 4 .
Figure 4. Example 2: Sparsity structure of the matrix obtained with the Wavelet Basis 3 of the Appendix for 256 dim = .The entries above the threshold -7 = 10 ε in absolute value are shown in black.

Figure 5
which is the famous Lena image.This image has 256  256 pixels.In Table 4 we report the compression coefficients comp C  obtained using the Wavelet Bases 1, 2, 3, 4 of the Appendix and two different values of the threshold  .The relative 2 L -errors made substituting the original image with the compressed image range between 3 10  and 2 10  .Figures 6-8 show the compressed Lena images obtained with the Wavelet Basis 3 of the Appendix for

2 L
Bases 1, 2, 3, 4 of the Appendix and two different values of the threshold  .The relative -errors made substituting the original image with the compressed image respectively.The comparisons between Tables 4-5 and Figures 5-15 suggest the following observations and comments.With respect to the compression property there is a different behaviour of the Wavelet Bases 1, 2, 3, 4 of the Appendix between operator compression and image
2 is one of the multi-wavelets bases introduced by Alpert in[9].Let us show now some functions 