HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Figures Only Full Text (PDF) Free Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Lovejoy, S. Articles by Schertzer, D. Search for Related Content PubMed Articles by Lovejoy, S. Articles by Schertzer, D. GeoRef GeoRef Citation Agricola Articles by Lovejoy, S. Articles by Schertzer, D. Related Collections Fractal Approaches Scaling Remote Sensing

SPECIAL SECTION: MULTISCALE MAPPING

Single- and Multiscale Remote Sensing Techniques, Multifractals, and MODIS-Derived Vegetation and Soil Moisture

^a Dep. of Physics, McGill Univ., 3600 University St., Montreal, QC H3A 2T8 Canada
^b Dep. of Applied Mathematics, E.T.S.I.A. Agricultural Engineering, Univ. Politécnica de Madrid, Ciudad Universitaria s.n. 28040, Madrid, Spain, and CEIGRAM–U.P.M., Madrid, Spain
^c GEOTOP, L'Université du Québec à Montréal, Montreal, Canada
^d Université de Paris Est, ENPC, CEREVE, 6-8, Avenue Blaise Pascal, Cité Descartes, 77455 Marne-la-Vallee Cedex, France
^* Corresponding author (anamaria.tarquis{at}upm.es).

All rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher.

Received 12 November 2007.

	ABSTRACT

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

 
Scaling processes are increasingly understood to be the result of nonlinear dynamic mechanisms repeating scale after scale from large to small scales leading to nonclassical resolution dependencies. This means that the statistical properties systematically vary in strong, power-law ways with the resolution. When present in geophysical and remotely sensed fields, it implies that when classical (single-scale) remote sensing algorithms are used to determine surrogates of various geophysical fields, they can at most be correct at the unique (and subjective) calibration resolution. Scaling analysis and modeling techniques were applied to MODIS TERRA Bands 1 through 7 and to the standard derived vegetation and soil moisture indices in order to quantitatively characterize the wide range of scaling of these fields. The scaling exponents we found are not so large; however, they act across wide scale ranges and imply large effects. For example, for the statistics near the mean, the MODIS (500-m) resolution would be biased by a factor of 1.52 when compared with similar results from an "ideal" sensor at 1-mm resolution. Applying the standard index algorithms on lower and lower resolution satellite data, we obtained indices with significantly different statistical properties than if the same algorithm was used at the finest resolution and then degraded to an intermediate value (a difference of a factor 1.54). This shows that the algorithms can, at best, be accurate at the unique calibration scale and this points to the need to develop resolution-independent algorithms based on the scaling exponents.

	INTRODUCTION

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

 
SOIL MOISTURE is an important hydrologic state variable critical to successful hydroclimatic and environmental predictions. Soil moisture varies in both space and time because of spatiotemporal variations in precipitation, soil properties, topographic features, and vegetation characteristics (Corwin et al., 2006). Monitoring all these processes across multiple scales is a fundamental component of many environmental and natural resource issues and remote sensing has long been identified as a technology capable of supporting such development (Viña et al., 2004; Vereecken et al., 2007).

The geophysical and biological processes that determine vegetation, soil moisture, and other surface characteristics are highly nonlinear; they involve interacting structures from planetary down to millimeter scales. A consequence is that the surface radiance fields are statistically interrelated across comparable ranges. Such statistical interrelations (correlations) are the basis of most remote sensing algorithms. A typical approach identifies wavebands that are particularly sensitive to the surface features of interest—such as the vegetation or soil moisture; then semiempirical algorithms are developed to relate them to the radiances. Typically, the algorithms assume that satellite data are homogenous or roughly homogenous at scales below the resolution.

The problem with this approach is that it ignores systematic scale or resolution dependencies due to strong subpixel variability and singles out a single scale for the calibration. This scale is essentially subjective; it is primarily determined by the available technology. A symptom of this scale problem is that when new, higher resolution satellite data are used in algorithms calibrated at earlier, lower resolutions, the algorithm typically must be recalibrated.

In the last 25 yr, an understanding of these strong resolution dependencies has started to emerge, and systematic techniques are now available for analyzing and modeling such behavior. Motivated by the fractal geometry of sets (Mandelbrot, 1983) and the development of (multifractal) cascade models in turbulence, the origin of this behavior has been traced to nonlinear dynamic mechanisms that repeat scale after scale from large to small scales. Paralleling this revolution in our understanding of the consequences of wide-range scaling is the generalization of the notion of scale itself to encompass systems that are scaling but highly anisotropic (generalized scale invariance, Schertzer and Lovejoy, 1985b; Pecknold et al., 1997; Lilley et al., 2004; Lovejoy and Schertzer, 2007). We illustrate some of these ideas using MODIS data and vegetation and soil moisture surrogates derived from them. Although we briefly review the basic multifractal notions, our aim is to perform systematic scale-by-scale analyses on the MODIS data and surrogates.

	Spectral Analysis

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

 
The data we used are calibrated MODIS TERRA radiances (product level B1), with wavelengths indicated in Table 1 ; some of the images are shown in Fig. 1 . The vegetation and soil surface moisture indices are standard products described in Rouse et al. (1973) and Lampkin and Yool (2004):

[1a]

[1b]

where _VI and _SM are the vegetation and soil surface moisture indices (for a complete list of symbols and definitions, see the appendix) and the Is are the radiances (the band number is indicated in the subscript); Eq. [1] is usually stated without reference to any particular scale. We may immediately note that, although the MODIS TERRA data have a resolution of 500 m (Bands 1 and 2 were degraded from 250 m) and the variability of the radiances and surface features continues to much smaller scales, the surrogates are defined at a single (subjective) resolution equal to that of the sensor. One of the applications of our analyses is to investigate how the relations between the surrogates and the bands used to define them change with scale: we anticipated that since the scaling properties are different, making the surrogates with data at different resolutions would produce fields with different properties.

View this table: [in this window] [in a new window]   TABLE 1. Characteristics of the MODIS data, taken over Guadalajara (central Spain) over 512 by 512 pixels, each 500 by 500 m (250 by 250 km2), on 29 July 2006.

View larger version (119K): [in this window] [in a new window]   FIG. 1. Grey-shade images of Guadalajara (central Spain) corresponding to an area of 250 by 250 km. Left column: (top) vegetation index, (middle) Band 1, (bottom) Band 2; right column: (top) soil moisture index, (middle) Band 6, (bottom) Band 7.

[2a]

[2b]

where r is a position vector and k is a wave vector. Since P is quadratic in I, it is a second-order statistic. In the definition, we have included the theoretically motivated ensemble average (denoted ) although, in fact, below we estimate P from a single realization using a fast Fourier algorithm and a standard Hanning window.

In Fig. 2 , we show P estimated for the vegetation and soil moisture indices. The indices were estimated from Eq. [1] at the highest resolution, i.e., they are single-scale surrogates (we discuss multiple-scale surrogates below). To make the contours clearer, we smoothed log₁₀P with a four-point Gaussian smoother. In the figure, we see that the contours are fairly roundish, indicating that, to some approximation, the second-order statistics are isotropic. Careful scale-by-scale analyses of the anisotropy would be rewarding (see, e.g., Pflug et al., 1993; Lewis et al., 1999; Beaulieu et al., 2007), but is outside our present scope. Analysis of the individual radiance bands revealed a degree of isotropy similar to the indices shown in Fig. 2.

