The compound product distribution; a solution to the distributional equation X=AX+1

TL;DR

This paper analyzes the solution of the distributional equation X=AX+1 for discrete A, presenting a fast algorithm to compute the heavy tail density, with detailed study of the compound product distribution for specific cases.

Contribution

It introduces a fast algorithm for calculating the heavy tail density of solutions to the distributional equation X=AX+1, with detailed analysis of the compound product distribution.

Findings

Developed a computationally efficient method for heavy tail density calculation.

Provided detailed analysis for specific distribution families.

Demonstrated the applicability to compound product distributions.

Abstract

The solution of $ X=AX+1 $ is analyzed for a discrete variable $ A $ with $ \mathbb{P}\left[A=0\right]>0 $. Accordingly, a fast algorithm is presented to calculate the obtained heavy tail density. To exemplify, the compound product distribution is studied in detail for some particular families of distributions.

Tables

Table 1. Table 1: Conditions for mean and variance to be finite for the CP distribution for different distributions for the variable A 𝐴 A .

$A\sim$	Conditions for	Conditions for
$\mathrm{Po}\left(\lambda \right)$
$\mathrm{Bin}(n,p)$	always for $n=1$	always for $n=1$
	, $n\ge 2$	, $n\ge 2$
$\mathrm{NB}(r,p)$		, $r=1$
		, $r\ge 2$
$\mathrm{Geo}\left(p\right)$

Table 2. Table 2: Mean and variance of the CP distribution for different distributions of A 𝐴 A . In each case the obtain formulas are equal to ∞ \infty if the denominator are equal or less to zero. Table   1 present the detailed conditions for the existence of the first and second moments.

$A\sim$	$𝔼\left[X\right]$	$\mathrm{Var}\left(X\right)$
$\mathrm{Po}\left(\lambda \right)$	$\frac{1}{1-\lambda }$	$-\frac{\lambda }{{\left(\lambda -1\right)}^2\left({\lambda }^2+\lambda -1\right)}$
$\mathrm{Bin}(n,p)$	$\frac{1}{1-np}$	$\frac{n\left(p-1\right)p}{{\left(np-1\right)}^2\left(n^2{p}^2-n\left(p-1\right)p-1\right)}$
$\mathrm{NB}(r,p)$	$\frac{p}{pr+p-r}$	$\frac{\left(p-1\right)p^{2}r}{{\left(pr+p-r\right)}^2\left(p^2\left(r^2-1\right)-pr\left(2r+1\right)+r\left(r+1\right)\right)}$
$\mathrm{Geo}\left(p\right)$	$\frac{p}{2p-1}$	$-\frac{\left(p-1\right)p^2}{{\left(1-2p\right)}^2\left(3p-2\right)}$

Table 3. Table 3: Moment estimators for the CP distribution. For simplicity of notation we denote by μ ′ = 1 n ​ ∑ i = 1 n X i superscript 𝜇 ′ 1 𝑛 superscript subscript 𝑖 1 𝑛 subscript 𝑋 𝑖 \mu^{\prime}=\frac{1}{n}\sum_{i=1}^{n}X_{i} and μ 2 ′ = 1 n ​ ∑ i = 1 n X i subscript superscript 𝜇 ′ 2 1 𝑛 superscript subscript 𝑖 1 𝑛 subscript 𝑋 𝑖 \mu^{\prime}_{2}=\frac{1}{n}\sum_{i=1}^{n}X_{i} , the first and second population moments, respectively.

$A\sim$	Moment matching estimators
$\mathrm{Po}\left(\lambda \right)$	$\hat{\lambda }=1-\frac{1}{{\mu }^{\prime }}$
$\mathrm{Bin}(n,p)$	$\hat{n}=\frac{{\left({\mu }^{\prime }-1\right)}^2{\mu }_2^{\prime }}{\left(2{\mu }^{\prime }-1\right){\mu }^{\prime 2}+\left({\mu }^{\prime 2}-3{\mu }^{\prime }+1\right){\mu }_2^{\prime }}$	$\hat{p}=\frac{\left(2{\mu }^{\prime }-1\right){\mu }^{\prime 2}+\left({\mu }^{\prime 2}-3{\mu }^{\prime }+1\right){\mu }_2^{\prime }}{\left({\mu }^{\prime }-1\right){\mu }^{\prime }{\mu }_2^{\prime }}$
$\mathrm{NB}(r,p)$	$\hat{r}=-\frac{{\left({\mu }^{\prime }-1\right)}^2{\mu }_2^{\prime }}{\left(2{\mu }^{\prime }-1\right){\mu }^{\prime 2}+\left({\mu }^{\prime 2}-3{\mu }^{\prime }+1\right){\mu }_2^{\prime }}$	$\hat{p}=-\frac{\left({\mu }^{\prime }-1\right){\mu }^{\prime }{\mu }_2^{\prime }}{\left(2{\mu }^{\prime }-1\right)\left({\mu }^{\prime 2}-{\mu }_2^{\prime }\right)}$
$\mathrm{Geo}\left(p\right)$	$\hat{p}=\frac{{\mu }^{\prime }}{2{\mu }^{\prime }-1}$

