Multiscaling Description of the Space-Time Structure of Rainfall

Abstract

Understanding the space-time structure of rainfall is one of the most complex problems faced in geophysics and atmospheric sciences. Rainfall has been studied in a similar fashion than turbulence phenomena and many of its mathematical tools have been applied to the description of rainfall patterns. Some of the tools used in turbulence phenomena and rainfall were born in the theory of fractal geometry. Such a theory has turned out to be the most suitable scientific toolbox for explaining complexity, nonlinear processes, and multiscaling dynamics. The aim of this chapter is to illustrate some of these ideas related to the statistical characterization of the space-time structure of rainfall, including the proposals for modeling rainfall fields and some achievements in the construction of a holistic multiscaling theory for rainfall field description.

1. Introduction

Rainfall can be described as irregular patterns that look random, and such a behavior has been mostly associated with the one found in turbulence phenomena. The statistical structure of rainfall is enriched by many self-similar fluctuations of different sizes, which possess a particular fractal dimension; besides, its structure is also characterized by a property known as scaling invariance that explains how rainfall patterns can be understood as self-similar objects. After many years of studying rainfall, the multifractal approach has been the most suitable tool for understanding how rainfall behaves in a statistical sense and also in a physical sense. From this approach, multifractal measures are widely used to study highly irregular distributions of physical quantities defined on a geometric support [1, 2].

Originally, the application of the multifractal approach for studying physical systems was born in the study of turbulence phenomena [3, 4, 5, 6]. The techniques and mathematical tools used in such phenomena were subsequently extended to the study of rainfall. On the study of fully developed turbulent flows, the starting point is in the Navier-Stokes equations:

$${\partial }_t\mathbf{u}+\mathbf{u}{\partial }_x\mathbf{u}={\partial }_{x}p+\nu {\partial }_{\mathit{xx}}\mathbf{u};{\partial }_x\mathbf{u}=0;\mathbf{u}\left(0,\mathbf{x}\right)=u_0$$E1

in which a broken symmetry property has been identified on them. One can identify on these equations that their scaling properties are represented by infinitely many scaling groups, particularly when the viscous term vanishes. Moreover, the Navier-Stokes equations are invariant under scaling transformations $g_{\lambda }$, that is:

$$g_{\lambda }:\mathbf{x}\to \lambda \mathbf{x};\mathbf{u}\to {\lambda }^H\mathbf{u};t\to {\lambda }^{1-H}t$$E2

where $H$ represents a self-similarity scaling exponent. As the energy production is given at a unique scale ${\lambda }_0$ (integral scale), the invariance property cannot be sustained for all scales. Kolmogorov’s theory [7] postulates that there is a finite rate of energy dissipation per unit mass in the turbulent flow, and that such a flow is statistically self-similar with a scaling exponent $H=1/3$.

Previous studies about intermittency in fully developed turbulent flows use the structure function to identify the self-similarity scaling exponent. For instance, the (longitudinal) structure function $S_q\left(\mathrm{\ell }\right)$ is defined as the $q$-th order moment of the velocity fluctuation on a distance $\mathrm{\ell }$, this is:

$$S_q\left(\mathrm{\ell }\right)\equiv ⟨{\left|\mathbf{u}\left(\mathbf{x}+\mathrm{\ell }\right)-\mathbf{u}\left(\mathbf{x}\right)\right|}^q⟩\equiv ⟨{\left|\mathrm{\Delta }\mathbf{u}\left(\mathrm{\ell }\right)\right|}^q⟩$$E3

where $q$ is the order of the moment, $\mathbf{u}\left(\mathbf{x}\right)$ is the velocity at location $\mathbf{x}$, and the operator $⟨\cdot ⟩$ represents an average. For homogeneous isotropic turbulence, the structure function follows a power law of the form:

$$S_q\left(\mathrm{\ell }\right)\propto {\mathrm{\ell }}^{\zeta \left(q\right)}$$E4

where $\zeta \left(q\right)$ represents a scaling function depending on the order of the moment $q$. In Kolmogorov’s theory, the energy spectrum of velocity signals for high-Reynolds flows follows a power law $\propto k^{-5/3}$ over a specific range of scales (inertial range) [7, 8]. Such a law can also be related by the structure function when $q=2$ and the value of $\zeta \left(2\right)$ corresponds to $2/3$, yielding a self-similarity scaling exponent equals to $H=\zeta \left(1\right)=1-\zeta \left(2\right)=1/3$.

However, measurements and numerical experiments have determined that there exists a conditional self-similarity property in data since the structure function scaling exponents change for higher order of the moments, which could possibly be associated with the complex nature of intermittency in turbulent flows [3]. While a theory suggests that the energy dissipation in fully developed turbulent flows is concentrated in a fractal set whose dissipative structures have a fractal dimension that can only be determined in the sixth-order structure function and its value is near to the embedding dimension [3], another theory suggested that the explanation of such a behavior in the structure function is given by considering that there is a hierarchy of singularities of rates of dissipation $\epsilon$ of turbulent kinetic energy, which are distributed on interwoven sets of different fractal dimensions [5, 6].

The description of such a distribution can be obtained by thinking that the dissipation energy flux occurs in a domain of size ${\lambda }_0$ (or eddies of size ${\lambda }_0$), which is subdivided into two subdomains of size ${\lambda }_1={\lambda }_2={\lambda }_0/2$, but the dissipation energy flux is distributed in proportions $p$ and $1-p$ (with $p<1$) for both resulting subdomains. The domain subdivision process is continued by applying the aforementioned distribution rule to each subdomain up to subdomains of size close to the Kolmogorov scale $\eta$, and the distribution of dissipation energy fluxes turns out by following the $p$-proportion rule to every sibling subdomains. At the end of the energy distribution process, there will be a structure statistically enriched by rates of dissipation $\epsilon$ of different sizes, which can be described by a hierarchy of singularities $\alpha$ with multiplicity of fractal dimensions $f\left(\alpha \right)$; thus, a multifractal description of the energy dissipation is suitable to explain the energy distribution and its statistical structure.

Turning to the study of rainfall, the use of the structure function exhibits some subtle differences. For rainfall time series, its structure function is defined in similar fashion as follows:

$$S_q\left(\tau \right)\equiv ⟨{\left|\mathbf{R}\left(t+\tau \right)-\mathbf{R}\left(t\right)\right|}^q⟩\equiv ⟨{\left|\mathrm{\Delta }\mathbf{R}\left(\tau \right)\right|}^q⟩$$E5

where $\mathbf{R}\left(t\right)$ represents the rainfall field at time $t$, $\tau =\mathit{nt}/T$ (for $n=1,2,\dots$) is the aggregation timescale, and $\mathrm{\Delta }\mathbf{R}\left(\tau \right)$ is the rainfall fluctuation for the aggregation timescale $\tau$, and $T$ is the observational period of time. Given the intrinsic relation between rainfall and turbulence, the expected rainfall’s structure function follows a power law given by:

$$S_q\left(\tau \right)\propto {\tau }^{\xi \left(q\right)}$$E6