View larger version (80K): [in this window] [in a new window]   FIG. 2. A gray shade rendition of the (log) spectral densities (P) of (a) the vegetation index, and (b) the soil moisture index.

[3]

where k is the modulus of the wave vector (the notation indicates angle integration in Fourier space). In Fig. 3a , we see E on a log–log plot showing that, with the exception of the single lowest wave number k = 1, E is quite accurately of the power-law form:

[4]

where β is the spectral exponent. Note that sometimes angle averaging (rather than integration) is performed so that in two dimensions, the corresponding exponent is β – 1. The advantage of using the present (turbulence-based) definition is that if the process is isotropic, then β is independent of the dimension of space (e.g., one-dimensional sections of two-dimensional data will have the same exponent). We can also see from the figure that Band 5 (with strong artifacts, "banding"), indicated by blue in Fig. 3a nevertheless has good scaling except at harmonics of the basic banding wave number. In Table 2 , we indicate the β values (estimated for k 2), finding for all the bands β 1.17 ± 0.08, where the uncertainty indicates the band-to-band spread in exponents.

View larger version (20K): [in this window] [in a new window]   FIG. 3. (a) Spectra of all five bands [E(k)] as a function of the modulus of the wave vector k (and using a Hanning window). In order from top to bottom at the point log10k = 0.7, the curves are: purple = Band 6, black = Band 1, magenta = Band 7, blue (with spikes) = Band 5, light green = Band 2, cyan = Band 4, dark green = Band 3; (b) log spectra compensated by k1.17 [i.e., log10 k1.17 E(k)] vegetation index (red, middle left) with Band 1 (black, top left) and Band 2 (green, bottom left) so that flat curves indicate k–1.17 spectra (note the greatly expanded scale, the fluctuations are small); and (c) spectrum compensated by k1.17 soil moisture (blue, bottom left) with Band 6 (purple, top left) and Band 7 (magenta, middle left).

View this table: [in this window] [in a new window]   TABLE 2. A comparison of the different exponents estimated in this study: β is the spectral exponent, H the basic fluctuation ("conservation") exponent, C1 the codimension characterizing the pure resolution dependence near the mean, is the Levy index, which determines how quickly the resolution exponents vary with statistical moment (intensity) of the fields. TM indicates that estimates are based on the trace moments, Struc indicates that they are based on the fluctuations via the structure functions, qD is the exponent of the probability distribution characterizing the extremes, VI is vegetation index, SM is soil moisture index, mss indicates statistics for the multiple-scale surrogate definition of index (all the other statistics are for single-scale surrogates [sss]). For all except (which is fit in q space), fits are across the smallest spatial scale to a factor of 4 from the lowest (except TM, which also excludes a small scale factor of 2).

	Multifractal Analysis

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

 
Generalized Structure Functions
As indicated in Eq. [2], the spectrum is a second-order statistic; if we assume statistical horizontal translational invariance (statistical homogeneity; this is the spatial counterpart of statistical stationarity, which refers to statistical translational invariance in time), then the Wiener–Khintchin theorem shows that the spectrum is the Fourier transform of the autocorrelation function. Before the advent of multifractals, it was believed that many turbulent and turbulent-like processes could be well modeled by Gaussian processes in which the statistics are completely determined by P(k) [indeed, isotropic Gaussian processes are completely determined by E(k)]. However, to explain the extreme variability associated with the phenomenon of intermittency, cascade models were developed in which a simple mechanism repeats scale after scale in a multiplicative manner while conserving the mean of the process. The resulting statistical behavior is given by

[5]

where is the conserved flux associated with an observable, such as an intensity, at resolution , which is the scale ratio across which the process has been developed, L is the outer scale, l is the inner scale, and K(q) is the (convex) moment scaling exponent (Kolmogorov, 1962; Mandelbrot, 1974). In remote sensing terms, Eq. [5] quantifies the rate at which spatial averaging (i.e., lower and lower resolution pixels, smaller ) smoothes the various powers of the field .

The symbol is used for the cascade process since the prototypical cascade quantity is the turbulent energy flux, denoted , which is conserved from one scale to another in hydrodynamic turbulence. Note that we are discussing systems far from equilibrium and it is not the energy that is conserved but rather the ensemble mean flux of energy from large scales to smaller scales (in fact, is only a true flux in Fourier space, i.e., a flux from low to high wave numbers). In terms of the moments in Eq. [5], this conservation implies that the (appropriately normalized) process respects = 1 so that K(1) = 0. Actually, although this ensemble "canonical" conservation is the most general one, it is not the only one possible. Indeed, more restrictive "microcanonical" conservation principles are quite popular—in spite of the fact that the variability of the corresponding multifractal processes is much lower and that such processes do not have simplifying universal behaviors (see Schertzer and Lovejoy, 1992, and below).

The typical inner scales of turbulent processes in the atmosphere are the viscous dissipation scales, which are typically on the order of millimeters or less. In remote sensing applications, while the actual inner scale of the radiances may be of the same order, the satellite averages them to the larger scale l (here 500 m) so that the relevant resolution scale in Eq. [5] is the resolution of the images. In this case, the small subsensor scales average out most of the small-scale variability. In general, however, they will not completely smooth it out; this is the origin of the interesting strong variability of the extremes, which goes variously under the names multifractal butterfly effect, nonclassical self-organized criticality, first-order multifractal phase transitions, and divergence of statistical moments, which are discussed below (for a review of this and other multifractal properties, see Schertzer et al., 1997).

The behavior described by Eq. [5] is called multiscaling because each statistical moment is scaling but with a different exponent. Continuing with the turbulence example, the famous Kolmogorov law for isotropic turbulence relates the energy flux to velocity fluctuations (v) across a distance x as follows:

[6]

where H is the scaling exponent of the fluctuations with respect to the flux and a is the exponent of the flux. The usual interpretation of this equation is that the equality is in the sense of scaling laws so that, taking the qth powers of both sides and ensemble averaging, we obtain:

obtain:

[7]

We can see that the fluctuations of typical observables have an extra linear scaling term Hq and we have introduced the (generalized) structure function exponent (q) (the usual structure function is for q = 2, for a single realization—no ensemble averaging—the latter is called a variogram). From Eq. [6], we can see that the statistics of the fluctuations across distances x has two components, the first because it depends of the flux at low resolution = L/x, the second because of the scaling relation between the fluctuations and the flux (x^H); H thus characterizes the distance from the (conserved) pure multiplicative process ; it is the degree of non- (scale-by-scale) conservation of the process.

To check this behavior on the MODIS data and surrogates, we need only estimate the gradients using:

[8]