Table 4. Table 4: Estimated values and Akaike information criterion, for the Moby Dick data.

$A\sim$	Estimators		AIC
$\mathrm{Po}\left(\lambda \right)$	$\hat{\lambda }=0.48$		174343.18
$\mathrm{Bin}(n,p)$	$\hat{n}=2$	$\hat{p}=0.37$	142566.88
$\mathrm{NB}(r,p)$	$\hat{r}=10$	$\hat{p}=0.94$	157936.02
$\mathrm{Geo}\left(p\right)$	$\hat{p}=0.67$		173737.18

Table 5. Table 5: Skewness of the infinite product distribution for different A 𝐴 A densities. In each case the obtain formulas are equal to ∞ \infty if the denominator are equal or less to zero.

$A\sim$	$\mathrm{Skewness}\left(X\right)$
$\mathrm{Po}\left(\lambda \right)$	$\frac{{\left(\lambda -1\right)}^3{\left({\lambda }^2+\lambda -1\right)}^2\left(5{\lambda }^2+2\lambda +1\right)}{{\lambda }^2\left({\lambda }^3+3{\lambda }^2+\lambda -1\right)}$
$\mathrm{Bin}(n,p)$	$-(({(np-1)}^3{(-n^2{p}^2+n(p-1)p+1)}^2$
	$4(n-1)n{p}^3+n(6-5n)p^2-2(n-1)p-1)/$
	$n^2{(p-1)}^2{p}^2(n^3{p}^3-3n^2(p-1)p^2+n(2p^2-3p+1)p-1))$
$\mathrm{NB}(r,p)$	$({(pr+p-r)}^3{(p^2(r^2-1)-pr(2r+1)+r(r+1))}^2$
	$(p{r}^2(-4r^2+r-5)+r(r^3+r+2)+p^3(-4r^4+3r^3-3r^2+3r+1)+$
	$p^4(r^4-r^3+r-1)+p^2(6r^4-3r^3+7r^2-6r+1)))/$
	${\left(p-1\right)}^2{p}^3{r}^2\left(p^3\left(r^2-r+1\right)+p^{2}r\left(2-3r\right)+3p{r}^2-r\left(r+1\right)\right)$
$\mathrm{Geo}\left(p\right)$	$\frac{{\left(2-3p\right)}^2{\left(2p-1\right)}^3\left(6p^2-13p+8\right)}{{\left(p-1\right)}^2{p}^3\left(2p^3-7p^2+12p-6\right)}$

Table 6. Table 6: Kurtosis of the infinite product distribution for different A 𝐴 A densities. In each case the obtain formulas are equal to ∞ \infty if the denominator are equal or less to zero.

$A\sim$	$\mathrm{Kurtosis}\left[X\right]$
$\mathrm{Po}\left(\lambda \right)$	$\frac{{\left(\lambda -1\right)}^4{\left({\lambda }^2+\lambda -1\right)}^3\left(3{\lambda }^6-25{\lambda }^5-55{\lambda }^4-32{\lambda }^3-47{\lambda }^2-11\lambda -1\right)}{{\lambda }^3{\left(\lambda +1\right)}^2\left({\lambda }^2+2\lambda -1\right)\left({\lambda }^3+5{\lambda }^2+2\lambda -1\right)}$
$\mathrm{Bin}(n,p)$	$({(np-1)}^4{(n^2{p}^2-n(p-1)p-1)}^3$
	$(1+(11n-6)p+(47n^2-65n+6)p^2+n(32n^2-165n+138)p^3+$
	$n\left(55n^3-177n^2+243n-120\right)p^4+n\left(25n^4-184n^3+339n^2-216n+36\right)p^5-$
	$3{(n-1)}^2{n}^2(n^2+4n-12)p^7+3n^2(n^4+10n^3-62n^2+93n-42)p^6))/$
	$(n^3{(p-1)}^3{p}^3(n^3{p}^3-3n^2(p-1)p^2+n(2p^2-3p+1)p-1)$
	$n^4{p}^4-6n^3(p-1)p^3+n^2(11p^2-18p+7)p^2+n(-6p^4+12p^3-7p^2+p)-1)$