where $\xi \left(q\right)$ represents a scaling function depending on the order of the moment $q$. Figure 1 shows an example of rainfall fluctuations time series for different aggregation timescales, whose data comes from IMERG satellite retrievals. One can easily appreciate in Figure 1 rough changes in magnitudes and intermittency from scale to scale. Figure 2 includes the corresponding structure function for the rainfall time series at aggregation timescales as shown in Figure 1. For $q>0$, the scale invariance property is sustained for a noncompact set of aggregation timescales (i.e., from $\tau =30$ min to $\tau =9$ h); therefore, we could suggest that there is also a conditional self-similarity property in the statistical structure of rainfall due to its intrinsic intermittency. In the technical literature, it is suggested for the first-order structure function $S_1\left(\tau \right)$ that the scaling exponent $\xi \left(1\right)$ specifies how far the rainfall field is from a conserved field, and therefore, it indicates its order of fractional integration [9, 10, 11, 12]. For the example shown in Figure 2, the scaling exponent values for the first four order structure functions are $\xi \left(1\right)\approx 0.8\approx 4/5$, $\xi \left(2\right)\approx 1.2\approx 6/5$, $\xi \left(3\right)\approx 1.4\approx 7/5$, and $\xi \left(4\right)\approx 1.6\approx 8/5$, whose scaling exponent values do not match the ones found in fully developed turbulence flows, that is, $\xi \left(1\right)=1/3$, $\xi \left(2\right)=2/3$, $\xi \left(3\right)=3/3$, and $\xi \left(4\right)=4/3$ [3, 4]. As seen in rainfall fluctuations, there seems to exist a conditional self-similarity property alike to the one identified in turbulence phenomena, which points out that the intermittency nature of rainfall patterns changes the symmetry property of the structure function.

Figure 1.

Rainfall fluctuations time series. Each panel depicts rainfall fluctuations $\mathrm{\Delta }R\left(\tau \right)$ for different timescales for the span period 2000 to 2020. The rainfall time series is derived from a spatial average of IMERG rainfall fields over the metropolitan area of Mexico valley on each time step. IMERG precipitation product has a time resolution of 30 minutes and spatial resolution of 0.1°. IMERG is the acronym of integrated multi-satellite retrievals for the global precipitation measurement mission.

Figure 2.

Rainfall’s structure function $S_q\left(\tau \right)$ estimations for several orders of the moments $q$. As seen, the data used for computing the structure functions corresponds to the rainfall fluctuations time series exhibited in Figure 1. The rainfall’s scale invariance mostly holds between the aggregation scales of $\tau =30$ minutes to $\tau =9$ hours. In this figure, the scaling exponent $\xi \left(q\right)$ of Eq. (6) is here obtained from the relation $S_q\left(\tau \right)\propto {\tau }^{\xi \left(q\right)}$ for the scale invariance range.

Based on turbulence theory and its advances in the identification and characterization of multifractality, one can have some special considerations for rainfall analysis: (i) the physics of rainfall appears to be as complex as turbulence and hence physical models designed to describe rainfall physics and its scaling properties can either be implicitly implied in Navier-Stokes equations as a subprocess of the turbulence physics or be an independent physical process that interacts with the turbulence, and whose interactions have as consequence a tantamount scaling behavior; (ii) the observational evidence for rainfall shows that there exists a local scale-invariant property, which determines the existence of a noncompact set of scaling exponents [9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19], yielding the property of infinitely many scaling groups that characterizes strictly multifractal sets may not necessarily be satisfied in the same way for rainfall, therefore, suggesting the definition of its patterns as multiscaling sets.

2. Rainfall multiscaling properties

In a general form, a rainfall field $\mathbf{R}\left(t,\mathbf{x}\right)$ at timescale ${\lambda }_t$ and spatial scale ${\lambda }_x$ will be formed by measures$\left\{{\mu }_i\right\}$ that possess multiscaling properties, that is, there exists a range of scaling (singularity) exponents $A_{\alpha }=\left({\alpha }_{\mathit{min}},{\alpha }_{\mathit{max}}\right)$ defined on fractal sets with dimensions $f\left(\alpha \right)$. As the space or timescale is diminished ($\lambda \to 0$), the measures $\left\{{\mu }_i\right\}$ defined on $\mathbf{R}\left(t,\mathbf{x}\right)$ should satisfy the following power-law relationships:

$${\mu }_{\alpha }\propto {\lambda }^{\alpha };N\left(\alpha \right)\propto {\lambda }^{-f\left(\alpha \right)}$$E7

where the $\alpha$-exponents are known as Hölder scaling singularity exponents of order $\alpha$, $N\left(\alpha \right)$ is the number of elements with singularity exponents of order $\alpha$, $f\left(\alpha \right)$ is the fractal dimension of the measure support with singularity exponents of order $\alpha$, and $\lambda$ represents the space (or time) scale of the rainfall field. The functional relationship between $\alpha$ and $f\left(\alpha \right)$ is known as the multifractal spectrum [1, 2, 20]. The former definition has a geometrical connotation as it implies the existence of a local singularity exponent $\alpha$ for every point of the measure function, which is not necessarily satisfied in observational data. Based on this restriction, the multiscaling properties of rainfall have also a probabilistic counterpart. If $g_X\left(X=\mu \left(\lambda \right)\right)\equiv g_X\left(\mu \right)$ is the probability density function (pdf) of the (random) measures $\left\{{\mu }_i\left(\lambda \right)\right\}$ defined on $\mathbf{R}\left(t,\mathbf{x}\right)$ at the scale $\lambda$, the statistical moments of measures $\left\{{\mu }_i\left(\lambda \right)\right\}$ are given by:

$$M_{\mu }^q\left(\lambda \right)=\mathbb{E}\left[{\mu }^q\left(\lambda \right)\right]={\int }_{\mathbb{R}}{\mu }^q\left(\lambda \right)g_X\left(\mu \right)\mathit{d\mu }$$E8

If Eq. 7 is satisfied for the rainfall measures set, the probability of having ${\mu }_{\alpha }\propto {\lambda }^{\alpha }$ should be proportional to ${\lambda }^{d-f\left(\alpha \right)}$. This is the probability of a ball at scale $\lambda$ corresponding to a set ${\mathcal{C}}_{\alpha }$ with dimension $f\left(\alpha \right)$ and Euclidean embedding dimension $d$. Therefore, the probability of having ${\mu }_{\alpha }\propto {\lambda }^{\alpha }$ can expressed as follows:

$$\underset{\lambda \to 0}{lim}\frac{log{\mathbb{P}}_r\left\{\mu >{\lambda }^{\alpha }\right\}}{log\lambda }=d-f\left(\alpha \right)$$E9

This means that the integral in Eq. (8) can be redefined as:

$$M_{{\mu }_{\alpha }}^q\left(\lambda \right)={\int }_{A_{\alpha }}{\lambda }^{q\alpha +d-f\left(\alpha \right)}d{\mu }_{\alpha }\equiv {\int }_{A_{\alpha }}{\lambda }^{q\alpha +c\left(\alpha \right)}d{\mu }_{\alpha }$$E10

where $c\left(\alpha \right)=d-f\left(\alpha \right)$ is the codimension function for the singularity exponents of order $\alpha$. As suggested by the steepest-decent method, in the limit $\lambda \to 0$, the smallest exponents of the power law into the integral of Eq. (10) dominate; therefore, there exists a scaling exponent function ${\theta }_q={inf}_{\alpha }\left[q\alpha +d-f\left(\alpha \right)\right]$ that explains how the probability changes, this is:

$$M_{\mu }^q\left(\lambda \right)\propto {\lambda }^{{\theta }_q};{\theta }_q=\underset{\lambda \to 0}{lim}\frac{log{M}_{\mu }^q\left(\lambda \right)}{log\lambda }$$E11

As such, the steepest-decent method points out that the integral of Eq. (10) is dominated by those terms where the integrand attains its maximum value, which also implies that:

$$\frac{d}{\mathit{d\alpha }}\left(f\left(\alpha \right)-\mathit{q\alpha }\right)=\frac{\mathit{d\tau }\left(q\right)}{\mathit{d\alpha }}=0$$E12

where $\tau \left(q\right)=f\left(\alpha \left(q\right)\right)-\mathit{q\alpha }\left(q\right)$ is defined as the mass exponent and represents a functional of three variables $q$, $\alpha \left(q\right)$, and $f\left(\alpha \left(q\right)\right)$. Under this definition of the mass exponent $\tau \left(q\right)$, the singularity exponents of order $\alpha$ and the fractal dimension associated with those singularity exponents $f\left(\alpha \right)$ can be computed as follows:

$$\alpha =-\frac{\mathit{d\tau }\left(q\right)}{\mathit{dq}}$$E13

$$f\left(\alpha \right)=\tau \left(q\right)+\mathit{q\alpha }$$E14

Eqs. (13) and (14) represent what is called a Legendre transformation. In a similar manner, the scaling exponent function ${\theta }_q$ is also associated with a Legendre transformation of the form:

$${\theta }_q=\underset{\alpha }{inf}\left\{q\alpha +d-f\left(\alpha \right)\right\};\alpha =\frac{d{\theta }_q}{\mathit{dq}};f\left(\alpha \right)=\underset{\alpha }{inf}\left\{q\alpha +d-{\theta }_q\right\};q=\frac{\mathit{df}\left(\alpha \right)}{\mathit{d\alpha }}$$E15

The function $f\left(\alpha \right)$ that describes the fractal dimension of singularities of order $\alpha$ turns out to be concave, and therefore, $d^{2}f\left(\alpha \right)/d{\alpha }^2\le 0$. Besides, for two singularity exponents $\left\{{\alpha }_i,{\alpha }_j\right\}$ with ${\alpha }_i<{\alpha }_j$, the concavity property suggests that for any value $\epsilon \in \left(0,1\right)$ the following inequality should be satisfied:

$$f\left(\left(1-\epsilon \right){\alpha }_i+\epsilon {\alpha }_j\right)\ge \left(1-\epsilon \right)f\left({\alpha }_i\right)+\epsilon f\left({\alpha }_j\right)$$E16

To summarize, ${\mathbb{P}}_r\left\{\mu >{\lambda }^{\alpha }\right\}$ as the probability distribution of the rainfall (random) measures $\left\{{\mu }_i\left(\lambda \right)\right\}$ at scale $\lambda$ and also the geometric distribution of measures with singularities of order $\alpha$ and fractal dimension $f\left(\alpha \right)$. The previous definitions suggest that there is a power-law description of the rainfall measures statistics and a multiscaling process that shapes the rainfall geometric patterns. The application of these concepts have been extended to the space-time description of rainfall and the statistical characterization of rainfall events [14, 19, 21, 22]. So far, the technical literature reports a good fit of this multiscaling theory for understanding the statistical structure of rainfall and consequently the formely descriptive models of rainfall have been redefined to include multiscaling properties.

3. Multifractal approaches for modeling rainfall

In the last 30 years, several multifractal techniques and theories have been developed for describing and modeling rainfall [10, 17, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34]. The first approaches for such a task apply the phenomenology of multiplicative random cascades (e.g., [10, 14, 16, 26]). First, the so-called binomial random cascade will be studied next, as an application example and to introduce some properties of multiplicative random cascades.

Initially, one can suppose that there is a dyadic interval $I^k=I_{{\beta }_1{\beta }_2\dots {\beta }_k}$ that is selected randomly and that such an interval represents the size of the simulated field and which helps to identify the location of simulated measures. Then a random sequence of digits ${\beta }_1{\beta }_2\dots {\beta }_k$ taking values of 0 or 1 with a uniform probability distribution given by:

$$Pr\left\{{\beta }_i\right\}=\left\{\begin{array}{cc}1/2,& \text{if}{\beta }_i=0\\ 1/2,& \text{if}{\beta }_i=1\end{array}\right.\text{for}i=\left\{\mathrm{0,1,2},\dots ,k-1\right\}\text{with}k=2^n$$E17

where $n=\mathrm{0,1,2},\dots$ may be used to determine the stages that define the random cascade with size determined by the dyadic interval $I^k$. The measures in such a dyadic interval are estimated following a multiplicative process, such that:

$$\mu \left(I^k\right)=\prod _{i=0}^{k-1}{W}_{{\beta }_i};\text{with}{W}_{{\beta }_i}=\left(\begin{array}{cc}{W}_0,& \text{if}{\beta }_i=0\\ W_1,& \text{if}{\beta }_i=1\end{array}\right.$$E18

where $W_{{\beta }_i}$ is a random variable whose possible values are $W_0$ and $W_1$ with probability 1/2, and the measure $\mu \left(I^k\right)$ is computed as products of $k$ independent and identically distributed random variables $W_{{\beta }_i}$. The singularity exponents of order $\alpha$ associated with the measures and resulting from the multiplicative process are estimated as follows:

$${\alpha }_k\left(I^k\right)=\frac{log{\prod }_{i=0}^{k-1}{W}_{{\beta }_i}}{log{2}^{-k}}=-\frac{1}{k}\sum _{i=0}^{k-1}{log}_2{W}_{{\beta }_i}=\frac{1}{k}\sum _{j=1}^k{V}_j$$E19

where $V_j$ is a random variable with probability distribution, which is given by:

$$Pr\left\{V_j\right\}=\left\{\begin{array}{cc}1/2,& \text{for}{V}_0=-{log}_2{W}_0\\ 1/2,& \text{for}{V}_1=-{log}_2{W}_1\end{array}\right.\text{for}j=\left\{1,2,\dots ,k\right\}$$E20

Hence, the singularity exponent ${\alpha }_k$ is also a random variable in which its expected value ${\alpha }_0$ can be computed as:

$$\mathbb{E}\left[V_j\right]=\frac{1}{2}{V}_0+\frac{1}{2}{V}_1=-\frac{1}{2}\left[{log}_2{W}_0+{log}_2{W}_1\right]=-\frac{1}{2}{log}_2{W}_0{W}_1$$E21

where $\mathbb{E}\left[V_j\right]\equiv {\alpha }_0$ represents the average singularity exponent to be found in a set of simulated measures. If one sets a restriction on the product of $W_0$ and $W_1$, such that its product is bounded in the range $0\le W_0{W}_1\le 1$, the exponent ${\alpha }_0$ will be singular (i.e., $0\le {\alpha }_0\le 1$) when the product is bounded in the range $1/4\le W_0{W}_1<1$ and nonsingular (i.e., ${\alpha }_0<1$) when the product is bounded in the range $W_0{W}_1>1/4$.

The singularity of a measures can also be understood through the cumulative measure function. At this stage, let us define the function $M\left(x\right)={\sum }_{i=0}^x\mu \left(i\right)$ as a descriptor of the local variability of the field measures $\mu \left(i\right)=M\left(x\left(i\right)+\lambda \right)-M\left(x\left(i\right)\right)$. With the help of Eq. (7), one can observe in the case when $\lambda =2^{-k}$ and $k\to \infty$, if $\alpha >1$ the measures $\mu \to 0$ and the function $M\left(x\right)$ becomes constant for all values of $x$ (nonsingular measures), but if $0\le \alpha \le 1$ the measures ${\mu }_i\to 1$ and the function $M\left(x\right)$ characterizes the measures as singular for all values of $x$. Therefore, the values to be selected for $W_0$ and $W_1$ will determine the singularity of the binomial random cascade measures.

Now, the variability of binomial random cascade measures can be approximated by the variance of ${\alpha }_k$, thus:

$$\mathrm{Var}\left[V_j\right]=\mathbb{E}\left[V_j^2\right]-{\left(\mathbb{E}\left[V_j\right]\right)}^2$$E22

$$=\left(\frac{1}{2}{V}_0^2+\frac{1}{2}{V}_1^2\right)-{\left(\frac{1}{2}{V}_0+\frac{1}{2}{V}_1\right)}^2=\frac{1}{4}{\left(V_0-V_1\right)}^2$$E23

$$=\frac{1}{4}{\left({log}_2\left(\frac{W_0}{W_1}\right)\right)}^2$$E24

and the expected limits of the singularity exponent $\alpha$ are defined by:

$${\alpha }_1=-\frac{1}{2}\left({log}_2{W}_0{W}_1+{log}_2\left(\frac{W_0}{W_1}\right)\right)=-\frac{1}{2}{log}_2{W}_0$$E25

$${\alpha }_2=-\frac{1}{2}\left({log}_2{W}_0{W}_1-{log}_2\left(\frac{W_0}{W_1}\right)\right)=-\frac{1}{2}{log}_2{W}_1$$E26

where ${\alpha }_1={\alpha }_{\mathit{min}}$ and ${\alpha }_2={\alpha }_{\mathit{max}}$ if $W_0>W_1$. On the other hand, if the strong law of large numbers (SLLN) is applied to Eq. (19), the expected value of the singularity exponent ${\alpha }_k$ will converge, almost surely, to the expected value ${\alpha }_0$ for $k\to \infty$, that is,

$${\mathbb{P}}_r\left\{\underset{k\to \infty }{lim}\frac{1}{k}\sum _{j=1}^k{V}_j={\alpha }_0\right\}=1$$E27

The previous result means that the binary expansion ${\beta }_1{\beta }_2\dots {\beta }_k$ will have the same frequency of zeros and ones with probability $1$. Moreover, for a large number of finite size $2^k$, there will be a range of singularity exponents spanning from ${\alpha }_{\mathit{min}}$ to ${\alpha }_{\mathit{max}}$; thus, deviations from the expected value ${\alpha }_0$ become significant depending on the value of $k$. As one can observe, the binomial random cascade is a model for describing multiscaling processes where its patterns are enriched with singularity exponents of order $\alpha$ but defined as sums of independent and identically distributed random variables $W_{{\beta }_i}$.

As mentioned above, the binomial random cascade approach has been used for describing the distribution of the energy dissipation fluxes in turbulent flows [5, 6]; however, its application for the description of rainfall field measures has not been completely successful. In observational rainfall data, there is not a well-defined compact and dense set of singularity exponents which possibly is due to measurements problems, insufficient length of rainfall records, or the intrinsic aperiodic structure of their patterns.

The next models to be presented exhibit considerations for modeling rainfall under a multiscaling approach. There are always limits for a complete description, and caution must be taken at the moment of selecting the model for particular applications. Let us begin with one pioneering multiscaling (or multifractal) rainfall model, which was introduced by Over and Gupta [16], and is set up on the theory of random cascades. In this model, rainfall fields are based on the subdivisions of a $d$-dimensional domain ${\left[0,L_0\right]}^d$ into $b^n$ subcubes. For instance, if the embedding dimension equals $d=2$ and the branching number equals $b=4$, there are $4^n$ subdivisions of the two-dimensional rainfall field at the $n$-th level of the cascade. As such, the $i$-th subdomain after $n$ partitions is denoted by ${\mathrm{\Delta }}_n^i$ and the length $L_n=L_0{b}^{-n/2}$ of each side of the subdomain ${\mathrm{\Delta }}_n^i$ at the $n$-th level represents the spatial scale.

The construction of the rainfall field begins when a domain ${\left[0,L_0\right]}^d$ adopts an initial rainfall mass $R_0{L}_0^d$, and the $i$-th subdomain of the first subdivision ${\mathrm{\Delta }}_1^i$ consequently takes a rainfall mass equals to $R_0{L}_0^d{W}_1^i{b}^{-1}$ for $i=\left\{1,\dots ,b\right\}$. Here, $W_1^i$ is a mutually independent random variable and the distribution of $W$ is called the cascade generator. Thus, for every subdivision ${\mathrm{\Delta }}_n^i$ of the cascade, the mass is estimated as follows:

$$\mu \left({\mathrm{\Delta }}_n^i\right)=\frac{R_0{L}_0^d}{b^n}\prod _{j=1}^n{W}_j^i$$E28

Following Eq. (28), the cascade limit mass ${\mu }_{\infty }$ is obtained by letting $n\to \infty$ (small-scale limit), so that the rainfall field is computed as:

$${\mu }_{\infty }=\mu \left({\mathrm{\Delta }}_n\right)Z_{\infty }$$E29

where,

$$Z_{\infty }=\frac{{\mu }_{\infty }\left({\left[0,L_0\right]}^d\right)}{R_0{L}_0^d}$$E30

is an independent and identically distributed set of random variables which is employed to designate the high-frequency (or small-scale) component of the rainfall field. $Z_{\infty }$ can also be understood as a normalized mass of the rainfall field when the multiplicative process is developed for the subdivision $n\to \infty$. The low-frequency (or large scale) component of the rainfall field is designated by $\mu \left({\mathrm{\Delta }}_n\right)$, and thus, the final construction of the random cascade is defined by the product of both high- and low-frequency components. Another important aspect in the construction of random cascade fields is their ability to maintain the ensemble mean of the random variable $W$ which is defined as $\mathbb{E}\left[W\right]=1$, thus preserving mass on the average at all scales for each partition of the random cascade. It ought be highlighted that this approach may yield degenerate cascades, that is, those with total mass ${\mu }_{\infty }$ of zero with probability of 1. For a nondegenerate cascade, a suitable cascade generator should be carefully selected and for all cases, and the canonical restriction $\mathbb{E}\left[W\right]=1$ and $\mathbb{E}\left[Z_{\infty }\right]=1$ should be satisfied.

An example of a possible outcome from this multifractal model is illustrated in Figure 3. An advantage of this random cascade model is that it is parsimonious, which means that it requires few parameters to describe the spatial pattern. However, the model foundations are based on the theory of random functions and, for this reason, there are some difficulties in validating the choice of a universal class of random generators for a particular case of study. On the other hand, it ought to be highlighted that a common drawback of random cascade models is that their simulations often exhibit what is called blockiness, that is, a concentration of the mass in blocks into the spatial patterns [26, 35]. Although there exist several procedures for vanishing the anomalous effect of blockiness in random cascade outputs, further research is required for improving the reduction of blockiness from the model itself.

Figure 3.

Random cascade field simulated with the algorithm suggested by reference [16]. The field was simulated for an unitary precipitation rate and a spatial domain of $128$ × $128$ spatial units. The model used a branching number of $4$, the subdivisions for the high-frequency component ($z_{\infty }$) were seven stages and the subdivisions for the low-frequency component ($\mu \left({\mathrm{\Delta }}_n\right)$) were three stages.

Having described a probabilistic approach to represent rainfall fields, here we devote our attention to another model, which defines multifractal measures but via a deterministic approach. This is the so-called Fractal-Multifractal (FM) approach that was introduced by Puente [32] and has been applied to rainfall modeling by several researchers in the last 30 years, for example, [23, 25, 29, 31, 33, 34, 36]. In this model, the observed patterns are built as illuminated canonical deterministic multifractal measures that are supported over the graphs of fractal interpolation function [37]. In the original FM approach, a combination of linear functions of the form $W_n\left(\mathbf{x}\right)=A_n\mathbf{x}+B_n$ (affine transformations) are employed to interpolate (in a fractal way) a set of points in ${\mathbb{R}}^N$ and then, by using projections (derived distributions) of canonical multifractal measures over the graphs of such functions, other, noncanonical, multifractal measures are constructed. For instance, at ${\mathbb{R}}^2$, the FM approach requires a set of at least two contractile affine maps for modeling patterns along a line, that is,

$$W_n\left(\begin{array}{c}x\\ y\end{array}\right)=\left(\begin{array}{cc}{a}_n& 0\\ c_n& d_n\end{array}\right)\left(\begin{array}{c}x\\ y\end{array}\right)+\left(\begin{array}{c}{e}_n\\ f_n\end{array}\right)=\mathit{Ax}+t,n=1,\dots ,N$$E31

subject to these restrains:

$$W_n\left(\begin{array}{c}{x}_0\\ y_0\end{array}\right)=\left(\begin{array}{c}{x}_{n-1}\\ y_{n-1}\end{array}\right),W_n\left(\begin{array}{c}{x}_N\\ y_N\end{array}\right)=\left(\begin{array}{c}{x}_n\\ y_n\end{array}\right),0\le d_n<1,n=1,\dots ,N$$E32

where $a_n$, $c_n$, $e_n$, and $f_n$ are parameters of the FM approach in ${\mathbb{R}}^2$, which are estimated directly from the restrains indicated above. $d_n$ is also a parameter of the FM approach, but this one defines the vertical scaling factor in the transformation $W_n$. Consequently, all affine transformations $W_n\left(\cdot \right)$ become contractile if for all $n$ is satisfied that $0\le d_n<1$. When all these conditions are satisfied, a unique fixed point exists, that is, a fractal interpolation function $f:x\to y$, such that $G=\left\{(x,f\left(x\right))|x\in [0,1]\right\}$, that satisfies $G={\cup }_{n=1}^N{W}_n\left(G\right)$ [37]. During the construction of the fractal interpolation function in ${\mathbb{R}}^2$, two probability measures $\mathit{dy}$ and $\mathit{dx}$ are generated by counting the frequency of points within coordinates $x$ and $y$. Under a suitable parameter configuration, both generated measures are multifractal objects: the one over $x$ is a canonical multinomial multifractal and the one over $y$ may resemble geometrically geophysical patterns [33]. In addition, if for all points at coordinates $x$ and $y$ a rotation $\theta$ is applied, the original FM approach increases its complexity, that is, enhances its space of possibilities for generating yet more interesting fractal patterns. The concept of projections of multifractal measures that are supported over a fractal interpolation function are indeed transformations $W_n\left(x\right):{\mathbb{R}}^2\to {\mathbb{R}}^2$, where the matrix $A_n$ in $W_n\left(x\right)=A_{n}x+t_n$ is rewritten as follows:

$$A_n=\left(\begin{array}{cc}{a}_n& 0\\ c_n& d_n\end{array}\right)\to A_n=\left(\begin{array}{cc}\cos \theta & -sin\theta \\ sin\theta & cos\theta \end{array}\right){\left(\begin{array}{cc}{a}_n& 0\\ c_n& d_n\end{array}\right)}^T$$E33

where the matrix $A_n$ transforms the relative space with rotations or contractions, and the vector $t_n$ specify a linear translational displacement of the map. The new measures under that comes out from the rotational transformation are now denotated as ${\left.\mathit{dy}\right|}_{\theta }$ and ${\left.\mathit{dx}\right|}_{\theta }$. Figure 4 illustrates some simulated rainfall fields via FM approach in ${\mathbb{R}}^2$. While the upper left graph shows a fractal interpolation function and over it the corresponding measures obtained from a binomial cascade, the one on the upper right shows the projection (or illumination) of those multifractal measures for a specific rotation angle $\theta$. The lower left graph of Figure 4 shows the comparison between a simulated rainfall field and the high-resolution rainfall pattern registered at Boston in Octuber 25, 1980 (the right plot) which was previously studied by Puente and Obregón [33]. Some applications of the FM approach for modeling spatial patterns can be found in Refs. [23, 33].

Figure 4.

Example of two simulated precipitation fields via the FM approach. In each subfigure, the left panel shows the construction of a fractal interpolation function in the plane $x$–$y$ and over such a function, multifractal measures $\mathit{dx}$ are exhibited. The projection concept of multifractal measures are represented by the model as a new derived measure ${\left.\mathit{dy}\right|}_{\theta }$ (second panel), which resemble complex geophysical patterns. In the lower subfigure is depicted an output of the FM approach that was originally reported in ref. [33], but adapted in this manuscript by the authors who used an own MATLAB algorithm. The simulated rainfall resembles graphically and statistically the rainfall pattern observed in Boston on October 25, 1980 (lower left panel).

4. Space: time multifractal modeling

There are many applications on meteorology and hydrology in which a good knowledge about the dynamics of the space-time precipitation field would be very useful to know in advance (e.g., early warning systems for flood control). Although there exists many climatic and weather models, trying to obtain prediction to be used as inputs of rainfall-runoff models, there are difficulties to be compatible with the scales among these models. If the output scale from precipitation forecasting models is higher than the necessary scale for running a hydrological model, then it is required to employ a downscaling procedure for parity among data scales, but such a procedure entails to information losses and a reduction of the forecasting likelihood. Clearly, the relationship between scales that are used for running weather models and those ones used for running hydrological models is still a problem to be solved.

In the multiscaling theory of the rainfall field, a self-similar process is assumed to explain the space-time dynamics of rainfall, that is, there are sequences of precipitation data that converge in law after a change of location and timescale. Formally, a process $\mathbf{R}\left(t,\mathbf{x}\right)$ on the interval ${\mathbb{R}}^{+}=\left[0,\infty \right)$ is self-similar with exponent $H$ if for any $c>0$ and $d>0$, the following property is satisfied:

$$\left\{\mathbf{R}\left(\mathit{ct},d\mathbf{x}\right)\right\}\stackrel{\mathrm{d}}{=}\left\{{\left(\mathit{cd}\right)}^H\mathbf{R}\left(t,\mathbf{x}\right)\right\}$$E34

In addition, the self-similarity property represents also the validity of the Taylor’s hypothesis of the “frozen turbulence” which says that the autocovariance at a time lag $\tau =t_2-t_1$ (with $t_2>t_1$), and at a fixed but arbitrary location $x_0$, is the same to the spatial covariance between two spatial points separated by the translation distance $\mathit{u\tau }$ for at a fixed but arbitrary time $t$ and a storm moving at a constant velocity $u$, this is:

$$\mathrm{Cov}\left[\mathbf{R}\left(t_1,x_0\right),\mathbf{R}\left(t_2,x_0\right)\right]=\mathrm{Cov}\left[\mathbf{R}\left(t,x_0\right),\mathbf{R}\left(t,\mathit{u\tau }\right)\right]$$E35

In general form, Taylor’s hypothesis reinterprets time variation at a fixed point as a spatial variation at some fixed time; besides, Taylor’s hypothesis represents a space-time process which is forced by a constant large-scale advection velocity $u=\lambda /\tau$ at the spatial scale $\lambda$. Such advection velocity $u$ satisfies the scaling transformation property $u\to {\lambda }^{H}u$, which was indicated previously in Eq. (2) and $H$ represents a self-similarity scaling exponent.

It is highlighted that there are some physical considerations in the application of Taylor’s hypothesis. Some empirical studies have reported that the large-scale advection velocity induces a space-time anisotropy (e.g., [29, 38]). As the advection velocity changes with time, the anisotropy also changes; therefore, there are limit values for the time lag $\tau$ (measured from the origin of the storm) for the Taylor’s hypothesis to be satisfied. Furthermore, there could be similar near-zero velocities in the rain cells or in the storm region that also break the hypothesis.

Several studies have proved that time and space scaling properties can be represented by power laws in their structure functions such that:

$$S_q\left(\tau \right)\propto {\tau }^{H_{\lambda }\left(q\right)}\left(\text{for}\mathrm{a}\text{local\hspace{0.17em}area\hspace{0.17em}of\hspace{0.17em}size}\lambda \right)$$E36

$$S_q\left(\lambda \right)\propto {\lambda }^{H_{\tau }\left(q\right)}\left(\text{for}\mathrm{a}\text{duration}\tau \right)$$E37

In practice, both space and time structure functions are restricted to some temporal and spatial scales; therefore, the scaling invariance property is bounded (i.e., there is a range in which the multifractal scaling properties are broken or change to represent a monofractal behavior) and such a consideration sets limits in the application of rainfall models. Another problem is the identification of scaling relationships between the exponents $H_{\lambda }\left(q\right)$ and $H_{\tau }\left(q\right)$ or to identify a new one that does not depend directly on either spatial or timescales. Deidda [38] suggests that rainfall can be modeled as a self-similar multifractal process (for a bounded spatial range ${\lambda }_0\le \lambda \le L$) with a constant advection velocity that satisfies the Taylor’s hypothesis. Moreover, the space-time structure function:

$$S_q\left(\tau ,\lambda \right)=⟨{\left|\mathbf{R}(t+\tau ,\mathbf{x}+\lambda )-\mathbf{R}(t,\mathbf{x})\right|}^q⟩$$E38

is given by the following space-time scaling law:

$$S_q\left(\tau ,\lambda \right)\propto {\lambda }^{K\left(q\right)}$$E39

where the space-time scaling exponent function $K\left(q\right)$ shown here does not dependent on either time or space scales in a direct form. If one considers only the space-time anisotropy due to the large-scale advection velocity, the time domain can be rescaled with the advection velocity in order to maintain the space-time self-similarity property.

In Deidda’s approach [38], there exists a scale-dependent velocity parameter $u_{\lambda }$ which is approximated as ${\lambda }^H$, where $H$ is the anisotropic scaling exponent. If $H=0$, the velocity parameter is constant and independent of the scale $\lambda$, if $H>0$ the velocity parameter increases at larger scales and conversely, if $H<0$ the velocity parameter increases at smaller scales. For a multiplicative random cascade process with branching number $b$, the spatial scale would be given by ${\lambda }_n={\mathit{Lb}}^{-n}$ (for $n=\left[0,\infty \right)$), the scale-dependent velocity parameter would be given by $u_{\lambda }=U_0{b}^{-\mathit{nH}}$, and the timescale would be ${\tau }_n={\lambda }_n{u}_{\lambda }^{-1}=L{U}_0^{-1}{b}^{-n\left(1-H\right)}$.

Deidda [38] also suggest that to estimate the anisotropic scaling, exponent is necessary to satisfy the self-affinity property which is given by:

$$\frac{P\left(b_n\omega \right)}{P\left(\omega \right)}=\frac{P\left(b_n{k}_x\right)}{P\left(k_x\right)}=\frac{P\left(b_n{k}_y\right)}{P\left(k_y\right)}$$E40

defined by the following power laws:

$$P\left(\omega \right)\propto {\omega }^{-s_t};P\left(k_x\right)\propto k_x^{-s_x};P\left(k_y\right)\propto k_y^{-s_y}$$E41

where $P\left(\cdot \right)$ represents the power spectrum, $\omega$ is the frequency, $k_x$ and $k_y$ are the wave numbers in the zonal and meridional directions, respectively, and $s_x$, $s_y$, and $s_t$ represent the slopes of the space-time power spectrum. The zonal anisotropic scaling exponent is then computed as $H_x=1-s_x/s_t$ and the meridional anisotropic scaling exponent as $H_y=1-s_y/s_t$.

On another point of view, the aforementioned self-affinity property can be understood as an space-time transformation such that:

$$\mathbf{x}\to \frac{\mathbf{x}}{b_s};t\to \frac{t}{b_t}$$E42

where $b_s$ identifies the branching number for space and $b_t$ is the branching number for time. Since the space and time branching numbers are related through the anisotropic scaling exponent $H$, the time branching number can be estimated as $b_t=b_s^{1-H}$.

In similar fashion to Deidda [38], Marzan et al. [28] have proposed a space-time multifractal model based on the concept of anisotropic scaling and the existence of a causal process in the space-time dynamics of rainfall (i.e., the past influences the future). In this model, rainfall intensity values $\mathbf{R}\left(t,\mathbf{x}\right)$ represent random spatial structures of scale $\lambda$ and density $\epsilon \left(\mathbf{x}\right)$. For a discrete random cascade model, densities evolve by computing products between random variables as follows:

$${\epsilon }_{n+1}\left(\mathbf{x}\right)={\epsilon }_n\left(\mathbf{x}\right)W_{n+1}\left(\mathbf{x}\right)$$E43

where $W$ is a random random variable that satisfies the scaling law $⟨W^q⟩\approx {\lambda }^{K\left(q\right)}$ and the positivity property ${\mathbb{P}}_r\left\{W<0\right\}=0$. At each step of the cascade, the structures have a scale given by ${\lambda }_n={\mathit{Lb}}^{-n}$, which become structures of scale ${\lambda }_{n+1}={\mathit{Lb}}^{-\left(n+1\right)}$ and so on. In the space-time domain, the evolution of scales is given by ${\lambda }_{n+1}={\mathit{Lb}}^{-\left(n+1\right)}$ for the spatial domain and ${\tau }_{n+1}={\mathit{Tb}}^{-\left(n+1\right)\left(1-H\right)}$ for the time domain. Under this conceptual model, if the spatial branching number equals to $b_s=8$ and the time branching number equals to $b_t=4$, the anisotropic scaling exponent should be equal to $H=1/3$, which is the value expected for passive scalar fields in turbulent environments. Clearly, one can observe that the lifetime of structures depend on the spatial scale and the degree of space-time anisotropy. For $H\to 0$, structures change at the same velocity either time or space, but for $H\to 1$, structures evolve faster in time than in the space. Figure 5 shows a discrete space-time simulation of the rainfall field by following the conceptual model suggested by Marzan et al. [28]. In this figure, the space is represented in a one-dimensional domain and so is the time; the rainfall intensity measures turn out from a random cascade simulation run, which evolves in time according to Eq. 43 and accounting an anisotropy scaling exponent $H\approx 1/3$.

Figure 5.

Example of a discrete space-time random cascade. This construction uses an anisotropic scaling exponent $H\approx 1/3$, a branching number $b=2$ and nine stages of the cascade development. The space and time domains were normalized between 0 and 1; however, the spatial domain is discretized by 1024 bins and the time domain by 128. The generator of the random variable $W$ used in this realization a standard normal distribution. The time changes of the 1-D spatial rainfall field evolve according to Eq. (43).

For generating continuous cascades, the multiplicative process can be replaced by an additive process of the logarithms of field values. Initially, the cascade generator is given by ${\mathrm{\Gamma }}_{\lambda }$ at the scale $\lambda$, then the densities of the field are computed as:

$${\epsilon }_{\lambda }\left(t,\mathbf{x}\right)=e^{{\mathrm{\Gamma }}_{\lambda }\left(t,\mathbf{x}\right)}$$E44

Once again, the field will be multifractal if the scaling law $⟨{\epsilon }_{\lambda }^q⟩\approx {\lambda }^{K\left(q\right)}$ or equivalently, $⟨e^{q{\mathrm{\Gamma }}_{\lambda }}⟩\approx e^{K\left(q\right)log\lambda }$. Besides, the multifractal field generator can be expressed as an impulse-response function between a scale function $g_{\lambda }\left(t,\mathbf{x}\right)\approx \parallel \left(t,\mathbf{x}\right){\parallel }^h$ and a random function ${\gamma }_{\lambda }\left(t,\mathbf{x}\right)$:

$${\mathrm{\Gamma }}_{\lambda }\left(t,\mathbf{x}\right)=g_{\lambda }\left(t,\mathbf{x}\right)\ast {\gamma }_{\lambda }\left(t,\mathbf{x}\right)={\int }_0^t{g}_{\lambda }\left(\tau ,\mathbf{x}\right){\gamma }_{\lambda }\left(t+\tau ,\mathbf{x}\right)\mathit{d\tau }$$E45

where this integration represents a sum of random variables at different scales. The random function ${\gamma }_{\lambda }\left(\mathbf{x},t\right)$ in Eq. (45) represents a causal filter, which can be expressed as follows:

$$g_{\lambda }^{-1}\left(\mathbf{x},t\right)\ast {\mathrm{\Gamma }}_{\lambda }\left(t,\mathbf{x}\right)={\gamma }_{\lambda }\left(\mathbf{x},t\right)$$E46

where $g_{\lambda }\left(\mathbf{x},t\right)$ represents a Green function and that corresponds to a fractional diffusion equation of the form:

$$\left({\partial }_t^{{\xi }_1}-{\partial }_{\mathbf{x}}^{{\xi }_2}\right){\mathrm{\Gamma }}_{\lambda }\left(\mathbf{x},t\right)={\gamma }_{\lambda }\left(\mathbf{x},t\right)$$E47

with,

$${\xi }_1=\frac{h\left(1+1/d\right)}{\left(1-H\right)};{\xi }_2=\frac{h\left(1+1/d\right)}{2};{\widehat{g}}_{\lambda }\left(\omega ,|k|\right)=\frac{1}{{\left|\omega \right|}^{2{\xi }_2}-{\left(i\omega \right)}^{{\xi }_1}}$$E48

where $d$ is the Euclidean spatial dimension and $H$ is the anisotropic scaling exponent. Eq. (47) suggests the necessity of of employing an anisotropic scaling to represent the space-time asymmetry observed in rainfall patterns. In order to describe mathematically the anisotropic scaling, a fractional differential equation is suggested for the multifractal field generator; furthermore, the well-known multiplicative process used for describing spatial fields is changed by an additive process which constitutes a mechanism of modeling random cascades as a fractional diffusive process.

Another approach for the space-time modeling, besides the ones mentioned above, and one which considers explicitly an anisotropy between time and space, was introduced by Over and Gupta [30]. In this approach, the spatial precipitation fields are modeled as random multiplicative cascades that change in time by considering that their cascade generators can be represented as a time-indexed stochastic process. Under such a consideration, spatial precipitation fields evolve by following a evolutionary causal stochastic process. Moreover, the approach is settled in a discrete theory instead that a continuous one since it eases the estimation of nonrainy areas and its application for solving some hydrologic problems.

In the previous section, the Over and Gupta’s approach [16] was presented for describing the construction of random cascade fields without a time indexation. The extension of this approach for describing space-time precipitation fields implies as mentioned above that the cascade generator of $W$ be replaced by an independent and identically distributed stochastic process $\left\{W_t\right\}$ indexed to time, which should satisfy for each time step:

$${\mathbb{P}}_r\left\{W_t\ge 0\right\}=1;\mathbb{E}\left[W_t\right]=1$$E49

The cascade limit mass at time $t$ is defined as ${\mu }_{\infty }\left({\mathrm{\Delta }}_n^i,t\right)$ with a cascade generator $W_t$, which possesses a temporal cross moment defined by:

$$M_n^q\left(t_1,t_2\right)=\sum _{i=1}^{b^n}{\mu }_{\infty }^q\left({\mathrm{\Delta }}_n^i,t_1\right){\mu }_{\infty }^q\left({\mathrm{\Delta }}_n^i,t_2\right)$$E50

Since the field masses are identically distributed and they obeyed the law of large numbers, the ensemble moments are computed as the expectation of the temporal cross moment, thus:

$$\mathbb{E}\left[M_n^q\left(t_1,t_2\right)\right]={\left(\frac{R_0\left(t_1\right)R_0\left(t_2\right)L_0^{2d}}{b^{2n}}\right)}^q{\left(b^2\mathbb{E}\left[W_{t_1}^q{W}_{t_2}^q\right]\right)}^n\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]$$E51

Taking $b$-base logarithms to this equation gives:

$$\begin{array}{l}{log}_b\mathbb{E}\left[M_n^q\left(t_1,t_2\right)\right]=q{log}_b\left(R_0\left(t_1\right)R_0\left(t_2\right)L_0^{2d}\right)+{log}_b\left(\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]\right)\\ +n\left({log}_b\left(\mathbb{E}\left[W_{t_1}^q{W}_{t_2}^q\right]\right)+2\left(1-q\right)\right)\end{array}$$E52

or equivalently:

$$\begin{array}{l}{log}_b\mathbb{E}\left[M_n^q\left(t_1,t_2\right)\right]=q{log}_b\left(R_0\left(t_1\right)R_0\left(t_2\right)L_0^{2d}\right)+{log}_b\left(\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]\right)\\ -2{log}_b\left(\frac{L_n}{L_0}\right)\left({log}_b\left(\mathbb{E}\left[W_{t_1}^q{W}_{t_2}^q\right]\right)+2\left(1-q\right)\right)\end{array}$$E53