where we have assumed statistical translational invariance and isotropy. It is worth noting that this definition of the fluctuation I is sometimes called the "poor man's wavelet"; other choices of definition are possible. Wavelets provide a systematic framework for this (see, e.g., Holschneider, 1995); in practice, however, the definition of Eq. [8] is usually adequate, the main restriction being that it is only appropriate when 0 H 1, a condition that is usually (although not always) satisfied in geophysical applications (here we will see that H 0.18). For example, when H > 1, one must measure fluctuations with respect to a local linear trend; this can be done either by fractionally differentiating the process (power-law filtering, Schertzer and Lovejoy, 1987), using appropriate wavelets (Bacry et al., 1989), or using the multifractal detrended fluctuation analysis technique (Kantelhardt et al., 2002). Note that by using the Wiener–Khintchin theorem, we obtain a simple relation between the second-order structure function exponent and the spectral exponent: β = 1 + (2); this is indeed approximately verified (see Table 2).

Using Eq. [8] to estimate the fluctuation I, we obtain the generalized structure functions shown in Fig. 4a for Bands 1, 2, 6, and 7 and in Fig. 4b for the vegetation and soil moisture indices. We can see that, except for the smallest factor of 4 to 8, the (multiple) scaling is excellent. To quantify the differences in the scaling, in Fig. 5a and 5b we show the slopes that are our estimates of (q). We can see that (as expected) (q) is concave downward. It is interesting to note that while the vegetation index has a (q) value intermediate between its defining band (q) values, the soil moisture (q) is somewhat lower, outside the range. It is clear from the figures that a direct quantitative intercomparison of these continuous concave curves (q) is not possible; we first need to reduce the problem to a finite number of manageable parameters. Although we discuss this in much more detail below, we may already make a first quantification using the exponent H. If we assume that, due to the scale-by-scale conservation of the underlying nonlinear cascade, K(1) = 0, we see from Eq. [6] that H = (1) [actually, according to Eq. [7], H = (1) + K(a), and, in general, we don't know a; however, the correction K(a) is generally small and will be ignored below]. In Table 2, we compare these estimates, finding that across all seven bands, H = 0.189 ± 0.034 while H_VI = 0.162 and H_SM = 0.144, so that the surrogates have H values that are a little lower than those of the individual bands. To gauge the significance of this effect across the range of scales in the MODIS images (a factor of 512), we see that in the mean [(q = 1) = H], the largest and smallest fluctuations differ by a factor 512^0.162 = 2.74.

View larger version (32K): [in this window] [in a new window]   FIG. 4. (a) Structure functions for Bands 1, 2, 6, and 7 (M = Iq = |I(x)|q, where x is the number of pixels of the "lag"; see Eq. [8]) with moments of order q = 0.1, 0.3, ..., 1.7, and 1.9 shown bottom to top, along with fits across the scaling range (eight pixels up to 256); (b) left, vegetation index and right, soil moisture index single-scale surrogate (sss) with moments of order q = 0.1, 0.3, ..., 1.7, and 1.9 shown bottom to top; and (c) left, vegetation index and right, soil moisture index scale by scale surrogate structure functions (mss) with moments of order q = 0.1, 0.3, ..., 1.7, and 1.9 shown bottom to top.

View larger version (9K): [in this window] [in a new window]   FIG. 5. Single-scale structure function exponents [(q)] for (a) the vegetation index (blue, middle) with Bands 1 (black, top) and 2 (green, bottom), and (b) the soil moisture index (purple, top) with Bands 6 and 7 (pink, bottom and middle, respectively); and (c) multiscale (q) for the vegetation index (blue, top) and the soil moisture index (red, bottom).

It remains to characterize the nonlinear K(q). Although we do this more fully below, we can nevertheless go a step further by characterizing the slope of the (q) curve near the mean, i.e., '(1) = H – K'(1). Defining C₁ = K'(1) as the codimension of the mean (see below for the justification of this terminology), we obtained the mean for all seven bands: C₁ = 0.032 ± 0.010 while C_1VI = 0.033 and C_1SM = 0.010. Although these values appear small, typical values in turbulence are only a bit larger, e.g., C₁ 0.07 for horizontal wind (Schmitt et al., 1992, 1994) and C₁ 0.04 for passive scalars in the horizontal (Lilley et al., 2004). The corresponding values in the vertical are about 1.8 times larger; this is as predicted by the 23/9D model of scaling stratification (Schertzer and Lovejoy, 1985a). Finally, the topography has C₁ 0.12 (pretty much the same for both continents and oceans). For an early review of these and other results, see Lovejoy and Schertzer (1995).

Since C₁ characterizes the statistics of near its mean, it (partially) quantifies the effect of pure changes in satellite (pixel) resolution. For example, using the mean band value C₁ = 0.032, we see that for MODIS we may expect that the bias due to a 500-m resolution with respect to an "ideal" 1-mm-resolution satellite (i.e., assuming that the smallest scale of variability in the radiance field is 1 mm) is (500/0.001)^0.032 = 1.52. This bias will be important when trying to calibrate satellite-derived quantities from in situ measurements.

Multiscale vs. Single-Scale Surrogates
We have already mentioned above that the vegetation and soil moisture surrogates are defined in an ad hoc way using the (subjective, finest available) resolution. The fact that the surrogates have different scaling [(q)s] implies that if the surrogates are defined at different resolutions, their statistical properties will be different; in other words, the single-scale surrogate can be (at most) correct at a single resolution. To see this more clearly, these single-scale surrogates [sss, at scale , ^(s)] can be contrasted with the corresponding multiscale surrogates [mss, at scale , ^(m)]. Mathematically, the difference can be expressed as

[9a]

[9b]

where the maximum scale ratio (satellite image scale/single pixel scale) L/l = and the notation (I_i,) and [^(s)] denotes averaging from this finest resolution up to the intermediate resolution < . The mss is the surrogate that would be obtained by applying an identical algorithm (Eq. [1]) to satellite data at the lower resolution , whereas the sss is the surrogate at the same resolution but based on the scale calibration. In Fig. 4c, we show that the resulting _mss(q) obtained across the same range of scales as the single-scale surrogate _sss (Fig. 4b) is quite different (it is statistically significantly larger; see below). To partially quantify this difference, we see from Table 2 that H_VImss = 0.231 and H_SMmss = 0.213, which are larger than the sss estimates above: H_VIsss = 0.162 and H_SMsss = 0.144. Similarly, determining the corresponding C₁ values we find the mss values C_1VImss = 0.049 and C_1SMmss = 0.059 (c.f., the sss values C_1VIsss = 0.033 and C_1SMsss = 0.010, which are substantially smaller, especially the soil moisture values). Note that whenever an exponent is quoted without the indication of sss or mss, the sss value is understood; we add it explicitly here only to contrast it with the mss estimates. These statistical differences imply that if the indices are calculated at the MODIS resolution and then degraded to an intermediate resolution , that the results would be quite different than if the same algorithms were used but directly on data at a lower resolution. The problem gets worse as the range of scales increases because the scaling exponents are different. For example, the exponent difference H_SMmss – H_SMsss = 0.213 – 0.144 = 0.069 so that the difference in the mean index fluctuations across the range of scales of the MODIS image is about a factor 500^0.069 1.54. Alternatively, considering the effect of the inadequate MODIS resolution [500 m rather than, say, a scale of 1 mm, where "ideal" indices might be defined, i.e., = (500/0.001)], then the ratio between the bias due to pixel size in the soil moisture statistics near the mean is ^{C_1SMmss–C_1SMsss}=1.90 (with C_1SMmss – C_1SMsss = 0.059 – 0.01 = 0.049). We therefore conclude that due to the scale dependence of the data, the surrogate indices cannot be valid across a very significant range of scales. This scale dependence must at the very least be taken into account during calibration with in situ data, which are typically at much smaller resolutions. This underscores the need to develop scale-invariant algorithms based on the scale-invariant exponents; see Lovejoy et al. (2001) for a discussion.

Trace Moment Analysis
We have seen that the generic statistical properties of nonlinear processes that are repeated scale after scale are characterized by a nonlinear exponent K(q), and that the observables will generally have an extra linear scaling term Hq. Since at least for low q the linear term Hq is often larger than the nonlinear K(q), in analyses it can mask the nonlinear term. It is therefore advantageous to first estimate the conserved flux from the observed v and then estimate K(q) directly. From Eq. [6], we see that, in principle, this can be done by removing the ^–H scaling. To do this, note that if we start with a field ^a and fractionally integrate it by H (a power-law filter k^–H), the resulting field will have the fluctuation statistics indicated by Eq. [6] (see Marsan et al., 1996). This suggests that, to obtain a flux from v (i.e., a conserved field with H = 0), it suffices to invert the power-law filter, i.e., to fractionally differentiate by an order H. It turns out that a finite difference approximation to this is obtained by an integer order differentiation H' > H followed by taking the absolute value of the result (the absolute value is necessary since the multiplicative cascade > 0; see, however, complex and vector cascades [Schertzer and Lovejoy, 1995]). Since, usually, 0 H 1, a first-order finite difference is sufficient. One simply takes the absolute differences at the finest available resolution and degrades them:

[10]

where again the notation () indicates the average of the finest resolution data across the intermediate resolution . The proportionality constant (l^–H) is unimportant (see [Eq. 6] with v replaced by I and with a = 1, x = l ). In one dimension (e.g., series), one can simply use I(x) = I(x + l) – I(x) where the inner scale is l = L/, where L is again the outer scale. In higher dimensions on data on regular grids (as here), if the data are fairly isotropic, then we can simply treat each line separately and then use the one-dimensional method; however, there are more possibilities. For example, we may use finite difference moduli of gradient vectors, or finite difference Laplacians (the latter being an isotropic order 2 derivative). In two dimensions, this yields, respectively,

[11a]

[11b]

(l = L/ is the interpixel distance). These and other variants have been extensively tested on numerical simulations and on data (see, e.g., Tessier et al., 1993; Lavallée et al., 1993); there is generally not much difference between the various choices (or between direct fractional differentiation using Fourier techniques). Due to analogies with quantum statistical mechanics, the methods based on degrading the resolution of the fluxes and averaging across the available samples is called the trace moment technique (Schertzer and Lovejoy, 1987). It is distinguished from similar partition function methods (Halsey et al., 1986) in that, in addition to spatial averaging, it also involves ensemble averaging. Since here we apply the method to individual realizations (images), the two are equivalent.

In Fig. 6a and 6b , we show the results on the sss indices. The plot shows various statistical moments of at various intermediate resolutions obtained by spatially averaging obtained from the modulus of the gradient vector across larger scales (smaller scale ratios) . Since the intermediate scales x are inversely proportional to the scale ratio ( = L/x), the smallest scales are on the right, the largest scales are on the left. As a first remark, we see that the overall range of variability is quite a bit smaller than for the corresponding structure functions; this is because for the moments shown (q < 2), the Hq term is larger than the nonlinear K(q) term. Notice that for q > 1, the effect of spatial averaging is to decrease the values (K > 0), whereas for q < 1, it increases them (K < 0). Next, we can see that, due to the nonlinear operations at the smallest scales (the definition in terms of the various bands, the modulus of the gradient operator), the scaling is not so good at the finest resolutions; the exponents K(q) (the slopes on the log–log figure) were estimated across the intermediate range shown. At the very largest scales, the scaling is also poor partly because of poor statistics (there aren't very many large-scale structures), but also because of the low wave-number atmospheric corrections mentioned above.

View larger version (19K): [in this window] [in a new window]   FIG. 6. Log–log plot of trace moments (q) using 20 values of scale ratio per order of magnitude in scale: (a) vegetation index; (b) soil moisture index. From bottom to top, q = 0.2,0.4, ..., 1.6, and 1.8.

View larger version (11K): [in this window] [in a new window]   FIG. 7. Trace moment determination of the scaling exponent of the qth order moment of the field [K(q)] for: (a) vegetation index (purple, top right) and soil moisture index (green, bottom right); (b) vegetation index (purple, top right) and Band 1 (orange, middle right) and Band 2 (green, bottom right); and (c) soil moisture index (green) and Band 6 (blue, bottom right) and Band 7 (purple, middle right).

[12]

which states that the probability (Pr) that a (resolution ) value exceeds the threshold is a power of the resolution with exponent –c() (Schertzer and Lovejoy, 1987). The approximation sign is to within constant or slowly varying factors (e.g., including subexponential factors such as log). Note that although, in the above, we used the expression probability distribution, it is not strictly appropriate since it represents the integral of the probability density from to rather than the usual – to (the latter would correspond to a true cumulative distribution function). The exponent is called a singularity since it quantifies the rate at which the field values diverge in the small-scale limit (if < 0, it is in fact a regularity). The function c() defined via Eq. [12] is called the (statistical) codimension function (Schertzer and Lovejoy, 1987) or sometimes the Cramer function (Mandelbrot, 1989). If it is low enough [c() < d, where d is the dimension of the embedding space; d = 2 for single images], then d() = d – c() is the (geometrical) fractal dimension of the set of points (at resolution ) with > . While the geometrical d() must satisfy 0 d() d, the only restriction on c() is that it is nonnegative.

The use of c in place of d is necessary in stochastic processes since events with c > d [i.e., those that would have an impossible negative geometrical dimension d() < 0] correspond to events that are almost surely absent in a d-dimensional space but on the contrary are almost surely present in a large enough sample. Indeed, if we have N_s independent realizations of the process, each across a range of scales , then the effective or sampling dimension of the sample is D_s = logN_s/log, and one would expect to find all levels of activity up to extreme events such that c(_s) = d + D_s, where _s is the sampling singularity (Schertzer and Lovejoy, 1989). In particular, there will be a largest moment q_s = c'(_s) that can be estimated from a finite sample of N_s realizations; for q > q_s, the K(q) becomes (spuriously) linear following the tangent K'(q_s).

As mentioned above, knowledge of K(q) is equivalent to knowledge of c() and visa versa; indeed, the two are conveniently related via Legendre transformations (Parisi and Frisch, 1985):

[13]

This establishes a one-to-one correspondence between orders of singularity and statistical moments q: = K'(q) and q = c'(). We can now see the justification for the term C₁ = K'(1) introduced above; in a precise sense, it corresponds to the singularity that gives the dominant contribution to the mean; in addition, since K(1) = 0, we find C₁ = c(C₁) is also the codimension of the = C₁ singularities.

Before using Eq. [12] to directly estimate c(), we recall that the approximation sign in Eq. [12] means that using the approximation c() –logPr/log is not so accurate. It is therefore better to calculate the histograms at a series of resolutions and then, for a fixed = log/log, regress logPr against log; the slope is –c(). This is the probability distribution multiple scaling technique (Lavallée et al., 1991). We show the result in Fig. 8a, 8b, and 8c . As can be seen at the larger (the extremes), the statistics become poor so that even though, in principle, there is information up until c = d = 2, it is very noisy. The extremes can be examined by determining the histograms directly at the highest resolution. This is especially important since, in general, as mentioned above, in general "canonical" cascade processes there is a fundamental distinction that must be made between the "bare" cascade, quantities _,b, i.e., that developed from the outer scale L down to a smaller (but finite) scale l, and the "dressed" process, _,d, obtained by averaging ("dressing") a bare process down to the small-scale limit (this is a good approximation to measured fields, which typically have variability much below that of their averaging resolution). While the former has all its positive moments converge [the corresponding bare K_b(q) is always finite], the dressed process is more variable; it generally has diverging moments for all moments q > q_D where q_D is a critical exponent dependent on the dimension D across which the process is averaged (Mandelbrot, 1974; Schertzer and Lovejoy, 1987) (K_d(q) = K_b(q) for q < q_D). In general,

[14]

In practice, the moments of a finite data set are always finite; nevertheless, the divergence of moments implies that the K(q) estimated for q > q_D will be dominated by a single (largest) singularity, the corresponding scaling will be spurious, and there will be discontinuity in the slope of K(q) at q = q_D. Since there is a mathematical analogy between classical thermodynamics and multifractals, the discontinuity K'(q) at q = q_D is called a first-order multifractal phase transition. Another term for the phenomenon is the multifractal butterfly effect (Lovejoy and Schertzer, 1998), a term that is justified because—just as for the classical deterministic chaos "butterfly effect"—the large-scale moments for q > q_D are in fact determined by the small-scale details. Finally, the divergence of moments occurs in combination with scaling fractal structures (indeed, it is the because of the buildup of such strong variability from large to small scales that the spatial averaging across finite sets cannot sufficiently "calm" the higher order statistical moments). Elsewhere in statistical physics, the combination of fractal structures with algebraic extreme probabilities has been called self-organized criticality (SOC; Bak et al., 1987); hence, we see that multifractals provide a "nonclassical" route to SOC (Schertzer and Lovejoy, 1994). Actually, this nonclassical SOC may be more generally physically relevant since it based on a (quasi-) constant flux, whereas the classical SOC is in fact only strictly valid in the (often unrealistic) zero-flux limit.

View larger version (13K): [in this window] [in a new window]   FIG. 8. The codimension function [c()] of the set of points with singularities for: (a) the seven bands, in order from top to bottom at the point = 0.23 (below the large spike in the blue curve): blue (with spikes) = Band 5, purple = Band 6, dark green = Band 3, cyan = Band 4, orange = Band 1, magenta = Band 7, and light green = Band 2; (b) the vegetation index (red) and Bands 1 and 3 (yellow and green, respectively); and (c) the soil moisture index (orange, bottom line) and Band 6 (purple, top) and Band 7 (magenta, middle). The straight reference lines (black) have slopes of 6.6.

View larger version (19K): [in this window] [in a new window]   FIG. 9. Log–log plot of the probability (Pr) that a value (indicated by the threshold s), with resolution , exceeds the threshold vs. the order of singularities () (see Eq. [12]): (a) vegetation index (blue, top right), Band 1 (orange, middle right) and Band 2 (green, bottom right); and (b) soil moisture index (bottom right), Bands 6 and 7 (middle and top right, respectively).

[15]

The above Legendre transformation pairs (Schertzer and Lovejoy, 1987) are valid for 0 < 1 and 1 < 2 [when = 1, K(q) = C₁qlogq]. Also, 0 C₁ d, since if C₁ > d, the process cannot be normalized on the subspace. The restriction to q 0 is necessary; since generally for < 2, the q < 0 moments diverge. This is why we did not plot them in Fig. 5 and 7 (the exception is the = 2 lognormal multifractal). Once again, the results for the dressed cascades apply, so that these formulae are only valid for q < q_D and < _D = K'(q_D). Because of this divergence, the terms log-Levy and log-Gaussian are misnomers; the probability tails are actually somewhat stronger than those terms would imply. It could be mentioned that in the literature, a weaker form of universality has been proposed leading to log-Poisson multifractals (She and Levesque, 1994). While log-Poisson multifractals share the infinite divisibility properties of the universal (Levy-generator) multifractals (necessary for cascades continuous in scale), the corresponding processes are neither stable nor attractive so that it would be surprising for them to emerge from complex real-world processes (see Schertzer et al., 1995).

The universal form of Eq. [15] shows that only two parameters are needed to specify the conserved flux of universal multifractals (); for the "observables" (v), we need the third parameter (H). We have seen above that two of these parameters (H and C₁) can be estimated from (1) and '(1); the above shows that only a single additional parameter () is needed to characterize the entire process. This "Levy index" can be estimated either from the radius of curvature of (q) [or K(q)] at q = 1 (equivalently, determined from the second derivative), or—in practice more accurately—by considering K(q) or (q) across a wider range of q values and exploiting the dependence on q. We tested one of these methods (Schmitt et al., 1995) on the MODIS data (for another popular method, double trace moments, see Lavallée et al., 1993). As long as '(0) is finite, we may remove the linear part of (q) by considering the power-law "residue" r(q):

[16]

Applied to universal multifractals (as long as > 1), we have

[17]

Hence we can estimate by a linear regression of log₁₀r(q) vs. log₁₀q. We show this in Fig. 10a and 10b using q values from 0.1 to 3 at intervals of 0.1 ['(0) is estimated from a local quadratic fit to (q) with q = 0, 0.1, and 0.2]. As can be seen, the residue is accurately a power law; for the bands, we find (Table 2) = 1.91 ± 0.03, while for the indices, we find _VI = 2.02 and _SM = 2.23, both of which are slightly outside the theoretically allowed range for universal processes (it can thus only be an approximate empirical characterization). This is presumably a further consequence of the poor scaling of the indices (but not the radiances). Finally, we compare these results to those of a nonlinear regression on the trace moment estimate K(q) curve. Due to the presence of the linear term in the universal K(q) (Eq. [15]), the nonlinear regression is not as accurate. We might also note that a potential problem with it is that one must only make a fit for q < min(q_s, q_D) [recall that for q > q_s or q > q_D, K_d(q) becomes linear and is spurious], and a priori, the appropriate values of q_s and q_D are unknown. Here we have seen that the most extreme q that can be used is on the order q 6 or greater and we made our nonlinear regressions using the relatively well estimated moments q 2.

View larger version (13K): [in this window] [in a new window]   FIG. 10. The nonlinear part of the scaling exponent of the qth order moment of the field [r(q) = q'(0) – (q)] for: (a) the vegetation index (blue, bottom left) and Band 1 (orange, top left) and Band 2 (green, middle left); and (b) soil moisture (red, bottom) with Bands 6 and 7 (pink and blue, respectively, top, nearly superposed).

www.physics.mcgill.ca/gang/multifrac/index.htm

View larger version (74K): [in this window] [in a new window]   FIG. 11. Two simulations to illustrate the effect of a small amount of anisotropy with universal multifractal parameters of basic fluctuation exponent (H) = 0.18, codimension characterizing pure resolution dependence near the mean (C1) = 0.05, and Levy index () = 1.9: (a) the same anisotropy at all scales, so that the scale-changing operator is a standard "zoom" but on an anisotropic initial structure, the degree and type of anisotropy is independent of scale; and (b) differential anisotropy with the direction and elongation of structures slowly changing with scale and the scale-changing operator itself anisotropic so that the degree of anisotropy changes with scale.

	Conclusions

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

In this study, we illustrated these ideas using soil moisture and vegetation indices and surrogates obtained from red and far red MODIS satellite imagery. We systematically reviewed the key notions and analysis methods that have emerged after more than 20 yr of study of cascades and other multifractal processes. First, using the classical technique of Fourier analysis, we showed that the basic isotropic scaling was reasonably well obeyed across most of the available range of scales (here, a factor of 512). We then considered the statistics of fluctuations, which can be systematically analyzed as functions of scale and intensity (by varying the order, q) using generalized structure functions (q). This showed that the basic radiance bands had good multiscaling [i.e., scaling with concave (q)], motivating the introduction of two different surrogates: the first, the classical sss, is defined by a nonlinear combination of radiances at the finest resolution; the second, the mss, is based on the same algorithm except that it is defined on successively lower resolution satellite radiances. Although the basic algorithm is identical (the surrogates are equal to the difference in radiances at two bands divided by their sums), their resolution dependencies were found to be different. This is a consequence of the scale dependence of the radiances and the nonlinear operation defining the surrogates. The fact that the scale-by-scale properties of the sss and mss surrogates are different (a fact quantified further with the parameters H and C₁) shows that unless the finest resolution of the satellite data just happens to coincide with the inner scale of the soil moisture and vegetation indices, the surrogates cannot be more than rough approximations to the indices, valid only near the resolutions at which they were calibrated. In particular, across the MODIS range of scales (a factor 512), we showed that the amplitude of the mean statistics fluctuation in mss and sss was 1.54.

We then went on to make a more complete scale-by-scale statistical characterization of the data and surrogates, showing how to obtain and analyze the underlying conserved fluxes using a technique called trace moments, as well as a scale-invariant characterization of the probability distributions, called the probability distribution multiple scaling technique. These analysis techniques were originally developed to quantify multifractal turbulent processes and we reviewed the relevant theory while explaining the various techniques. At this level, the intercomparison of the scaling involves two exponent functions, one for the moments [K(q)], one for the probabilities [c()]. Although these are related by a Legendre transform—so that they contain essentially the same information—they still effectively involve an infinite number of parameters (e.g., a different exponent for each K or c value). To make the problem manageable, we argued that in the real world, multifractal processes are likely to be in the basin of attraction of the three parameter "universal" multifractal processes so that analysis requires the estimation of just three parameters: in addition to H and C₁ mentioned above, the Levy index , which characterizes the degree of multifractality. We then estimated the remaining index. Overall (Table 2), we found the scaling parameters (H, C₁, and ) were not too different for the different bands but were significantly different for the surrogates (which also had scaling that was less accurately followed, a consequence of the nonlinear transformation at the smallest scale). Even though the exponents H and C₁ seem small, because they act across wide ranges of scale they lead to large effects. For example, we pointed out that the difference in the mean intensity fluctuation across the MODIS image scale range (factor of 512) was about 1.47 depending on whether the indices are defined scale by scale or at a fixed scale. In contrast, there is another multifractal approach to remote sensing (Lévy-Vehel and Mignot, 1994; Lévy-Vehel, 1995) that is agnostic about the existence of wide-range scaling; in practice it only considers very small ranges of scale (e.g., factors of 2–4) so that it leads to applications that are not too different from classical remote sensing techniques.

The similitude of the scaling properties of the different radiance bands does not mean that overall they have the same statistics; the analysis techniques applied in this study are essentially isotropic. Even though by spectral analysis we showed that the anisotropy is not too large, the (unanalyzed) deviations from isotropic scaling can nevertheless lead to significant variations in the morphologies of structures. This is perhaps most graphically shown in the simulations (Fig. 11). Systematic analysis of the scaling anisotropy is not easy, however, and is outside our present scope.

Since our conclusions are based on an analysis of a single scene (albeit at seven wavelengths with nearly 2 x 10⁶ pixels in all), we expect to find differences when data from other sites are analyzed. Since the qualitative (multiscaling) behavior that we found is theoretically expected, and the corresponding exponents are very fundamental, we anticipate that differences will be primarily in the (possibly large, but "normal") sample-to-sample variability of multifractal fields coupled with possible differences in the scale-by-scale anisotropy mentioned above.

Although we criticized classical scale-bound remote sensing algorithms, we did not propose a specific alternative. The general problem is to obtain a scale-invariant characterization of the statistical interrelation between the various radiances or fields. The basic approach is to consider a "state vector," with each relevant field represented by a different component. We are thus led to scale-invariant vectors (Lie cascades, Schertzer and Lovejoy, 1995), an approach that has been discussed to some extent in Lovejoy et al. (2001); however, very little is known about vector multifractals or about the appropriate analysis techniques.

	Appendix: Symbol Definitions

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

 

Ii Intensity in band i (Eq. [1]) VI, SM Vegetation index (VI), soil moisture index (SM) (Eq. [2]) P Spectral density (Eq. [2]) E Spectrum (angle integral of P; Eq. [3]) k Wave vector (Eq. [2]) r Position vector (Eq. [2]) β Spectral exponent (Eq. [4]) L Largest scale in the system, here typically the image scale l Smallest scale in the system, here typically the pixel scale Scale ratio: largest/resolution scale Maximum scale ratio: L/l, Conserved flux associated with an observable, such as an intensity, at resolution (Eq. [5]) , b The "bare" , i.e., the result of a cascade over a ratio that stops at scale , d The "dressed" , i.e., the result of a cascade down to infinitely small scales that is then averaged to resolution K(q) Scaling exponent of the qth order moment of field quantifying how the qth order statistics change with resolution (Eq. [5]) q Order of statistical moments (e.g., Eq. [5]) qs Highest order statistic that can be reliably estimated with Ns samples at dimension d. qD Critical dimension for the divergence of statistical moments, equal to the asymptotic (absolute) probability distribution exponent v The velocity radiance across a distance x =L/ (Eq. [6]) I The velocity radiance across a distance x =L/ x Horizontal distance a Exponent of the flux (Eq. [6]) H Scaling exponent of the fluctuations with respect to the flux (Eq. [6]), scale-by-scale conservation exponent (q) Exponent of the qth order structure function (Eq. [7]) c() (Statistical) codimension of the set of points with singularities (Eq. [12]) Order of singularity corresponding to a field level ( = log/log; Eq. [12]) d Dimension of space, here (for images) d = 2 d() Dimension of the set with singularity [ = d – c()] Ds Sampling dimension = logNs/log Ns Number of independent fields of dimension d, resolution , here it = 1 Pr Probability s Dummy (threshold) variable (Eq. [14]) C1 Codimension corresponding to the mean (i.e., q = 1) Levy index quantifying how quickly the higher and lower order statistics depart from the mean (q = 1) behavior ' Auxiliary variable (Eq. [15]) D Critical singularity corresponding to qD s Critical singularity corresponding to qs r(q) Pure nonlinear part of K(q) (Eq. [17])

	ACKNOWLEDGMENTS

 
We are grateful to Carmelo Alonso for providing the images. This work has been partially supported by CICYT under Project MTM2006-15533-C02-02.

	REFERENCES

TOP ABSTRACT INTRODUCTION Spectral Analysis Multifractal Analysis Conclusions Appendix: Symbol Definitions REFERENCES

 

• Bacry, A., A. Arneodo, U. Frisch, Y. Gagne, and E. Hopfinger. 1989. Wavelet analysis of fully developed turbulence data and a measurement of scaling exponents. p. 703–718. In M. Lessieur and O. Metais (ed.) Turbulence and coherent structures. Kluwer Acad. Press, Dordrecht, The Netherlands.
• Bak, P., C. Tang, and K. Wiesenfeld. 1987. Self-organized criticality: An explanation of the 1/f noise. Phys. Rev. Lett. 59:381–384.[CrossRef][Web of Science][Medline]
• Beaulieu, A., H. Gaonac'h, and S. Lovejoy. 2007. Anisotropic scaling of remotely sensed drainage basins: The differential anisotropy scaling method. Nonlin. Processes Geophys. 14:337–350.
• Corwin, C.L., J. Hopmans, and G.H. de Rooij. 2006. From field- to landscape-scale vadose zone processes: Scale issues, modeling, and monitoring. Vadose Zone J. 5:129–139.[Abstract/Free Full Text]
• Gagnon, J., S. Lovejoy, and D. Schertzer. 2006. Multifractal earth topography. Nonlin. Processes Geophys. 13:541–570.
• Gaonac'h, H., S. Lovejoy, and D Schertzer. 2003. Resolution dependence of infrared imagery of active thermal features at Kilauea volcano. Int. J. Remote Sens. 24: 2323–2344.[CrossRef]
• Halsey, T.C., M.H. Jensen, L.P. Kadanoff, I. Procaccia, and B.I. Shraiman. 1986. Fractal measures and their singularities: The characterization of strange sets. Phys. Rev. A: Gen. Phys. 33:1141–1151.[CrossRef]
• Harvey, D. A., H. Gaonac'h, S Lovejoy, and D Schertzer. 2002. Multifractal characterization of remotely sensed volcanic features: A case study from Kiluaea volcano, Hawaii. Fractals 10:265–274.[CrossRef][Web of Science]
• Holschneider, M. 1995. Wavelets: An analysis tool. Oxford Univ. Press, New York.
• Kantelhardt, J.W., S.A. Zschiegner, K. Koscielny-Bunde, S. Havlin, A. Bunde, and H.E. Stanley. 2002. Multifractal detrended fluctuation analysis of nonstationary time series. Physica A 316:87–114.[CrossRef][Web of Science]
• Kolmogorov, A.N. 1962. A refinement of previous hypotheses concerning the local structure of turbulence in viscous incompressible fluid at high Reynolds number. J. Fluid Mech. 83:349–353.[CrossRef]
• Laferrière, A., and H. Gaonac'h. 1999. Multifractal properties of visible reflectance fields from basaltic volcanoes. J. Geophys. Res. 104:5115–5126.[CrossRef]
• Lampkin, D.J., and S.R. Yool. 2004. Monitoring mountain snowpack evolution using near-surface optical and thermal properties. Hydrol. Processes 18:3527–3542.[CrossRef]
• Lavallée, D., S. Lovejoy, and D. Schertzer. 1991. On the determination of the codimension function. p. 99–110. In D. Schertzer and S. Lovejoy (ed.) Nonlinear variability in geophysics: Scaling and fractals. Kluwer Acad. Publ., Dordrecht, The Netherlands.
• Lavallée, D., S. Lovejoy, D. Schertzer, and P. Ladoy. 1993. Nonlinear variability and landscape topography: Analysis and simulation. p. 171–205. In L. De Cola and N. Lam (ed.) Fractals in geography. Prentice Hall, Englewood Cliffs, NJ.
• Levy, P. 1925. Calcul des probabilites. Gautier-Villars, Paris.
• Lévy-Vehel, J. 1995. Multifractal analysis of remotely sensed images. p. 85–101. In G.G. Wilkinson (ed.) Fractals in geoscience and remote sensing. Official Publ. of the Eur. Commun., Luxembourg.
• Lévy-Vehel, J., and P. Mignot. 1994. Multifractal segmentation of images. Fractals 2:371–377.[CrossRef]
• Lewis, G.M., S. Lovejoy, D. Schertzer, and S. Pecknold. 1999. The scale invariant generator technique for parameter estimates in generalized scale invariance. Comput. Geosci. 25:963–978.[CrossRef]
• Lilley, M., S. Lovejoy, K. Strawbridge, and D. Schertzer. 2004. 23/9 Dimensional anisotropic scaling of passive admixtures using lidar aerosol data. Phys. Rev. E 70(3):036307, doi:10.1103/PhysRevE.70.036307.
• Lovejoy, S., and D. Schertzer. 1995. How bright is the coast of Brittany? p. 102–151. In G. Wilkinson (ed.) Fractals in geoscience and remote sensing. Office for Official Publ. of the Eur. Commun., Luxembourg.
• Lovejoy, S., and D. Schertzer. 1998. Stochastic chaos and multifractal geophysics. p. 38–52. In F.M. Guindani and G. Salvadori (ed.) Chaos, fractals, and models '96: Proc. Int. Conf. Univ. of Pavia, Pavia, Italy.
• Lovejoy, S., and D. Schertzer. 2006. Multifractals, cloud radiances, and rain. J. Hydrol. 322:59–88.[CrossRef]
• Lovejoy, S., and D. Schertzer. 2007. Scale, scaling, and multifractals in geosciences: Twenty years on. p. 311–337. In A.A. Tsonis and J.B. Elsner (ed.) Nonlinear dynamics in geosciences. Springer, New York.
• Lovejoy, S., D. Schertzer, and V. Allaire. 2008. The remarkable wide range spatial scaling of TRIMM precipitation. Atmos. Res. doi:10.1016/j.atmosres.2008.02.016.
• Lovejoy, S., D. Schertzer, Y. Tessier, and H. Gaonac'h. 2001. Multifractals and resolution-independent remote sensing algorithms: The example of ocean colour. Int. J. Remote Sens. 22:1191–1234.[CrossRef]
• Mandelbrot, B.B. 1974. Intermittent turbulence in self-similar cascades: Divergence of high moments and dimension of the carrier. J. Fluid Mech. 62:331–350.[CrossRef]
• Mandelbrot, B.B. 1983. The fractal geometry of nature. Freeman, San Francisco, CA.
• Mandelbrot, B.B. 1989. Fractal geometry: What is it and what does it do? Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, Volume 423, Issue 1864p. 3–16. In M. Fleischmann et al. (ed.) Fractals in the natural sciences. Princeton Univ. Press, Princeton, NJ.
• Marsan, D., D. Schertzer, and S. Lovejoy. 1996. Causal space–time multifractal processes: Predictability and forecasting of rain fields J. Geophys. Res. 101(D21):26333–26346.[CrossRef]
• Parisi, G., and U. Frisch. 1985. A multifractal model of intermittency. p. 84–88. In M. Ghil et al. (ed.) Turbulence and predictability in geophysical fluid dynamics and climate dynamics. North Holland, Amsterdam.
• Pecknold, S., S. Lovejoy, D. Schertzer, and C. Hooge. 1997. Multifractals and the resolution dependence of remotely sensed data: Generalized scale invariance and geographical information systems. p. 361–394. In M.G.D. Quattrochi (ed.) Scaling in remote sensing and geographical information systems. Lewis Publ., Boca Raton, FL.
• Pecknold, S., S. Lovejoy, D. Schertzer, C. Hooge, and J.F. Malouin. 1993. The simulation of universal multifractals. p. 228–267. In J.M. Perdang and A. Lejeune (ed.) Cellular automata: Prospects in astronomy and astrophysics. World Scientific, Singapore.
• Pflug, K., S. Lovejoy, and D. Schertzer. 1993. Generalized scale invariance, differential rotation and cloud texture. J. Atmos. Sci. 50:538–553.[CrossRef]
• Rouse, J.W., Jr., R.H. Haas, J.A. Schell, and D.W. Deering. 1973. Monitoring the vernal advancement and retrogradation (green wave effect) of natural vegetation. Prog. Rep. RSC 1978-1. Remote Sensing Ctr., Texas A&M Univ., College Station.
• Schertzer, D., and S. Lovejoy. 1985a. The dimension and intermittency of atmospheric dynamics. p. 7–33. In B. Launder (ed.) Turbulent shear flow 4. Springer-Verlag, New York.
• Schertzer, D., and S. Lovejoy. 1985b. Generalised scale invariance in turbulent phenomena. Physicochem. Hydrodyn. 6:623–635.
• Schertzer, D., and S. Lovejoy. 1987. Physical modeling and analysis of rain and clouds by anisotropic scaling of multiplicative processes. J. Geophys. Res. 92:9693–9714.[CrossRef]
• Schertzer, D., and S. Lovejoy. 1989. Nonlinear variability in geophysics: Multifractal analysis and simulation. p. 49. In L. Pietronero (ed.) Fractals: Physical origin and consequences. Plenum Press, New York.
• Schertzer, D., and S. Lovejoy. 1992. Hard and soft multifractal processes. Physica A 185:187–194.[CrossRef][Web of Science]
• Schertzer, D., and S. Lovejoy. 1994. Multifractal generation of self-organized criticality. p. 325–339. In M.M. Novak (ed.) Fractals in the natural and applied sciences. Elsevier, Amsterdam.
• Schertzer, D., and S. Lovejoy. 1995. From scalar cascades to Lie cascades: Joint multifractal analysis of rain and cloud processes. p. 153–173. In R.A. Feddes (ed.) Space/time variability and interdependance for various hydrological processes. Cambridge Univ. Press, New York.
• Schertzer, D., and S. Lovejoy. 1997. Universal multifractals do exist! J. Appl. Meteor. 36:1296–1303.[CrossRef]
• Schertzer, D., S. Lovejoy, and F. Schmitt. 1995. Structures in turbulence and multifractal universality. p. 137–144. In M. Meneguzzi (ed.) Small-scale structures in three-dimensional hydrodynamic and magnetohydrodynamic turbulence: Proc. Worksh., Nice, France. 10–13 Jan. 1995. Springer-Verlag, New York.
• Schertzer, D., S. Lovejoy, F. Schmitt, Y. Chigirinskaya, and D. Marsan. 1997. Multifractal cascade dynamics and turbulent intermittency. Fractals 5:427–471.[CrossRef][Web of Science]
• Schmitt, F., D. Lavallee, D. Schertzer, and S. Lovejoy. 1992. Empirical determination of universal multifractal exponents in turbulent velocity fields. Phys. Rev. Lett. 68:305–308.[CrossRef][Web of Science][Medline]
• Schmitt, F., S. Lovejoy, and D. Schertzer. 1995. Multifractal analysis of the Greenland Ice-Core Project climate data. Geophys. Res. Lett. 22:1689–1692.[CrossRef][Web of Science]
• Schmitt, F., D. Schertzer, S. Lovejoy, and Y. Brunet. 1994. Empirical study of multifractal phase transitions in atmospheric turbulence. Nonlin. Processes Geophys. 1:95–104.
• She, Z.S., and E. Levesque. 1994. Universal scaling laws in fully developed turbulence. Phys. Rev. Lett. 72:336–339.[CrossRef][Web of Science][Medline]
• Tessier, Y., S. Lovejoy, and D. Schertzer. 1993. Universal multifractals: Theory and observations for rain and clouds. J. Appl. Meteorol. 32:223–250.[CrossRef]
• Vereecken, H., R. Kasteel, J. Vanderborght, and T. Harter. 2007. Upscaling hydraulic properties and soil water flow processes in heterogeneous soils. Vadose Zone J. 6:1–28.[Abstract/Free Full Text]
• Viña, A., A.A. Gitelson, D.C. Rundquist, G. Keydan, B. Leavitt, and J. Schepers. 2004. Monitoring maize (Zea mays L.) phenology with remote sensing. Agron. J. 96:1139–1147.[Abstract/Free Full Text]
• Wilson, J., S. Lovejoy, and D. Schertzer. 1991. Physically based modelling by multiplicative cascade processes. p. 185–208. In D. Schertzer and S. Lovejoy (ed.) Nonlinear variability in geophysics: Scaling and fractals. Kluwer Acad. Press, Dordrecht, The Netherlands.

This article has been cited by other articles:

				S. D. Logsdon, E. Perfect, and A. M. Tarquis Multiscale Soil Investigations: Physical Concepts and Mathematical Techniques Vadose Zone J., May 27, 2008; 7(2): 453 - 455. [Full Text] [PDF]

This Article Abstract Figures Only Full Text (PDF) Free Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Lovejoy, S. Articles by Schertzer, D. Search for Related Content PubMed Articles by Lovejoy, S. Articles by Schertzer, D. GeoRef GeoRef Citation Agricola Articles by Lovejoy, S. Articles by Schertzer, D. Related Collections Fractal Approaches Scaling Remote Sensing