$\mathrm{NB}(r,p)$	$-(({(pr+p-r)}^4{(p^2(r^2-1)-pr(2r+1)+r(r+1))}^3$
	$(r^2{(r+1)}^2(3r^3+5r^2-3r+6)-2p{r}^2(12r^5+32r^4+24r^3+16r^2+15r+3)+$
	$p^{3}r\left(-168r^6-106r^5-140r^4-29r^3+76r^2+35r+10\right)+$
	$p^{6}r\left(84r^6-211r^5+238r^4-126r^3+98r^2-38r-24\right)+$
	$p^{2}r\left(84r^6+139r^5+98r^4+68r^3-2r^2-18r+6\right)+$
	$p^{4}r\left(210r^6-85r^5+210r^4-103r^3-9r^2-66r-19\right)+$
	$p^7\left(-24r^7+86r^6-108r^5+39r^4+6r^3-33r^2+30r+1\right)+$
	$p^8\left(3r^7-14r^6+20r^5-4r^4-16r^3+20r^2-7r-2\right)-$
	$2p^5(84r^7-122r^6+140r^5-91r^4+66r^3-50r^2-2r-1)))/$
	$({(p-1)}^3{p}^4{r}^3(p^3(r^2-r+1)+p^{2}r(2-3r)+3p{r}^2-r(r+1))$
	$(p^{3}r(-4r^2+6r-3)+r(r^2+3r+2)-$
	$p^4(r^3-3r^2+2r-1)+p^2(6r^3+r)+2pr(2r^2+3r+1))))$
$\mathrm{Geo}\left(p\right)$	$\frac{{\left(1-2p\right)}^4{\left(3p-2\right)}^3\left(42p^6-173p^5+105p^4+435p^3-872p^2+642p-180\right)}{{\left(p-1\right)}^3{p}^4\left(2p^3-7p^2+12p-6\right)\left(15p^3-50p^2+60p-24\right)}$

Equations

$$X_t=A_t{X}_{t-1}+B_t,t\in \mathbb{Z},$$

$$X_t=Z_t(\beta +\lambda X_{t-1}^2{)}^{1/2},t\in \mathbb{Z},$$

$$X_t=A_t{X}_{t-1}+1,$$

$$X\stackrel{d}{=}AX+1.$$

$$X\stackrel{d}{=}1+A_1+A_1{A}_2+A_1{A}_2{A}_3+\dots ,$$

$$\mathbb{E}\left[X^m\right]=\frac{{\sum }_{i=0}^{m-1}\left(\genfrac{}{}{0ex}{}{m}{i}\right)\mathbb{E}\left[A^i\right]\mathbb{E}\left[X^i\right]}{1-\mathbb{E}\left[A^m\right]}.$$

$$\mathbb{E}\left[X^m\right]=\sum _{i=0}^m\left(\genfrac{}{}{0ex}{}{m}{i}\right)\mathbb{E}\left[A^i\right]\mathbb{E}\left[X^i\right],$$

$$X=1+\sum _{i=1}^N{A}_1^{\prime }{A}_2^{\prime }\dots A_i^{\prime },$$

$$\mathbb{P}\left[X=n\right]$$

$$=\sum _{k=0}^{n-1}\mathbb{P}\left[\sum _{i=1}^N{A}_1^{\prime }\dots A_i^{\prime }=n-1|N=k\right]p_0(1-p_0{)}^{k-1},$$

$$\mathbb{P}\left[\sum _{i=1}^N{A}_1^{\prime }\dots A_i^{\prime }=n-1|N=k\right]=0,k=n,n+1,\dots .$$

$$n=1+{\xi }_1+{\xi }_1{\xi }_2+{\xi }_1{\xi }_2{\xi }_3+\dots +{\xi }_1\cdots {\xi }_{n-1}.$$