From which the temporal cross moment can rewritten as follows:

$${log}_b{M}_n^q\left(t_1,t_2\right)=\frac{{log}_b\left(M_n^q\left(t_1,t_2\right)\right)}{{log}_b\left(\mathbb{E}\left[M_n^q\left(t_1,t_2\right)\right]\right)}+{log}_b\left(\mathbb{E}\left[M_n^q\left(t_1,t_2\right)\right]\right)$$E54

$$={log}_b\left(Y_n^q\left(t_1,t_2\right)\right)+{log}_b\left(\mathbb{E}\left[M_n^q\left(t_1,t_2\right)\right]\right)$$E55

where $Y_n^q\left(t_1,t_2\right)$ is a positive martingale sequence allowing the temporal cross moment to be approximated as follows:

$$\begin{array}{l}{log}_b{M}_n^q\left(t_1,t_2\right)\approx {log}_b\left(Y_n^q\left(t_1,t_2\right)\right)+{log}_b\left(\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]\right)-\\ 2{log}_b\left(\frac{L_n}{L_0}\right){log}_b\left(\mathbb{E}\left[W_{t_1}^q{W}_{t_2}^q\right]\right)\end{array}$$E56

This equation shows that for $L_n/L_0\to 0$, the temporal cross moments of a single space-time realization are scale-invariant. Moreover, the stochastic definition of the random variable $W_t$ determines the structure and variability of the temporal cross moments of the space-time rainfall field. For instance, if $W_t$ is a time stationary stochastic process, then:

$$\begin{array}{l}{log}_b{M}_n^q\left(t_1,t_2\right)\approx {log}_b\left(Y_n^q\left(t_1,t_2\right)\right)+{log}_b\left(\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]\right)-\\ 2{log}_b\left(\frac{L_n}{L_0}\right){log}_b\left(\mathbb{E}\left[W_{t_1}^q\right]\right){log}_b\left(\mathbb{E}\left[W_{t_2}^q\right]\right)\end{array}$$E57

Now, if $\mathbb{E}\left[W_t^q\right]$ possess a power law distribution given by $\mathbb{E}\left[W_t^q\right]\approx {\left(L_n/L_0\right)}^{{\theta }_q}$, then:

$${log}_b{M}_n^q\left(t_1,t_2\right)\approx {log}_b\left(Y_n^q\left(t_1,t_2\right)\right)+{log}_b\left(\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]\right)-2{\theta }_q^2{\left({log}_b\left(\frac{L_n}{L_0}\right)\right)}^3$$E58

where ${\theta }_q$ is a space-time moment scaling exponent. If in addition the time is related to the spatial scale via an anisotropic scaling exponent $H$, as suggested by [28]:

$$\frac{L_n}{L_0}={\left(t_2-t_1\right)}^{1/\left(1-H\right)}$$E59

