Wavelet-based multifractal analysis

Themultifractal formalism was introduced in the context of fully-developedturbulencedata analysis and modeling to account for the experimental observation of some deviationtoKolmogorov theory (K41) of homogenous and isotropic turbulence (Frisch, 1995). Thepredictions of various multiplicative cascade models, including the weighted curdling (binomial)model proposed by Mandelbrot (1974), were tested using box-counting (BC) estimatesof the so-called \(f(\alpha)\) singularity spectrum of the dissipation field (Meneveau & Sreenivasan,1991). Alternatively, the intermittent nature of the velocityfluctuations were investigatedvia the computation of the \(D(h)\) singularity spectrum using the structure function (SF)method (Parisi & Frisch, 1985). Unfortunately, both types of studies suffered from severe insufficiencies.On the one hand, they were mostly limited by one point probe measurements tothe analysis of one (longitudinal) velocity component and to some 1D surrogate approximationof the dissipation (Aurell et al., 1992). On the other hand, both the BC and SF methodologieshave intrinsic limitations and fail to fully characterize the corresponding singularity spectrumsince only the strongest singularities are a priori amenable to these techniques (Arneodo et al.,1995b; Bacry et al., 1993; Muzy et al., 1993, 1994). In the early nineties, a statistical approachbased on thecontinuous wavelet transform was proposed as a unified multifractal descriptionof singular measures and multi-affine functions (Arneodo et al., 1995b; Bacry et al., 1993;Muzy et al., 1993, 1994). Applications of the so-calledwavelet transform modulus maxima(WTMM) method have already provided insight into a wide variety of problems, e.g., fullydeveloped turbulence, econophysics, meteorology, physiology andDNA sequences (Arneodoet al., 2002). Let us note that alternative approaches to the multifractal description havebeen developed using discretewavelet bases (Abry et al., 2000, 2002a,b; Veitch & Abry, 1999)including the recent use of wavelet leaders (Jaffard et al., 2006; Wendt & Abry, 2007; Wendtet al., 2007). Later on, the WTMM method was generalized to 2D for multifractal analysisof rough surfaces (Arneodo et al., 2000; Decoster et al., 2000), with very promising resultsin the context of the geophysical study of the intermittent nature of satellite images of thecloud structure (Arneodo et al., 1999a, 2003; Roux et al., 2000) and the medical assist in thediagnosis in digitized mammograms (Arneodo et al., 2003; Kestener et al., 2001). Recently,the WTMM method has been further extended to 3D scalar as well as 3Dvector field analysisand applied to 3D numerical data issue from isotropic turbulence direct numerical simulations(DNS) (Kestener & Arneodo, 2003, 2004, 2007).

The continuous wavelet transform

Introduction

The continuous wavelet transform (WT) is a mathematical technique introduced in signalanalysis in the early 1980s (Goupillaud et al., 1984; Grossmann & Morlet, 1984).Since then, it has been the subject of considerable theoretical developmentsand practical applications in a wide variety of fields.The WT has been early recognized as a mathematical microscope thatis well adapted to reveal the hierarchy that governs the spatialdistribution of singularities of multifractalmeasures (Arneodo et al., 1988, 1989, 1992).What makes the WT of fundamental use in the present study is that itssingularity scanning ability equally applies to singular functionsthan to singular measures (Arneodo et al., 1988, 1989, 1992; Holschneider, 1988; Holschneider & Tchamitchian, 1990; Jaffard, 1989, 1991; Mallat & Hwang, 1992; Mallat& Zhong, 1992).This has led Alain Arneodo and his collaborators (Arneodo et al., 1995b; Bacry et al., 1993; Muzy et al., 1991, 1993, 1994)to elaborate a unified thermodynamic description of multifractaldistributions including measures and functions, the so-called WaveletTransform Modulus Maxima (WTMM) method.By using wavelets instead of boxes, one can take advantage of thefreedom in the choice of these "generalized oscillating boxes"to get rid of possible (smooth) polynomial behavior that mighteither mask singularities or perturb the estimation of theirstrength \(h\) (Hölder exponent), remedying in this way for oneof the main failures of the classical multifractal methods(e.g. the box-countingalgorithms in the case of measures and thestructure function method in the case of functions (Arneodo et al., 1995b; Bacry et al., 1993; Muzy et al., 1993, 1994)).The other fundamental advantage of using wavelets is thatthe skeleton defined by the WTMM (Mallat & Hwang, 1992; Mallat & Zhong, 1992),provides an adaptative space-scale partitioning from which onecan extract the \(D(h)\) singularity spectrum via the Legendretransform of the scaling exponents \(\tau(q)\) (\(q\) real, positiveas well as negative) of some partition functions defined from theWT skeleton.We refer the reader to Bacry et al. (1993), Jaffard (1997a,b) for rigorous mathematical resultsand to Hentschel (1994) for the theoretical treatmentof random multifractal functions.

Definition

Figure 1: Set of analyzing wavelets \(\psi(x)\) that can be used in Eq. (1). Subfigures (a) \(g^{(1)}\) and (b) \(g^{(2)}\) as defined in Eq. (3).

The WT is a space-scale analysis which consists inexpanding signals in terms ofwavelets which areconstructed from a single function, theanalyzing wavelet\(\psi\ ,\) by means of translations and dilations.The WT of a real-valued function \(f\) is definedas (Goupillaud et al., 1984; Grossmann& Morlet, 1984)\[\tag{1}T_\psi[f](x_0,a)= \frac{1}{a}\int_{-\infty}^{+\infty} f(x)\psi(\frac{x-x_0}{a})dx \; ,\]

where \(x_0\) is the space parameter and \(a\) (\(>0\)) the scale parameter.The analyzing wavelet \(\psi\) is generallychosen to be well localized in both space and frequency.Usually \(\psi\) is required to be of zero mean for the WTto be invertible. But for the particular purposeof singularity tracking that is of interest here, we will furtherrequire \(\psi\) to be orthogonal to low-orderpolynomials (Arneodo et al., 1995b; Bacry et al., 1993; Holschneider & Tchamitchian, 1990; Jaffard, 1989, 1991; Mallat & Hwang, 1992; Mallat & Zhong, 1992; Muzy et al., 1991, 1993, 1994):

\[\tag{2}\int_{-\infty}^{+\infty}x^m \psi(x) dx = 0 \;, \;\; \;0 \leq m < n_\psi\; .\]

As originally pointed out by Mallat and collaborators (Mallat & Hwang, 1992; Mallat & Zhong,1992), for the specific purpose of analyzing the regularity of a function, one can get rid of theredundancy of the WT by concentrating on the WT skeleton defined by its modulus maxima only.These maxima are defined, at each scale \(a\ ,\) as the local maxima of \(|T_\psi[f](x,a)|\) considered as a function of \(x\ .\)As illustrated in Figure2(e,f),these WTMM are disposed on connected curves in the space-scale (or time-scale) half-plane,calledmaxima lines.Let us define \(\mathcal{L}(a_0)\) as the setof all the maxima lines that exist at the scale \(a_0\)and which contain maxima at any scale \(a \leq a_0\ .\)An important feature of these maxima lines, when analyzing singular functions, is that there is at leastone maxima line pointing towards each singularity (Mallat & Hwang, 1992; Mallat & Zhong,1992; Muzy et al., 1994).

Analyzing wavelets

There are almost as many analyzing wavelets as applications of the continuous WT (Arneodoet al., 1988, 1989, 1992, 1995b; Bacry et al., 1993; Muzy et al., 1991, 1993, 1994). A commonlyused class of analyzing wavelets is defined by the successive derivatives of the Gaussian function:

\[\tag{3}g^{(N)}(x)=\frac{d^N}{dx^N}e^{-x^2/2} \; ,\]

for which \(n_\psi=N\) and more specifically \(g^{(1)}\) and \(g^{(2)}\) that areillustrated in Figure1(a,b).Note that the WT of a signal \(f\) with \(g^{(N)}\) (Eq. (3)) takes the following simple expression:

\[\tag{4}\begin{array}{lcl} T_{g^{(N)}}[f](x,a) &=& \frac{1}{a} \int_{-\infty}^{+\infty} f(y) g^{(N)} \left( \frac{y-x}{a} \right) dy, \\ &=& a^N \frac{d^N}{dx^N}T_{g^{(0)}}[f](x,a).\end{array}\]

Equation (4) shows that the WT computed with \(g^{(N)}\) at scale \(a\) is nothing but the Nthderivative of the signal \(f(x)\) smoothed by a dilated version \(g^{(0)}(x/a)\) of theGaussian function. This property is at the heart of various applications of the WT microscopeas a very efficient multi-scale singularity tracking technique (Arneodo et al., 2002).

Scanning singularities with the wavelet transform modulus maxima

Figure 2: WT of monofractal and multifractal functions.
DNA walk:
(a) DNA walk \(f(x)\) of the largest intron of human dystrophin gene (\(L=109574\)) using the “G” mononucleotide coding
(c) WT of DNA walk landscape color coded, independently at each scale \(a\ ,\) using 256 colors from black (\(|T_\psi|=0\)) to red (\(\max|T_\psi|\))
(e) WT skeleton defined by the set of all maxima lines
Fully developed turbulence:
(b) Longitudinal velocity signal recorded in the Modane wind-tunnel experiment (\(R_\lambda\simeq2000\)) over about two integral scales
(d) WT of the velocity signal in (b) using the same color coding as in (c)
(f) corresponding WT skeleton.
The analyzing wavelet is the Mexican hat \(g^{(2)}\) (Eq. (3)).

The strength of the singularity of a function \(f\) at point \(x_0\) is givenby theHölder exponent,i.e., the largest exponent such that thereexists a polynomial \(P_n(x-x_0)\) of order \(n < h(x_0)\) anda constant \(C > 0\ ,\) so that for any point \(x\) in a neighborhood of \(x_0\ ,\) one has(Bacry et al., 1993; Holschneider & Tchamitchian, 1990; Jaffard, 1989, 1991; Mallat & Hwang, 1992; Mallat & Zhong,1992; Muzy et al., 1994)\[\tag{5}|f(x) - P_n(x-x_0)| \leq C|x-x_0|^h.\]

If \(f\) is \(n\) times continuously differentiable at the point\(x_0\ ,\) then one can use for the polynomial \(P_n(x-x_0)\ ,\)the order-\(n\) Taylor series of \(f\) at \(x_0\) and thus provethat \(h(x_0) > n\ .\) Thus \(h(x_0)\) measures how irregularthe function \(f\) is at the point \(x_0\ .\) The higher theexponent \(h(x_0)\ ,\) the more regular the function \(f\ .\)

The main interest in using the WT for analyzingthe regularity of a function lies in its ability to beblind to polynomial behavior by an appropriate choiceof the analyzing wavelet \(\psi\ .\) Indeed, let usassume that according to Eq.(5),\(f\) has, at the point \(x_0\ ,\) a local scaling (Hölder)exponent \(h(x_0)\ ;\) then, assuming that the singularity is notoscillating (Arneodo et al., 1997a, 1998b; Mallat& Zhong, 1992),one can easily prove thatthe local behavior of \(f\) is mirrored by the WTwhich locally behaves like (Arneodo et al., 1995b; Bacry et al., 1993; Holschneider & Tchamitchian, 1990;Jaffard, 1989, 1991, 1997a,b; Mallat & Hwang, 1992; Mallat & Zhong, 1992; Muzy et al., 1991,1993, 1994)\[\tag{6}T_\psi[f](x_0,a) \sim a^{h(x_0)}\; , a \rightarrow 0^+\; ,\]

provided \(n_\psi>h(x_0)\ ,\) where \(n_\psi\) is the number of vanishing moments of \(\psi\) (Eq.(2)).Therefore one can extract the exponent\(h(x_0)\) as the slope of a log-log plot of the WT amplitude versus thescale \(a\ .\) On the contrary, if one chooses \(n_\psi < h(x_0)\ ,\) the WT still behavesas apower-law but with a scaling exponent which is \(n_\psi\ :\)\(\tag{7}T_\psi[f](x_0,a) \sim a^{n_\psi}\; , a \rightarrow 0^+\; .\)

Thus, around a given point \(x_0\ ,\) the faster the WTdecreases when the scale goes to zero, the more regular \(f\) isaround that point. In particular,if \(f\in C^\infty\) at \(x_0\) (\(h(x_0) = +\infty\)), then the WT scalingexponent is given by \(n_\psi\ ,\)i.e. a value which is dependent on the shape of the analyzing wavelet.According to this observation, one can hope to detect the points where\(f\) is smooth by just checkingthe scaling behavior of the WT when increasing the order \(n_\psi\)of the analyzing wavelet (Arneodo et al., 1995b; Bacry et al., 1993; Muzy et al., 1991, 1993, 1994).

Remark

A very important point (at least for practical purpose) raisedby Mallat and Hwang (Mallat& Hwang, 1992) is that the local scaling exponent \(h(x_0)\)can be equally estimated by looking at the valueof the WT modulus along a maxima line converging towards the point \(x_0\ .\)Indeed one can prove that both Eqs.(6) and (7) still holdwhen following a maxima line from largedown to small scales (Mallat & Hwang, 1992; Mallat & Zhong, 1992).

The wavelet transform modulus maxima method for multifractal analysis

Singularity spectrum

As originally defined by Parisi & Frisch (1985), the multifractal formalism of multi-affine functionsamounts to compute the so-calledsingularity spectrum \(D(h)\) defined as theHausdorff dimension of the set where the Hölder exponent is equal to\(h\) (Arneodo et al., 1995b; Bacry et al., 1993; Muzyet al., 1994):

\[\tag{8}D(h)=dim_H\{ x\;,\; h(x)=h \} \; ,\]

where \(h\) can take,a priori, positive as well asnegative real values (e.g., the Dirac distribution\(\delta(x)\) corresponds to the Hölder exponent \(h(0)=-1\)) (Jaffard, 1997a).

The WTMM method

Figure 3: Determination of the \(\tau(q)\) and \(D(h)\) multifractal spectra of the “G” DNA walks of a set of 2184 human introns (\(L\geq800\)bp) (\({\color{Red}\circ}\)) and of the turbulent velocity signals (\(\color{Blue}\bullet\)) using the WTMM method. \(\log_2Z(q,a)\) vs \(\log_2a\) for, from bottom to top, \(q=\)-2, -1, 0, 1, 2, 3 and 4 in (a) and for \(q=\)-3, 0, 3 and 6 in (b). (c) \(\tau(q)\) vs \(q\ ;\) the solid lines correspond respectively to a monofractal spectrum \(\tau(q)=qH-1\) with \(H=0.60\pm0.02\) (red) and to a quadratic multifractal spectrum \(\tau(q)=c_1q-\frac{c_2}{2}q^2-1\) with \(c_1=0.36\pm0.02\) and a non zero value for the intermittent coefficient \(c_2=0.028\pm0.003\) (blue). (d) \(D(h)\) vs \(h\ ;\) while the singularity spectrum of the monofractal DNA walk landscape reduces to a point (\(D(h=0.6)=1\) and 0 elsewhere, \(\color{Red}\circ\)), it has a parabolic shape and extends from \(h_{min}=0.12\) to \(h_{max}=0.60\) as the signature of the “intermittent” multifractal nature of Eulerian turbulence (\(\color{Blue}\bullet\)). The analyzing wavelet is \(g^{(2)}\) (Eq. (3)).