$$\mathbb{P}\left[X=n\right]=\sum _{(i,j)\in {\Delta }_n}{p}_i\mathbb{P}\left[X=j\right],n\ge 1,$$

$$\omega (n)\sim \log \log n\text{and}\Omega (n)\sim \log \log n.$$

$$\hat{p}=\frac{1}{\hat{n}{\mu }^{\prime }}$$

$$\hat{p}=\hat{r}\frac{\mu \left((2\mu -1){\mu }^2+\left({\mu }^2-3\mu +1\right){\mu }_2\right)}{(\mu -1)(2\mu -1)\left({\mu }^2-{\mu }_2\right)}.$$

Full text

K****eywords Power laws $\cdot$ Random coefficient autoregressive model $\cdot$ Stochastic difference equation $\cdot$ Heavy tail distributions

1 Introduction

The Kesten’s stochastic recurrent equation

[TABLE]

where $(A_{t},B_{t})$ are independent of $X$, has many interesting properties that makes it appealing for applications, Kesten (1973). This framework allows flexibility in describing and estimating the conditional distribution including, for instance, fat tails and stationary behaviors. For example, de Haan et al. (1989) analyze the ARCH model

[TABLE]

as a solution of (1) with $A_{t}=\lambda Z_{t}^{2}$, $B_{t}=\beta Z_{t}^{2}$. Other uses of this structure range from processes with renewal epochs to model growth-collapse behaviors (see Boxma et al. (2011)) to limit theorems to analyze TCP congestion (see Dumas et al. (2002) and Guillemin et al. (2004)). A recent monograph that unites the results of this equation is Buraczewski et al. (2016). Probably the most remarkable property about the solutions to Kesten’s equation is that, under general assumptions, the tails of $X_{t}$ are of power-law-type even for light-tailed input variables $(A_{t},B_{t})$ (see Kesten (1973)). This implies that these models had many interesting applications on regular variation, weak convergence of probability measures and point processes, as is pointed out in Buraczewski et al. (2016).

The main contribution of the present work is a recursive algorithm to calculate the general density of the solution of (1), when $B=1$ and the support of $A$ is restricted to the non-negative integers with $\mathbb{P}\left[A=0\right]>0$. Here and subsequently these restrictions over the density of $A$ are assumed. In other words, in this manuscript it is analyzed the probability structure of $X_{t}$ for the recursive equation

[TABLE]

which is equivalent to the analysis of the distributional equality

[TABLE]

Under this framework the fractal structure that reside in (3) is rich enough to allow many probability structures and, at the same time, does not complicate their use to model the behavior of real heavy tailed data. Moreover, the methodology presented here could be used as a new method to obtain heavy tail distributions.

The remaining part of this document is organized as follows: In Section 2 it is considered some general properties of the density of $X$ and also it is studied its probability structure under four different families of densities for $A$. Section 3 establishes a computational algorithm to calculate its density. The behavior of the method of moments and maximum likelihood estimators are presented in Section 4. Furthermore, in this section a fitting exercise with real data is presented. Final points and conclusions are deferred to Section 5. Additionally, in the Appendix are presented the expressions of skewness and kurtosis, for the four studied distributions.

2 The compound product distribution

The existence of a solution of (3) is given in Corollary 2.1.2 and Theorem 2.1.3 of Buraczewski et al. (2016). These results allows rewriting $X$ as

[TABLE]

for a sequence $\left\{A_{n}\right\}_{n\in\mathbb{N}}$ of independent identically distributed (i.i.d.) random variables. This expression is the reason to say that $X$ has a compound product (CP) distribution. Mena and Walker (2012) also consider the behavior of sums of products to define a probability model and, as in this contribution, they obtain a recursive formula to evaluate the model density. In order to limit the summation in (4) it is necessary to put some restrictions on the density of $A$. In particular, the $\mathbb{P}\left[A=0\right]>0$ assumption implies that with probability one only a finite number of the terms of the series $\{A_{1}\cdots A_{n}\}_{n\geq 1}$ are non-zero. This condition is not particularly restrictive and allows to set many distributional behaviors for $X$, by changing the distribution of $A$. The next result implies that a subtle shift over the moments of $A$ gets amplified in the density of $X$.

Theorem 1**.**