The temporal cross moment would be given by the following expression:

$${log}_b{M}_n^q\left(t_1,t_2\right)\approx {log}_b\left(Y_n^q\left(t_1,t_2\right)\right)+{log}_b\left(\mathbb{E}\left[Z_{\infty }^q\left(t_1\right)Z_{\infty }^q\left(t_2\right)\right]\right)-2{\theta }_q^2{\left(\frac{{log}_b\left(t_2-t_1\right)}{1-H}\right)}^3$$E60

The last equation suggests that the divergence of temporal cross moments is set for an anisotropic scaling exponent close to unity (i.e., $H\to 1$); hence, a well-defined algebraic structure of the space-time scaling moments requires that $H>1$.

Clearly, Over and Gupta’s approach [30] and the others mentioned above have some limitations to be used as a general theory of the space-time structure of rainfall fields. Certainly, a key aspect in the definition of a suitable rainfall model for a specific region, besides of having an anisotropic scaling exponent, is the elucidation of a Markovian process for the random cascade generator in order to arrive at a causal process. Further works should be addressed to structure an overall Markovian process to represents the space-time nature of rainfall dynamics, hopefully, in connection with a physical atmospheric theory.

5. Summary and conclusions

This chapter highlighted that the multiscaling approach is a landmark theory for the description of the space-time structure of rainfall, crucial for improving the current research and applications in meteorology and hydrological sciences. Clearly, the usage of such schemes in the world depends on the advances in the interpretation and prediction of upcoming climatic events, even more in those regions with scarcity of meteorological observations or where the climatic data are necessarily reconstructions of models outcomes.