A natural way of performing a multifractalanalysis offractal functions consists in generalizing the “classical”multifractal formalism (Collet et al., 1987; Grassberger et al., 1988; Halsey et al., 1986;Paladin & Vulpiani, 1987; Rand, 1989)using wavelets instead of boxes.By taking advantage of the freedom in the choice of the“generalized oscillating boxes” that arethe wavelets, one can hope to get rid of possible smooth behaviorthat could mask singularities or perturb the estimationof their strength \(h\ .\) But the major difficulty withrespect to box-counting techniques (Argoul et al.,1990; Farmer et al., 1983; Grassberger & Procaccia, 1983; Grassberger et al., 1988; Meneveau& Sreenivasan, 1991) for singular measures,consists in defining a covering of the support ofthe singular part of the function with our set ofwavelets of different sizes.As emphasized in (Arneodo et al., 1995b; Bacry et al., 1993; Muzy et al., 1991, 1993, 1994),the branching structure of the WT skeletons of fractal functions inthe \((x,a)\) half-plane enlightens the hierarchical organizationof their singularities ( Figure2(e,f)).The WT skeleton can thus be used as a guide to position, at a consideredscale \(a\ ,\) the oscillating boxes in order to obtain a partition of thesingularities of \(f\ .\)The wavelet transform modulus maxima (WTMM) method amounts to computethe following partition function in terms of WTMMcoefficients (Arneodo et al., 1995b; Bacry et al.,1993; Muzy et al., 1991, 1993, 1994):