Let $X$ be the solution of (3). If $\mathbb{E}\left[A^{m}\right]<1$ for $m\in\mathbb{N}$, then

[TABLE]

Proof.

The binomial expansion applied to (3) implies

[TABLE]

and (5) follows. ∎

Table 1 presents the conditions for finitude of the first and second moments and Table 2 presents the mean and variance, for the CP distribution under the four families of distributions. The important point to note from these two tables is the heavy tail of the distribution of $X$. For instance, to have a finite variance under $A\sim\mathrm{Po}(\lambda)$ we need $\lambda<(\sqrt{5}-1)/2\approx 0.62$. These tables exemplify how heavy tailed the CP distribution could be.

2.1 Examples of families of compound product densities

This section is focus on four different families of distributions for $A$: Poisson ($\mathrm{Po}(\lambda)$, with mean $\lambda$), binomial ($\mathrm{Bin}(n,p)$, with mean $np$), negative binomial ($\mathrm{NB}(r,p)$, with mean $pr/(1-p)$) and geometric ($\mathrm{Geo}(p)$, with mean $p/(1-p)$) . Although the geometric distribution is a particular case of the negative binomial family, its analysis expose some properties that could be obscure under the negative binomial framework. An example of the versatility of the negative binomial distribution is presented in Leisen et al. (2019), to model the stationary time series behavior. On the contrary, an analysis of the Bernoulli distribution it is not worthy since in this case $X$ is a geometric random variable.

It is possible to study the CP density as the result of a geometric compound process with dependent terms. An equivalent representation of the CP distribution is given as the distribution of the random variable $X$, given as

[TABLE]

where $N\sim\mathrm{geo}(p_{0})$ and $\mathbb{P}\left[A_{j}^{\prime}=n\right]=p_{n}/(1-p_{0})$, with $p_{n}:=\mathbb{P}\left[A=n\right]$. This expression is obtain by conditioning (4) to the random variable $N=\inf\{i\in\mathbb{N}:A_{i}=0\}$. One advantage of using this representation is that it gives an straight way to calculate the density of $X$. By conditioning over $N$ and applying (6), it is obtained that $\mathbb{P}\left[X=0\right]=p_{0}$, and for $n\geq 1$,

[TABLE]

where the last equality is given by the fact that

[TABLE]

The problem with (7) is that the number of combinations of $A_{1}^{\prime},A_{2}^{\prime},\ldots,A_{N}^{\prime}$ that fulfills $\sum_{i=1}^{N}A_{1}^{\prime}\ldots A_{i}^{\prime}=n-1$ is computational expensive to calculate, since is not available a closed form expression. To solve this issue an recursive algorithm to compute $\mathbb{P}\left[X=n\right]$ is presented in Section 3. Figures 1 and 2 present the changes of the CP density for different values of the parameters. Figure 1(a) shows the change in the Poisson density for values $\lambda<1$ (finite mean case) and Figure 1(b) for values $\lambda>1$ (infinite mean case). An equivalent analysis is presented in Figure 2 for the geometric distribution.

3 An algorithm to calculate the compound product density

In this section is presented a computational algorithm to calculate $\mathbb{P}\left[X=n\right]$ for $n\geq 1$. The problem to calculate these probabilities centers on the next question: Which are the different vectors $(\xi_{1},\ldots,\xi_{n-1})\in\mathbb{N}^{n-1}$ that fulfill the next equation?

[TABLE]

To confront this problem the algorithm here presented focuses around two strategies: prime decomposition and recursion. To explain how the algorithm works, let us first focus on the decomposition of a natural number into two factors. Let us assume that $n-1$ has the prime decomposition $n-1=z_{1}^{a_{1}}\cdots z_{r}^{a_{r}}$. This implies that the set $\Delta_{n}=\{(i,j)\in\mathbb{N}^{2}:i*j=n-1\}$ has cardinality equal to the number of ways to choose how many factors goes to the number $i$; this cardinality is then equal to $(a_{1}+1)(a_{2}+1)\cdots(a_{r}+1)$. By conditioning on the values of $A$, we obtain the recursive equality

[TABLE]

since $X$ and $A$ are independent. Under these observations is obtained the Algorithm 1, which returns all the values $\mathbb{P}\left[X=i\right]$ for $i=1,\ldots,n$, for a fixed $n$. Let us mention some particularities of this algorithm. Using the notation $\omega(n)$ for the number of different prime factors of $n$, and $\Omega(n)$ for the total number of prime factors of $n$, it is well known (see Hardy et al. (1927)) that these quantities behave asymptotically as