One of the ideas highlighted in this chapter is that some complex natural phenomena including turbulence can be described by multiscaling properties. In a similar fashion, the space-time structure of rainfall can be described by multiscaling properties, which enlighten the complex composition of the rainfall field in terms of singular measures, self-similar groups, and power laws. As explained in the chapter, the identification of multiscaling properties in the rainfall fields relies on adopting and adapting the mathematical techniques that were originally applied in turbulence theory. Although rainfall and turbulence are interwoven by similar physical processes, the statistical and geometrical properties vary between them; moreover, the concept of universality in turbulence does not seem to hold on rainfall, because many studies have reported changes in rainfall scaling properties depending on regional or local conditions. The aforementioned observation suggests that rainfall is an inherently complex natural process that requires being understood locally; hence, rainfall models should be designed to specific regional conditions and attached to the available measurement techniques.

As shown here, fractal theories represent the cornerstone in the formulation of a theory for the description of space-time rainfall processes. The predominant characteristics of rainfall statistics are indeed given by power laws, singular scaling exponents, and fractal dimensions for the set of rainfall measures. Our understanding of such characteristics have allowed the birth of several models that are settled in stochastic and geometrical theories. The multiscaling models for rainfall that were mentioned herein can describe space or time patterns with a good resemblance of those ones found in observational fields; however, these models exhibit limitations that prevent their usage under general conditions, specially in their adaptation to the space-time analysis.

In rainfall, time and space models are supposed to be linked via an isotropic scaling exponent, which in turn represents a symmetry breaking of the time-space domain. On the other hand, the most acknowledged space-time multifractal (or multiscaling) conceptualization recognizes the necessity of including a causal process for the description of rainfall field temporal changes. Under the basis of the linear theory, Marzan et al. [28] proposed a fractional diffusive approach for getting a causal process, but the Over and Gupta’s approach [30] suggests the usage of a time-indexed stochastic process. Although these approaches exhibit some conceptual parallelism on their theoretical formulation, there still exits differences that preclude defining a unified theory. Moreover, their applications in general regions require an appropriate knowledge of the space-time scaling structure of the rainfall field in order to correctly define their parameters.

Acknowledgments

The authors thank all the technical support provided by the Institute of Atmospheric Science and Climate Change (ICAyCC)—Universidad Nacional Autonóma de México (UNAM). The first author thanks the General Direction of Academic Staff Affairs (DGAPA—UNAM) for its postdoctoral scholarship.

Conflict of interest

The authors declare no conflict of interest.

Nomenclature

$\alpha$	singularity scaling exponent
${\alpha }_0$	average singularity scaling exponent of a multiplicative cascade
$\beta$	binary number identifier
$b$	branching number
$b_s$	spatial branching number
$b_t$	temporal branching number
${\mathcal{C}}_{\alpha }$	set of $\alpha$-elements
$c\left(\alpha \right)$	codimension function of measures with singularity scaling exponent of order $\alpha$
Cov	covariance operator
$d$	embedding dimension of the Euclidean space
$\mathit{dx}$	probability measure which is generated by the FM approach on the $x$ axis
${\left.\mathit{dx}\right|}_{\theta }$	rotated probability measure resulting from applying a transformation to the measure $\mathit{dx}$ with a rotation angle $\theta$
$\mathit{dy}$	probability measure which is generated by the FM approach on the $y$ axis
${\left.\mathit{dy}\right|}_{\theta }$	rotated probability measure resulting from applying a transformation to the measure $\mathit{dy}$ with a rotation angle $\theta$
$\mathrm{\Delta }$	increment
${\mathrm{\Delta }}_n^i$	$i$-th subdomain after $n$ partitions of a random cascade
${\epsilon }_n$	discrete random variable
$\mathbb{E}$	expectation value
$\epsilon$	rate of kinetic energy dissipation
$f\left(\alpha \right)$	fractal dimension of measures with singularity scaling exponent of order $\alpha$
$g_{\lambda }\left(\mathbf{x},t\right)$	green function/scaling transformation group
${\mathrm{\Gamma }}_{\lambda }$	random cascade generator at the scale $\lambda$
$G$	set identifier
$h$	power spectrum exponent
$H$	self-similarity scaling exponent/anisotropic scaling exponent
$H_x$	zonal anisotropic scaling exponent
$H_y$	meridional anisotropic scaling exponent
$I^k$	numeric interval of order $k$
$\lambda$	spatial or timescale
$K\left(q\right)$	moment scaling exponent function for the moments of order $q$
$\mathrm{\ell }$	spatial longitude
$L_0$	length of the spatial domain
$L_n$	length of the spatial domain
${\lambda }_0$	integral scale
$S_q$	structure function of order $q$
$k$	wave number
$\eta$	Kolmogorov scale
$\mu$	measures (normalized) of rainfall
${\mu }_{\infty }$	limit mass (or measure) of a random cascade
$\mu \left({\mathrm{\Delta }}_n\right)$	low frequency component of a random cascade outcome
$M_{\mu }^q$	moment scaling function of order $q$ for the measures $\mu$
$M$	cumulative measure function
$N$	number of elements in a set/number of contractile affine maps
$g_{\mu }\left(\mu \right)$	probability density function of measures $\mu$
$\xi$	power spectrum exponent
$\xi \left(q\right)$	scaling exponent function for the structure function of order $q$ for the rainfall field
${\theta }_q$	scaling exponent function for the structure function of order $q$ for the rainfall field
$\zeta \left(q\right)$	scaling exponent function for the structure function of order $q$ for the velocity field
$q$	order of the statistical moment
${\mathbb{P}}_r$	probability
$P$	power spectrum
$\mathrm{R}$	Rainfall field
$R_0$	initial rainfall intensity in a random cascade
$\Re$	real numbers
${\Re }^{+}$	real positive numbers
$\theta$	rotation angle for the multifractal measures projection
$\tau$	time scale/temporal lag
$\tau \left(q\right)$	mass exponent scaling function for the moment $q$
$t$	continuous time position
$T$	transpose operator
$x$	spatial position for the standard basis $e_1$
$\mathbf{x}$	spatial vector $\mathbf{x}={\sum }_i^d{x}_i{e}_i$
$\mathbf{u}$	velocity vector $\mathbf{u}={\sum }_i^d{u}_i{e}_i$
$V$	random variable
$W$	random variable/random cascade generator
$W_n$	contractile map
$y$	spatial position for the standard basis $e_2$
$Z_{\infty }$	high-frequency component of a random cascade outcome

Abbreviations

FM	fractal multifractal
pdf	probability density function
SLLN	strong law of large numbers