\[\tag{9}Z(q,a)=\sum_{l\in \mathcal{L}(a)} {\left( \sup_{\stackrel{\scriptstyle{(x,a')\in l}}{\scriptstyle{a'\leq a}}} \vert T_\psi[f](x,a')\vert \right)}^q \; ,\]

where \(q\in\mathbb{R}\) and the \(\sup\) can be regarded as a way to define a scaleadaptive “Hausdorff-like” partition.Now from the deep analogy that links the multifractal formalism tothermodynamics (Arneodo et al., 1995b; Bohr & Tel, 1988), onecan define the exponent\(\tau(q)\) from the power-law behavior of the partition function:

\[\tag{10}Z(q,a) \; \sim \; a^{\tau(q)} \; ,\;\; a \rightarrow 0^+\; ,\]

where \(q\) and \(\tau(q)\) play respectively the role of the inversetemperature and the free energy.The main result of this wavelet-based multifractalformalism is that in place of the energy and theentropy(i.e. the variables conjugated to \(q\) and \(\tau\)),one has \(h\ ,\) the Hölder exponent, and \(D(h)\ ,\) thesingularity spectrum. This means that thesingularity spectrum of \(f\) can be determined from the Legendre transformof the partition function scaling exponent \(\tau(q)\)(Bacry et al., 1993; Jaffard,1997a,b):

\[\tag{11}D(h)=\min_q(qh-\tau(q)) \; .\]

Monofractal versus multifractal functions

From the properties of the Legendre transform, it is easy to see thathomogeneous monofractal functions that involve singularities of unique Hölder exponent\(h=\partial \tau/\partial q\ ,\) are characterized by a \(\tau(q)\) spectrum which is alinearfunction of \(q\) (Figure3(c)).On the contrary, anonlinear \(\tau(q)\) curve is the signature of nonhomogeneous functions thatexhibitmultifractal properties, in the sense that the Hölder exponent \(h(x)\) is afluctuating quantity that depends upon the spatial position \(x\) (Figure3(c)).As illustrated in Figure3(d), the \(D(h)\) singularity spectrum of a multifractalfunction displays a single humped shape that characterizes intermittent fluctuations correspondingto Hölder exponent values spanning a whole interval \([h_{min}, h_{max}]\ ,\) where \(h_{min}\) and\(h_{max}\) are the Hölder exponents of the strongest and weakest singularities respectively.

Applications of the WTMM method

DNA sequences: monofractality of DNA walks

Figure 4: Probability density functions of the wavelet coefficient values of the “G” DNA walk of 2834 human introns (\(L\geq600\)) (red open symbols) and of the turbulent velocity signal (blue filled symbols). (a) and (b) \(\rho_a(T)\) for the sets of scales \(a=\) 4 (\(\color{Red}\vartriangle\)), 16 (\(\color{Red}\square\)), 64 (\(\color{Red}\circ\)) and \(a=\) 27 (\(\color{Blue}\blacktriangle\)), 144 (\(\color{Blue}\blacksquare\)), 760 (\(\color{Blue}\blacktriangledown\)), 3993 (\(\color{Blue}\bullet\)). (c) and (d) \(a^H\rho_a(a^HT)\) for the same sets of scales; \(H=0.60\) in (c) and \(H=1/3\) in (d). The analyzing wavelet is \(g^{(2)}\ .\)

A DNA sequence is a four-letter (A, C, G, T) text where A, C, G and T stand for the basesadenine, cytosine, guanine and thymine respectively.A popular method to graphically portray the genetic information stored in DNA sequences isto used the so-called “DNA walk” representation (Peng et al., 1992).It consists first in converting the DNA text into a binary sequence by coding for examplewith \(\chi(i)=1\) at a given nucleotide positions and \(\chi(i)=-1/3\) at otherpositions (Voss, 1992), and then in defining the graph of the DNA walk bythe cumulative variables \(f(n)=\sum_{i=1}^n\chi(i)\ .\)The DNA walk obtained with the “G” mononucleotide coding for the largest intron of thehuman dystrophin gene is shown in Figure2(a) for illustration.Figure2(c) illustrates the WT when using an analyzing wavelet of sufficientlyhigh order, namely \(g^{(2)}\) (\(n_\psi=2\)), to get rid of the linear trends in the DNA walklandscape inherent to the heterogeneity of composition of genomicsequences (Arneodo et al., 1995a, 1996).Figure3(a) displays some plots of \(\log_2Z(q,a)\) vs \(\log_2(a)\) for differentvalues of \(q\ ,\) where the partition function \(Z(q,a)\) has been computedon the WTMM skeleton (Figure2(e)), according to the definition(Eq. (9)) for a set of 2184 human introns of size \(L\geq800\)bp.Using a linear regression fit, we then obtain the slopes \(\tau(q)\) of these graphs.As shown in Figure3(c), when plotted versus \(q\ ,\) the data for the exponents \(\tau(q)\)consistently fall on a straight line that is remarkably fitted by

\[\tag{12}\tau(q) = qH-1 \; ,\]

with \(H=0.60\pm0.02\ .\) From the Legendre transform of this linear \(\tau(q)\) (Eq. (11)),one gets a \(D(h)\) singularity spectrum that reduces to a single point:

\[\tag{13}\begin{array}{lcl} D(h) = 1 \; & \text{if} & \; h = H \; ,\\ = -\infty \; & \text{if} & \; h \neq H \; ,\end{array}\]

as the signature of a nowhere differentiable homogeneous fractal signalwith a unique Hölder exponent \(h=H=0.60\ .\)Note that similar good estimates are obtained when using analyzing wavelets of differentorder (e.g. \(g^{(3)}\)).

Within the perspective of confirming the monofractality of DNA walks, we have studied theprobability density function (pdf) of wavelet coefficient values \(\rho_a(T_{g^{(2)}}(.,a))\ ,\)as computed at a fixed scale \(a\) in the fractal scalingrange. According to the monofractal scaling properties,one expects these pdfs to satisfy the self-similarityrelationship (Arneodo et al., 1995a, 1996, 2002):

\[\tag{14}a^H \rho_a(a^H T) = \rho(T) \; ,\]

where \(\rho(T)\) is a “universal” pdf (actually the pdf obtained at scale \(a=1\))that does not depend on the scale parameter \(a\ .\)As shown in Figure4(a,c), when plotting\(a^H \rho_a(a^H T)\) vs \(T\ ,\) all the \(\rho_a\) curves corresponding todifferent scales (Figure4(a)) remarkably collapse ona unique curve when using a unique exponent \(H=0.60\) (Figure4(c)).Furthermore the so-obtained universal curve cannot be distinguished froma parabola in semi-log representation as the signature of monofractalGaussian statistics.Therefore, the fluctuations of DNA walks about the composition induced lineartrends cannot be distinguished from persistent fractional Brownian motion(fBm) \(B_{H=0.60}\) that display long-range correlations (LRC)(\(H>0.5\)) (Arneodo et al., 1996, 2002; Muzy et al., 1994).Similar LRC were found in non-coding sequences as well as in coding regions(e.g. coding exons) in eukaryotic genomes (but not for eubacterial sequencesfor which \(H=0.5\)) as the signature of nucleosomal structure, the first step ofcompaction of DNA in eukaryotic nuclei (Audit et al., 2001, 2002).

Fully developed turbulence

It is now well accepted (Frisch, 1995) that in the fully developed regime, a turbulentflow is likely to be in a universal state that can be experimentally characterized bystatistical quantities such as the multifractal spectra \(\tau(q)\) and \(D(h)\ .\)For more than thirty years, one of the main features recognized experimentally is theintermittency of small scales (Frisch, 1995; Meneveau & Sreenivasan, 1991;Monin & Yaglom, 1975) whichmanifests in a significant departureof the experimental velocity data from the monofractal prediction \(\tau(q)=q/3-1\)of Kolmogorov (K41) (Kolmogorov, 1941) based on the homogeneity assumption that, at eachpoint of the fluid, the longitudinal velocity increments have the same scaling behavior\(\delta v_l(x)\sim l^{1/3}\ ,\) which yields the well known \(E(k)\sim k^{-5/3}\) energyspectrum (Frisch, 1995).The pioneering studies (Anselmet et al., 1984; Frisch, 1995) were performed using thestructure functions method which, as discussed in Muzy et al. (1993), intrinsicallyfails to fully characterize the \(D(h)\) singularity spectrum.

In Figure2(b,d,f), Figure3(b,c,d) and Figure4(b,d) arereported the results of a multifractal analysis of single point longitudinal velocity datafrom high Reynolds 3D turbulence using the WTMMmethod (Arneodo et al., 1998a,c, 1999b; Delour et al., 2001).The data were obtained by Gagne and collaborators in the large wind tunnel S1 of ONERAat Modane.The Taylor scale based Reynolds number is \(R_\lambda\simeq 2000\) and the extent of theinertial range following approximately the Kolmogorov \(k^{-5/3}\) law is about four decades(integral scale \(L\simeq7\)m, dissipation scale \(\eta\simeq0.27\)mm).In Figure2(b) is illustrated a sample of the longitudinal velocity signalof length of about two integral scales, when using the Taylor hypothesis (Frisch, 1995).The corresponding WT and WT skeleton as computed with \(g^{(2)}\) are shown inFigure2(d) and Figure2(f) respectively.As shown in Figure3(b), when plotted versus the scale parameter \(a\)in a logarithmic representation, the annealed average over 28000 integral scales ofthe partition functions \(Z(q,a)\) displays a well defined scaling behavior in the inertialrange for a rather wide range of \(q\) values\[-4\leq q\leq 7\ .\]When processing to a linear regression fit of the data, one gets a non-linear \(\tau(q)\)spectrum, the hallmark of multifractal scaling, that is well approximated by the quadraticspectrum of log-normal processes:

\[\tag{15}\tau(q) =c_1q-\frac{c_2}{2}q^2-d\; ,\]

with \(c_1=0.36\pm0.02\ ,\) \(c_2=0.028\pm0.003\) and where \(d=1\) isthe spatial dimension (1D cut of the 3D velocity field).Similar, quantitative agreement is observed for the \(D(h)\) singularity spectrumin Figure3(d) which displays a remarkable parabolic shape:

\[\tag{16}D(h) = d - \frac{(h+c_1)^2}{2c_2}\; ,\]

that characterizes intermittent fluctuations corresponding to Hölder exponentvalues ranging from \(h_{min}=0.12\) to \(h_{max}=0.60\ ,\) the largest dimension\(D(h(q=0))=-\tau(0)=0.999\pm0.001=d\) being attained for \(h=c_1=0.36\pm0.02\ ,\)i.e., a value which is slightly larger than the K41 prediction \(h=1/3\ .\)This multifractal diagnosis is confirmed in Figure4(b) where the pdfof WT coefficients has a shape which evolves across scales from Gaussian at largescales to more intermittent profiles with stretched exponential-like tails at smallerscales.As illustrated in Figure4(d), there is no way to collapse all the WT pdfson a single curve with a unique exponent \(H\) as expected from the self-similarityrelationship (14).Instead, this can be done (Arneodo et al., 1997b, 1998c, 1999b) by usinga Gaussian kernel that strongly supports the log-normal cascadephenomenology (Castainget al., 1990; Delour et al., 2001; Kolmogorov, 1962; Oboukhov, 1962) offully developed turbulence.

Generalizing the WTMM method to d-dimensional image analysis

The generalization of the WTMM method in higher dimension is directlyinspired from Mallatet al. (Mallat & Hwang, 1992; Mallat & Zhong, 1992)reformulation of Canny multiscale edge detector (Canny, 1986) in terms of 2DWT. The general idea is to start smoothing the discrete image data byconvolving it with a filter and then compute the gradient of thesmoothed image. This method has been implemented, tested and appliedto 2D (Arneodo et al., 1999a, 2000, 2003; Decoster et al., 2000; Roux et al., 2000) and3D (Kestener & Arneodo, 2003) scalar field.

The d-dimensional WTMM method

Let us define d analyzing wavelet \(\psi_i(\mathbf{x}=(x_1,x_2,\ldots,x_d))\) thatare respectively, the partial derivatives of a smoothing scalarfunction \(\phi(\mathbf{x})\ :\)\[\tag{17}\psi_i(\mathbf{x}=(x_1,x_2,\ldots,x_d))=\partial \phi(\mathbf{x}=(x_1,x_2,\ldots,x_d))/\partial x_i\; , \;\; i=1,2,\ldots,d. \]

\(\phi(\mathbf{x})\) is supposed to be an isotropic function that depends on\(|\mathbf{x}|\) only and that is well localized around \(|\mathbf{x}|=0\ .\) Commonly usedsmoothing functions are the Gaussian function :

\[\tag{18}\phi(\mathbf{x}=(x_1,x_2,\ldots,x_d))=e^{-|\mathbf{x}|^2/2},\]

and the isotropic Mexican hat :

\[\tag{19}\phi(\mathbf{x}=(x_1,x_2,\ldots,x_d))=(d-x^2) e^{-|\mathbf{x}|^2/2},\]

that correspond to a first-order (\(n_{\boldsymbol{\psi}}=1\)) and a third-order(\(n_{\boldsymbol{\psi}}=3\)) analyzing wavelet respectively.

For any scalar function \(f(x_1,x_2,\ldots,x_d)\in L^2\)(\(\mathbb{R}^d\)), the WTat point \(\mathbf{b}\) and scale \(a\) can be expressed in a vectorialform (Arneodo et al., 2000; Decoster et al., 2000; Kestener & Arneodo,2003):

\[\tag{20}{\mathbf T}_{\boldsymbol{\psi}} [f] ({\mathbf b},a)= \left\{\begin{array}{l} T_{\psi_1}[f] = a^{-d} \int d^d {\mathbf x} \; \psi_1 \bigl( a^{-1} ( {\mathbf x} - {\mathbf b} ) \bigr) f( {\mathbf x} ) \\ T_{\psi_2}[f] = a^{-d} \int d^d {\mathbf x} \; \psi_2 \bigl( a^{-1} ( {\mathbf x} - {\mathbf b} ) \bigr) f( {\mathbf x} )\\ \qquad\quad\vdots \\ T_{\psi_d}[f] = a^{-d} \int d^d {\mathbf x} \; \psi_d \bigl( a^{-1} ( {\mathbf x} - {\mathbf b} ) \bigr) f( {\mathbf x} ) \end{array}\right\}\]

Figure 5: 2D wavelet transform analysis of a Landsat image of marine Sc clouds captured at \(l=30\)m resolution on July 7,1987, off the coast of San Diego (CA)(Arneodo et al., 1999a; Roux et al., 2000). (a) 256 grey-scale coding of a (\(1024\times 1024\)) portion of the original radiance image. In (b) \(a=2^{2.9}\sigma_W\ ,\) (c) \(a=2^{1.9}\sigma_W\) and (d) \(a=2^{3.9}\sigma_W\) (where \(\sigma_W=13\) pixels \(\simeq 390\)m), are shown the maxima chains; the local maxima of \(\mathcal{M}_{\boldsymbol{\psi}}\) along these chains are indicated by (\(\bullet\)) from which originates an arrow whose length is proportional to \(\mathcal{M}_{\boldsymbol{\psi}}\) and its direction (with respect to the \(x\)-axis) is given by \(\mathcal{A}_{\boldsymbol{\psi}}\ ;\) only the central (\(512\times 512\)) part delimited by a dashed square in (a) is taken into account to define the WT skeleton. In (b), the smoothed image \(\boldsymbol{\phi}_{\mathbf{b},a} * I\) is shown as a grey-scale coded background from white (min) to black (max). \(\boldsymbol{\psi}(\mathbf{x})\) is the first-order radially symmetric analyzing wavelet.

Figure 6: The WT skeleton defined by the maxima lines obtained after linking the WTMMM detected at different scales.

Then, after a straightforward integration by parts, \({\mathbfT}_{\boldsymbol{\psi}}\) can be expressed as the gradient vector field of \(f(\mathbf{x})\)smoothed by dilated versions \(\phi(\mathbf{x}/a)\) of the smoothing filter. Ata given scale \(a\) the WTMM are defined by the positions \(\mathbf{b}\) wherethe modulus \(\mathcal{M}_{\boldsymbol{\psi}}[f]({\mathbf b},a)=|{\mathbf T}_{\boldsymbol{\psi}}[f](\mathbf{b},a)|\) is locally maximum along the direction of the WTvector. These WTMM lie on connected (d-1) hypersurfaces called maximahypersurfaces (see Figure5 and Figure8).In theory, at each scale \(a\ ,\) one only needs to record the position ofthe local maxima of \(\mathcal{M}_{\boldsymbol{\psi}}\) along the maximahypersurfaces together with the value of \(\mathcal{M}_{\boldsymbol{\psi}}[f]\) and thedirection \(\mathcal{A}_{\boldsymbol{\psi}}[f]({\mathbf b},a)\) of \({\mathbf T}_{\boldsymbol{\psi}}[f]\ .\)These WTMMM are disposed along connected curves across scalescalled maxima lines living in a (d+1)-space\((x_1,x_2,\ldots,x_d,a)\ .\) The WT skeleton is then defined as the setof maxima lines that converge to the \((x_1,x_2,\ldots,x_d)\) hyperplanein the limit \(a\rightarrow 0^+\) (see Figure6). As originallydemonstrated in Arneodo et al. (1999a),Decoster et al. (2000) and Kestener & Arneodo (2003),the local Hölderregularity of \(f(\mathbf{x})\) can be estimated from the power-law behavior of\(\mathcal{M}_{\boldsymbol{\psi}}[f] \bigl(\mathcal{L}_{{\mathbf x}_0}(a) \bigr)\sim a^{h({\mathbf x}_0)}\)along the maxima line \(\mathcal{L}_{{\mathbf x}_0}(a)\) pointing to thepoint \({\mathbf x}_0\) in the limit \(a\rightarrow 0^+\ ,\) provided\(h({\mathbf x}_0)\) be smaller than the number \(n_{\boldsymbol{\psi}}\)(\(=\min_j n_{\psi_j}\)) of zero moments of the analyzing wavelet \(\boldsymbol{\psi}\ .\) Then,very much like in 1D (Eq. (9)), one can use thescale-partitioning given by the WT skeleton to define the followingpartition functions :

\[\tag{21}{\mathcal Z}(q,a)=\sum_{{\mathcal L}\in {\mathcal L}(a)} \left ( \mathcal{M}_{\boldsymbol{\psi}}[f]({\mathbf x},a)\right)^q \; ,\]

where \(q\in\mathbb{R}\) and \({\mathcal L}(a)\) is the set of maxima lines that existat scale \(a\) in the WT skeleton. As before, the \(\tau(q)\) spectrumwill be extracted from the scaling behavior of \({\mathcal Z}(q,a)\)(Eq. (10)) and in turn the \(D(h)\) singularity spectrum will beobtained from the Legendre transform of \(\tau(q)\)(Eq. (11)) (Arneodo et al., 1999a; Decoster et al., 2000; Kestener & Arneodo, 2003).

Application of the 2D WTMM method to high-resolution satellite images of cloud structure

Figure 7: Determination of the multifractal spectra of radiance Landsat images of marine Sc using the 2D WTMM method with either a first order (\(\color{Magenta}\bullet\)) or a third-order (\(\color{Turquoise}\boldsymbol{\circ}\)) radially symmetric analyzing wavelet (Eqs (17),(18) and (17),(19) respectively). \(\tau(q)\)vs \(q\ ;\) (b) \(D(h)\)vs \(h\ .\) In (a) and (b), the solid lines correspond to the theoretical multifractal spectra for log-normal cascade processes namely, Eqs (15) and (16) with parameter values \(c_1=0.38\) and \(c_2=0.07\) and \(d=2\ .\) The \(D(h)\) singularity spectrum of longitudinal velocity (dotted line) and temperature (dashed line) fluctuations in fully developed turbulence are shown for comparison in (b).

Stratocumulus are one of the most studied clouds types (Davis et al., 1996).Being at once persistent and horizontal extended,marine Sc layers are responsible for a large portion of the earth'sglobal albedo, hence, its overall energy balance.Figure5(a) shows a typical 1024x1024 pixels portion among 14overlapping subscenes of the original Sc Landsat images wherequasi-nadir viewing radiance at satellite level is digitized on aneight-bit grey scale. The different steps of the 2D WTMM methodologyare illustrated in Figure5 (b,c,d) where the WTMM chains and thelocal maxima of \(\mathcal{M}_{\boldsymbol{\psi}}\) along these chains computed with thefirst order (\(n_{\boldsymbol{\psi}}=1\)) analyzing wavelet, are shown at differentscales. In Figure7 are reported the \(\tau(q)\) and \(D(h)\)multifractal spectra obtained from the scaling behavior of \({\mathcalZ}(q,a)\) over the range of scales \(390\) m \(\lesssim a \lesssim 3120\) m(Arneodo et al., 1999a; Roux et al., 2000).Both spectra are clearly non linear and very well fitted by thetheoretical quadratic spectra of log-normal cascade processes(Eqs (15) and (16) with \(d=2\)). However, with thefirst-order analyzing wavelet, the best fit is obtained with theparameter values \(c_1=0.38\pm 0.02\) and \(c_2=0.070\pm 0.005\ ,\) whilefor the third-order wavelet these parameters take slightly differentvalues, namely \(c_1=0.37\pm 0.02\) and \(c_2=0.060\pm 0.005\ .\) Theintermittency coefficient \(c_2\) is therefore somehow reduced whengoing from \(n_{\boldsymbol{\psi}}=1\) to \(n_{\boldsymbol{\psi}}=3\ .\) Actually, it is a lack ofstatistical convergence because of insufficient sampling which is themain reason for this uncertainty in the estimate of \(c_2\ .\)

In Figure7(b) are shown for comparison the \(D(h)\) singularityspectra of turbulent longitudinal velocity data recorded at the Modanewind tunnel (\(R_\lambda\simeq 2000\)) and of temperature fluctuationsrecorded in a \(R_\lambda=400\) turbulent flow (Ruiz-Chavarria et al., 1996).The \(D(h)\)curve for marine Sc clouds is much wider than the velocity \(D(h)\)spectrum (the intermittency coefficient \(c_2\) being almost three timelarger) and it is rather close to the temperature \(D(h)\) spectrum. Ifit is well recognized that liquid water is not really passive, theresults derived with the 2D WTMM method in Figure7 show thatfrom a multifractal point of view, the intermittency captured by theLandsat satellite looks statistically equivalent to the intermittencyof a passive scalar in fully developed 3D turbulence. The fact thatthe internal structure of Sc cloud somehow reflects some statisticalproperties of atmospheric turbulence is not such a surprise in thishighly turbulent environment (Arneodo et al., 1999a; Roux et al., 2000).

Application of the 3D WTMM method to 3D isotropic turbulence simulations

A central quantity in the K41 theory of fully developed turbulence isthe mean energy dissipation \(\epsilon = \frac{\nu}{2}\sum_{i,j}(\partial_j v_i +\partial_i v_j)^2\) which is supposed to beconstant. Indeed, \(\epsilon\) is not spatially homogeneous but undergoeslocal intermittent fluctuations (Frisch, 1995; Meneveau & Sreenivasan, 1991). There havebeen early numerical and experimental attempts to measure themultifractal spectra of \(\epsilon\) or of its 1D surrogateapproximation \(\epsilon'=15\nu (\partial u/\partial x)^2\) (where \(u\) isthe longitudinal velocity) using classical box counting techniques (Meneveau & Sreenivasan, 1991).In Figure8 and Figure9 are reported the results of theapplication of the 3D WTMM method (Kestener & Arneodo, 2003)to isotropic turbulence direct numerical simulations (DNS) dataobtained by Meneguzzi with the same numerical code as previously usedby Vincent & Meneguzzi (1991), but at a\((512)^3\) resolution and a viscosity of \(5.10^{-4}\) corresponding to aTaylor Reynolds number \(R_\lambda = 216\) (one snapshot of thedissipation 3D spatial field). The main steps of the 3D WT computationare illustrated in Figure8. Note that the WTMMM points thatdefine the WT skeleton, now lie on WTMM 2D surfaces at a given scale.The multifractal spectra obtained from this WT skeleton are shown inFigure9. The \(\tau_\epsilon(q)\) spectrum in Figure9(a)significantly deviates from a straight line the hallmark ofmultifractality. But surprisingly, the data obtained from the 3D WTMMmethod \(\tau_\epsilon^{\mathrm{WT}}(q)\) significantly differ from the spectrum\(\tau_\epsilon(q)=\tau_\epsilon^{\mathrm{BC}}(q)-3q\) estimated withbox-counting technique (Kestener & Arneodo, 2003). Actuallythe WT estimate of the cancellation exponent\(\tau_\epsilon^{\mathrm{WT}}(q=1)+3=-0.19\pm 0.03 < 0\ ,\) the signature ofa signed measure. Indeed, as shown in Figure9(a), the\(\tau_\epsilon^{\mathrm{WT}}(q)\) data are rather nicely fitted by thetheoretical spectrum \(\tau_\mu(q)\) of the nonconservative binomial \(p\)model (Mandelbrot, 1974; Meneveau & Sreenivasan, 1991) with weights\(p_1=0.36\) and \(p_2=0.78\) (\(p_1+p_2=1.14>1\)). By construction, the BCalgorithms systematically provide a misleading conservative\(\tau_\epsilon(q)\) spectrum diagnostic with\(p=p_1/(p_1+p_2)=0.32\) and \(1-p=p_2/(p_1+p_2)=0.68\ .\) The differencebetween the two spectra is nothing but a fractional integration ofexponent \(H^*=\log_2(p_1+p_2)\sim 0.19\ .\) This result is confirmed inFigure9(b) where the singularity spectrum\(D_\epsilon^{\mathrm{BC}}(h)\) is misleadingly shifted to the right by\(H^*\) (\(= -\) the cancellation exponent), without any change of shape ascompared to \(D_\epsilon^{\mathrm{WT}}(h)\) (Kestener & Arneodo, 2003).This observation seriouslyquestions the validity of most of the experimental and numerical BCestimates of \(\tau_\epsilon^{\mathrm{BC}}(q)\) and\(f_\epsilon^{\mathrm{BC}}(\alpha)=D_\epsilon^{\mathrm{BC}}(h+3)\) spectrareported so far in the literature. Besides the fact that the\(\tau_\epsilon^{\mathrm{WT}}(q)\) and \(D_\epsilon^{\mathrm{WT}}(h)\) spectraseem to be even better fitted by a parabola, as predicted fornon-conservative log-normal cascade processes, these results raise thefundamental question of the possible asymptotic decrease to zero of thecancellation exponent in the infinite Reynolds number limit.

Figure 8: 3D WT analysis of the DNS dissipation field \(\epsilon\ .\) \(\boldsymbol{\psi}\) is the first-order analyzing wavelet (\(\phi(\mathbf{x})\) is the Gaussian). (a) Isosurface plot of \(\epsilon\) in a \((128)^3\) subcube. (b) \(\log(\epsilon)\) in the \((128)^3\) central part as coded using 64 grey levels. (c) Field lines of \({\mathbf T}_{\boldsymbol{\psi}} [f] ({\mathbf b},a)\) for \(a=\sigma_W=7\) pixels. (d) Same as (c) for \(a=2\sigma_W\ .\) (e) WTMM surfaces at scale \(a=2\sigma_W\ ;\) from the WTMMM along these surfaces originates a black segment whose length is proportional to \(\mathcal{M}_{\boldsymbol{\psi}}\) and direction is along the WT vector; the colors on the WTMM surfaces are proportional to \(\mathcal{M}_{\boldsymbol{\psi}}\ .\) (f) same as (e) for \(a=4\sigma_W\ .\)

Figure 9: Multifractal analysis of Meneguzzi DNS simulation data (\(d=3\)). using the 3D WTMM method (\(\circ\)) and BC techniques (\(\bullet\)). (a) \(\tau_\epsilon(q)= \tau_\epsilon^{\mathrm{WT}}(q)\) or \(\tau_\epsilon^{\mathrm{BC}}(q)-3q\)vs \(q\ ;\) the solid (dashed) lines correspond to the \(p\)-model with weights \(p_1=0.36\ ,\) \(p_2=0.78\) (\(p_1=0.32\ ,\) \(p_2=0.68\)). (b) \(D_\epsilon(h)=D_\epsilon^{\mathrm{WT}}(h)\) or \(D_\epsilon^{\mathrm{BC}}(h)=f_\epsilon^{\mathrm{WT}}(h+3)\ ;\) the solid and dashed lines have the same meaning as in (a); the (\(\triangle\)) correspond to some average \(D_{\epsilon'}^{\mathrm{BC}}(h)\) spectrum of experimental (\(d=1\)) surrogate dissipation data (Meneveau & Sreenivasan, 1991).

Perspectives

For many years, the multifractal description has been mainly devotedto scalar measures and functions. However, in physics as well as inother fundamental and applied sciences, fractals appear not only asdeterministic or random scalar fields but also as vector-valueddeterministic or random fields.Very recently, Kestener & Arneodo (2004, 2007)have combined singular value decomposition techniques andWT analysis to generalize the multifractal formalism to vector-valuedrandom fields. The so-called Tensorial Wavelet Transform ModulusMaxima (TWTMM) method has been applied to turbulent velocity andvorticity fields generated in \((256)^3\) DNS of the incompressibleNavier-Stokes equations. This study reveals the existence of anintimate relationship \(D_V(h+1)=D_\omega(h)\) between the singularityspectra of these two vector fields that are found significantly moreintermittent that previously estimated from longitudinal andtransverse velocity increment statistics. Furthermore, thanks to thesingular value decomposition, the TWTMM method looks very promisingfor future simultaneous multifractal and structural (vorticity sheets,vorticity filaments) analysis of turbulent flows (Kestener & Arneodo,2004, 2007).

References

P. Abry, P. Flandrin, M. S. Taqqu & D. Veitch (2000). Wavelets for the analysis, estimation and synthesis of scaling data. InSelf Similar Network Traffic Analysis and Performance Evaluation (K. Park & W. Willinger, eds.), pp. 39–88. Wiley.

P. Abry, R. Baraniuk, P. Flandrin, R. Riedi & D. Veitch (2002a). Multiscale nature of network traffic.IEEE Signal Process. Mag.19, 28–46.

P. Abry, P. Flandrin, M. S. Taqqu & D. Veitch (2002b). Self-similarity and long-range dependence through the wavelet lens. InTheory and Applications of Long Range Dependence (P. Doukhan, G. Oppenheim & M. S. Taqqu, eds.). Birkhauser, Boston.

F. Anselmet, Y. Gagne, E. J. Hopfinger & R. A. Antonia (1984). High-order velocity structure functions in turbulent shear flows.J. Fluid Mech.140, 63-89.

F. Argoul, A. Arneodo, J. Elezgaray & G. Grasseau (1990). Wavelet analysis of the self-similarity of diffusion-limited aggregates and electrodeposition clusters.Phys. Rev. A41, 5537-5560.

A. Arneodo, G. Grasseau & M. Holschneider (1988). Wavelet transform of multifractals.Phys. Rev. Lett.61, 2281-2284.

A. Arneodo, F. Argoul, J. Elezgaray & G. Grasseau (1989). Wavelet transform analysis of fractals: application to nonequilibrium phase transitions. InNonlinear Dynamics (G. Turchetti, ed.), pp. 130-180. World Scientific, Singapore.

A. Arneodo, F. Argoul, E. Bacry, J. Elezgaray, E. Freysz, G. Grasseau, J.-F. Muzy & B. Pouligny (1992). Wavelet transform of fractals. InWavelets and Applications (Y. Meyer, ed.), pp. 286 - 352. Springer, Berlin.

A. Arneodo, E. Bacry, P. V. Graves & J.-F. Muzy (1995a). Characterizing long-range correlations in DNA sequences from wavelet analysis.Phys. Rev. Lett.74, 3293-3296.

A. Arneodo, E. Bacry & J.-F. Muzy (1995b). The thermodynamics of fractals revisited with wavelets.Physica A213, 232-275.

A. Arneodo, Y. d'Aubenton-Carafa, E. Bacry, P. V. Graves, J.-F. Muzy & C. Thermes (1996). Wavelet based fractal analysis of DNA sequences.Physica D96, 291-320.

A. Arneodo, E. Bacry, S. Jaffard & J.-F. Muzy (1997a). Oscillating singularities on Cantor sets: A grand-canonical multifractal formalism.J. Stat. Phys.87, 179-209.

A. Arneodo, J. F. Muzy & S. G. Roux (1997b). Experimental analysis of self-similarity and random cascade processes: Application to fully developed turbulence data.J. Physique II7, 363-370.

A. Arneodo, B. Audit, E. Bacry, S. Manneville, J.-F. Muzy & S. G. Roux (1998a). Thermodynamics of fractal signals based on wavelet analysis: application to fully developed turbulence data and DNA sequences.Physica A254, 24-45.

A. Arneodo, E. Bacry, S. Jaffard & J.-F. Muzy (1998b). Singularity spectrum of multifractal functions involving oscillating singularities.J. of Fourier Anal. & Appl.4, 159-174.

A. Arneodo, S. Manneville & J.-F. Muzy (1998c). Towards log-normal statistics in high reynolds number turbulence.Eur. Phys. J. B1, 129-140.

A. Arneodo, N. Decoster & S. G. Roux (1999a). Intermittency, log-normal statistics, and multifractal cascade process in high-resolution satellite images of cloud structure.Phys. Rev. Lett.83, 1255-1258.

A. Arneodo, S. Manneville, J.-F. Muzy & S. G. Roux (1999b). Revealing a lognormal cascading process in turbulent velocity statistics with wavelet analysis.Phil. Trans. R. Soc. Lond. A357, 2415-2438.

A. Arneodo, N. Decoster & S. G. Roux (2000). A wavelet-based method for multifractal image analysis. I. Methodology and test applications on isotropic and anisotropic random rough surfaces.Eur. Phys. J. B15, 567-600.

A. Arneodo, B. Audit, N. Decoster, J.-F. Muzy & C. Vaillant (2002). Wavelet based multifractal formalism: Application to DNA sequences, satellite images of the cloud structure and stock market data. InThe Science of Disasters: Climate Disruptions, Heart Attacks, and Market Crashes (A. Bunde, J. Kropp & H. J. Schellnhuber, eds.), pp. 26-102. Springer Verlag, Berlin.

A. Arneodo, N. Decoster, P. Kestener & S. G. Roux (2003). A wavelet-based method for multi- fractal image analysis: From theoretical concepts to experimental applications.Adv. Imaging Electr. Phys.126, 1-92.

B. Audit, C. Thermes, C. Vaillant, Y. d'Aubenton Carafa, J.-F. Muzy & A. Arneodo (2001). Long-range correlations in genomic DNA: a signature of the nucleosomal structure.Phys. Rev. Lett.86, 2471-2474.

B. Audit, C. Vaillant, A. Arneodo, Y. d'Aubenton Carafa & C. Thermes (2002). Long-range correlations between DNA bending sites: relation to the structure and dynamics of nucleosomes.J. Mol. Biol.316, 903-918.

E. Bacry, J.-F. Muzy & A. Arneodo (1993). Singularity spectrum of fractal signals from wavelet analysis: exact results.J. Stat. Phys.70, 635-674.

T. Bohr & T. Tèl (1988). The thermodynamics of fractals. InDirection in Chaos, Vol. 2 (B. L. Hao, ed.), pp. 194-237. World Scientific, Singapore.

J. Canny (1986). A computational approach to edge-detection.IEEE Trans. Patt. Anal. Mach. Intell.8, 679-698.

B. Castaing, Y. Gagne & E. J. Hopfinger (1990). Velocity probability density-functions of high Reynolds-number turbulence.Physica D46, 177-200.

P. Collet, J. Lebowitz & A. Porzio (1987). The dimension spectrum of some dynamical systems.J. Stat. Phys.47, 609-644.

A. Davis, A. Marshak, W. J. Wiscombe & R. F. Cahalan (1996). Multifractal characterizations of intermittency in nonstationary geophysical signals and fields - A model-based perspective on ergodicity issues illustrated with cloud data. InCurrent Topics in Nonstationary Analysis (G. Treviño, ed.), pp. 97-158. World Scientific, Singapore.

N. Decoster, S. G. Roux & A. Arneodo (2000). A wavelet-based method for multifractal image analysis. II. Applications to synthetic multifractal rough surfaces.Eur. Phys. J. B15, 739-764.

J. Delour, J.-F. Muzy & A. Arneodo (2001). Intermittency of 1D velocity spatial profiles in turbulence: a magnitudecumulant analysis.Eur. Phys. J. B23, 243-248.

J. D. Farmer, E. Ott & J. A. Yorke (1983). The dimension of chaoticattractors.Physica D7, 153-180.

U. Frisch (1995).Turbulence. Cambridge University Press, Cambridge.

P. Goupillaud, A. Grossmann & J. Morlet (1984). Cycle-octave and related transforms in seismic signal analysis.Geoexploration23, 85-102.

P. Grassberger & I. Procaccia (1983). Measuring the strangeness of strange attractors.Physica D9, 189-208.

P. Grassberger, R. Badii & A. Politi (1988).Scaling laws for invariant measures on hyperbolic and non hyperbolic attractors.J. Stat. Phys.51, 135 -178.

A. Grossmann & J. Morlet (1984). Decomposition of Hardy functions into square integrable wavelets of constant shape.S.I.A.M. J. of Math. Anal.15, 723-736.

T. C. Halsey, M. H. Jensen, L. P. Kadanoff, I. Procaccia & B. I. Shraiman (1986). Fractal measures and their singularities: The characterization of strange sets.Phys. Rev. A33, 1141- 1151.

H. G. E. Hentschel (1994). Stochastic multifractality and universal scaling distributions.Phys. Rev. E50, 243-261.

M. Holschneider (1988). On the wavelet transform of fractal objects.J. Stat. Phys.50,963-993.

M. Holschneider & P. Tchamitchian (1990). Régularité locale de la fonction non-différentiable de Riemann. InLes Ondelettes en 1989 (P. G. Lemarié, ed.), pp. 102-124. Springer, Berlin.

S. Jaffard (1989). Hölder exponents at given points and wavelet coefficients.C. R. Acad. Sci. Paris Sér. I308, 79-81.

S. Jaffard (1991). Pointwise smoothness, two-microlocalization and wavelet coefficients.Publ. Mat.35, 155-168.

S. Jaffard (1997a). Multifractal formalism for functions part I: results valid for all functions.SIAM J. Math. Anal.28, 944-970.

S. Jaffard (1997b). Multifractal formalism for functions part II: self-similar functions.SIAM J. Math. Anal.28, 971-998.

S. Jaffard, B. Lashermes & P. Abry (2006). Wavelet leaders in multifractal analysis. InWavelet Analysis and Applications (T. Qian, M. I. Vai & Y. Xu, eds.), pp. 219–264. Birkhauser Verlag, Basel, Switzerland.

P. Kestener & A. Arneodo (2003). Three-dimensional wavelet-based multifractal method: The need for revisiting the multifractal description of turbulence dissipation data.Phys. Rev. Lett.91, 194501.

P. Kestener & A. Arneodo (2004). Generalizing the wavelet-based multifractal formalism to random vector fields: Application to three-dimensional turbulence velocity and vorticity data.Phys. Rev. Lett.93, 044501.

P. Kestener & A. Arneodo (2007). A multifractal formalism for vector-valued random fields based on wavelet analysis: application to turbulent velocity and vorticity 3D numerical data.Stoch. Environ. Res. Risk Assess. DOI: 10.1007/s00477-007-0121-6 .

A. N. Kolmogorov (1941). The local structure of turbulence in incompressible viscous fluid for very large Reynolds numbers.C. R. Acad. Sci. URSS30, 301-305.

A. N. Kolmogorov (1962). A refinement of previous hypotheses concerning the local structure of turbulence in a viscous incompressible fluid at high Reynolds number.J. Fluid Mech.13, 82-85.

S. Mallat & W. Hwang (1992). Singularity detection and processing with wavelets.IEEE Trans. Info. Theory38, 617-643.

S. Mallat & S. Zhong (1992). Characterization of signals from multiscale edges.IEEE Trans. Patt. Recog. Mach. Intell.14, 710-732.

B. B. Mandelbrot (1974). Intermitttent turbulence in self-similar cascades: divergence of high moments and dimension of the carrier.J. Fluid Mech.62, 331-358.

C. Meneveau & K. R. Sreenivasan (1991). The multifractal nature of turbulent energy-dissipation.J. Fluid Mech.224, 429-484.

A. S. Monin & A. M. Yaglom (1975).Statistical Fluid Mechanics: Mechanics of turbulence, vol. 2. MIT Press, Cambridge, MA.

J.-F. Muzy, E. Bacry & A. Arneodo (1991). Wavelets and multifractal formalism for singular signals: Application to turbulence data.Phys. Rev. Lett.67, 3515-3518.

J.-F. Muzy, E. Bacry & A. Arneodo (1993). Multifractal formalism for fractal signals: The structure-function approach versus the wavelet-transform modulus-maxima method.Phys. Rev. E47, 875-884.

J.-F. Muzy, E. Bacry & A. Arneodo (1994). The multifractal formalism revisited with wavelets.Int. J. Bifurc. Chaos4, 245-302.

A. M. Oboukhov (1962). Some specific features of atmospheric turbulence.J. Fluid Mech.13, 77-81.

G. Paladin & A. Vulpiani (1987). Anomalous scaling laws in multifractal objects.Phys. Rep.156, 147-225.

G. Parisi & U. Frisch (1985). Fully developed turbulence and intermittency. InTurbulence and Predictability in Geophysical Fluid Dynamics and Climate Dynamics (M. Ghil, R. Benzi & G. Parisi, eds.), Proc. of Int. School, pp. 84-88. North-Holland, Amsterdam.

C.-K. Peng, S. V. Buldyrev, A. L. Goldberger, S. Havlin, F. Sciortino, M. Simons & H. E. Stanley (1992). Long-range correlations in nucleotide sequences.Nature356, 168-170.

D. Rand (1989). The singularity spectrum for hyperbolic Cantor sets and attractors.Ergod. Th. Dyn. Sys.9, 527-541.

S. G. Roux, A. Arneodo & N. Decoster (2000). A wavelet-based method for multifractal image analysis. III. Applications to high-resolution satellite images of cloud structure.Eur. Phys. J. B15, 765-786.

G. Ruiz-Chavarria, C. Baudet & S. Ciliberto (1996). Scaling laws and dissipation scale of a passive scalar in fully developed turbulence.Physica D99, 369-380.