[TABLE]

This two asymptotic behaviors contributes to the speed of the algorithm. Even though Algorithm 1 has a loop nested on other loop, it is fast since (10) give each loop a few terms. These observations are depicted in Figure 3. This figure compares the computational time of $\mathbb{P}\left[X=n\right]$ by direct calculation and by Algorithm 1. To obtain the direct calculation it is applied (7). As a consequence of (10), the time consumption of this algorithm is so suitable that it is almost lineal. For instance, the algorithm takes one second for each thousand values added to $n$ from $1$ to $10^{7}$, computed on Intel i7 processor. Also, it is important to mention that with the current computer advances the calculation of the prime decomposition is not an issue of time or memory consumption for this algorithm.

4 Estimation of parameters

4.1 Estimators with the method of moments

Table 3 presents the estimators of method of moments for the CP distribution under the Poisson, binomial, negative binomial and geometric distributions. The existence of these estimators depend on the existence of their respective moments. Consequently, these estimators behave badly when the true parameters are near to conditions that make the correspondent moments infinite (see Table 1).

The confidence intervals for the moment estimators of the Poisson and geometric distributions are presented in Figures 4(a) and 4(b), respectively. In these figures the black line represents the true value of the estimators and the gray area is its confidence interval from 5% to 95%, obtained through simulation of 1000 variables for each parameter value. In both cases the confidence interval is slim since their relative error is small.

In the case of the binomial and negative binomial densities, it is important to notice that extra restrictions must be applied to confine the values of $\hat{n}$ and $\hat{r}$ to the positive integers. For instance, to set these restrictions in the binomial case, one could solve the method of moments equations for the closest integer next $\hat{n}$ and then substitute this value in $\hat{p}$ using the relation

[TABLE]

Similarly, in the negative binomial case

[TABLE]

However, these extra restrictions hamper the accuracy of these estimators.

4.2 Maximum likelihood estimators

Using Algorithm 1 one could obtain the MLE estimators, since it is computationally fast enough to compute them. Like with the method of moments, the estimators become inferior as the value of the true parameters tends to the frontiers of a finite mean distribution. In contrast to the method of moments, the MLE method do not have problems to estimate $\hat{n}$ and $\hat{r}$, for the binomial and negative binomial distributions. Moreover, the restricting of the searching region to the positive integers speeds-up its computation.

Figures 5(a) and 5(b) display a sample density of the MLEs for the binomial and negative binomial distributions, respectively; these densities were obtained through the simulation of 1000 samples of sizes 100 for each point estimator for the fixed true value. In both cases, the densities center around the true value of the distribution. Furthermore, the errors of the estimators $\hat{n}$ and $\hat{r}$ were at most one unit apart for the true values 90% of the time, for $n,r\in\{2,3,4,5,6,7,8\}$ and $p\in(0,1)$.

4.3 Application to real data

To present an example of how to fit the CP distribution to real data, here are fitted the four distributional models to the Moby Dick word count data (see Alsmeyer and Löwe (2013)). This data is the words occur in the novel Moby Dick, which is well-known to have a heavy-tail behavior. The Poisson and geometric distributions are fitted using the method of moments, and the binomial and negative binomial using the MLE method. To select among the four model it is calculated the Akaike information criterion. The results of the estimation for the Moby Dick data are presented in Table 4. Under this criterion the best model is the $\mathrm{Bin}(2,.37)$, which has finite mean but infinite variance.

5 Conclusions

In this manuscript is presented a fast algorithm to calculate the general density of the solution of the distributional equation $X=AX+1$, for a discrete random variable $A$ with $\mathbb{P}\left[A=0\right]>0$. To exhibit the behavior of this distribution, four commonly used parametric families of distributions for $A$ are particularly studied. The results here presented show that practitioners can use this methodology to generate new heavy tailed distributions. The evidence of the real data experiment suggests that the estimation method works correctly for the these four families of distributions.

6 Acknowledgments

The author gratefully acknowledges the helpful suggestions of Prof. Begoña Fernández during the preparation of the manuscript. This research was supported by a DGAPA-UNAM Posdoctoral Scholarship.

Appendix