Abstract
This paper is devoted to a fractional generalization of the Dirichlet distribution. The form of the multivariate distribution is derived assuming that the n partitions of the interval [0, Wn] are independent and identically distributed random variables following the generalized Mittag-Leffler distribution. The expected value and variance of the one-dimensional marginal are derived as well as the form of its probability density function. A related generalized Dirichlet distribution is studied that provides a reasonable approximation for some values of the parameters. The relation between this distribution and other generalizations of the Dirichlet distribution is discussed. Monte Carlo simulations of the one-dimensional marginals for both distributions are presented.
1 Introduction
Let us consider a finite sequence of n positive random variables Z1, ..., Zn. For instance, these variables can represent the wealth of n economic agents if indebtedness is not allowed. Let us denote the sum as Wn = Z1 + … + Zn. In the wealth interpretation this is the total wealth. If we define the wealth fraction of the i-th agent as Qi = Zi/Wn, we get a partition of the interval [0, 1] represented by the sequence Q = (Q1, …, Qn) such that Q1 + … + Qn = 1 almost surely. We are particularly interested in multivariate distributions for the sequence Q whose one-dimensional marginals have heavy tails. If we further assume that the random variables Z1, …, Zn are independent and identically distributed, there is a nice and immediate relationship with point processes of renewal type. In this case, the variables Zi can be interpreted as inter-event intervals and the partial sums Wk =
In order to clarify the relationship, we start by recalling some basic facts on the time-fractional Poisson process as we are going to use and generalize it in the next section. From [1, 14] we know that the time-fractional Poisson process Nν = (Nν(t))t≥0, ν ∈ (0, 1], can be defined as a renewal process with independent and identically distributed inter-event waiting times 𝓣j, j ∈ ℕ* = {1, 2, …}, with probability density function (pdf)
where
is the two-parameter Mittag–Leffler function. Note that for ν = 1, the waiting times 𝓣j are exponentially distributed and N1 is the homogeneous Poisson process. Moreover, the Laplace transform of the pdf (1.1) takes a very compact form. Indeed, we have
Let us now indicate with Tk, k ∈ ℕ*, the random occurrence time of the k-th event of the stream of events defining Nν. From the renewal structure of Nν we readily obtain that the Laplace transform of Tk reads
which in turn corresponds to the Laplace transform of a function involving the three-parameter Mittag–Leffler function (also known as the Prabhakar function – see [7]). In particular, the three-parameter Mittag–Leffler function is defined as
and we know by direct calculation that (see e.g. [15], formula (2.3.24))
where ℜ(α) > 0, ℜ(β) > 0, ℜ(z) > 0, z > |ζ|1/ℜ(α). Using (1.6), we obtain
Remark 1.1
Note that, for ν = 1, the above density reduces to that of an Erlang(λ, k) distributed random variable. This can be seen by simply noticing that
Remark 1.2
The Erlang(λ, k) distribution is a special case of the Gamma(a, c) distribution. Consider a sequence Z1, …, Zn of independent random variables each following a Gamma distribution of parameter (a1, c), …, (an, c). It is well known that their sum Wn is still a Gamma of parameter (a1 + … + an, c). Then the sequence of fractions Q1, …, Qn has a joint (N – 1)-dimensional Dirichlet distribution of parameters a1, …, an with density
with q1 + … + qn = 1 and is independent of Wn.
The proof of the results in Remark 1.2 can be found in several textbooks and lecture notes (see e.g. [2], Lemma 1.5).
Remark 1.3
The random variables 𝓣j have the following asymptotic behaviour for t → ∞ [8]:
therefore, their sums Tk belong to the basin of attraction of the ν-stable subordinator.
Remark 1.4
The distributions considered in the present paper belong to the class of distributions on the simplex discussed in [5] (see (2.2) below and [3]).
This paper contains the following material. Section 2 concerns the definition and properties of the fractional Dirichlet distribution. Section 3 mirrors Section 2 and is devoted to the generalized Dirichlet distribution. Section 4 explains how to simulate the fractional Dirichlet distribution and presents the results of Monte Carlo simulations in order to illustrate the relation between the fractional Dirichlet distribution and the generalized Dirichlet distribution.
2 Construction of the fractional Dirichlet distribution
Based on Remark 1.1 and Remark 1.2, we now define a generalization of the Gamma distribution and we immediately present a fractional generalization of the Dirichlet distribution.
Definition 2.1
(Fractional Gamma distribution). Let X be a positive real valued random variable with distribution
where λ > 0, x > 0, β > 0, ν ∈ (0, 1]. Then X is said to be distributed as a fractional Gamma of parameters λ, β, ν (we write X ∼ FG(λ, β, ν)) (see [19]; for applications to renewal processes see [4, 16, 17, 18]).
Remark 2.1
The Laplace transform of μ reads
By means of (2.1), we will construct a generalization of the Dirichlet distribution. We consider n independent random variables Zi, i = 1, …, n, distributed as fractional Gamma random variables of parameters (1, βi, ν), ν ∈ (0, 1], βi > 0, i = 1, …, n, respectively. Furthermore, define the sum W = Z1 + … + Zn, set Qi = Zi/W, i = 1, …, n, and consider the transformation
vskip-2pt Note that, from (2.2), the distribution of W is fractional Gamma as well, i.e. W ∼ FG(1, β̄, ν), where β̄ =
The joint pdf of Q = (Q1, …, Qn–1) is then obtained by marginalization. Hence,
Remark 2.2
On the n-dimensional simplex Δn the probability density of the random vector (Q1, …, Qn), where ∑nQn = 1 a.s., writes
Notice that for ν = 1 the integral in the rhs of (2.6) can be easily solved and the Dirichlet(β = (β1, …, βn)) is obtained. In this case (Q1, …, Qn) is uniformly distributed on Δn for βi = 1, i ∈ ℕ*.
If ν ∈ (0, 1) with βi = 1, we have
which is symmetric but not uniform.
If we let instead βi = 1/ν (again symmetric), we obtain
2.1 Properties
The derivation of the marginal moments can be done explicitly using the formulas in Section 2.2 of [5].
Proposition 2.1
LetQ = (Q1, …, Qn–1) be a random vector distributed with pdf(2.5). For eachj = 1, …, n – 1, we have,
Proof
By Proposition 2 of [5] we have
Similarly, the second moment writes
and hence after some computation
□
Remark 2.3
Notice that the first factor of the variance (2.10) is in fact the variance of a one-dimensional marginal of a Dirichlet(β) distribution. It follows that the marginals are overdispersed with respect to those of a Dirichlet(β) distribution.
We now proceed by analyzing the aggregation property and therefore the marginal distributions.
Proposition 2.2
(Aggregation property). Consider the pdf defined in equation(2.5)and the random variable 𝓠 =
Proof
The proof is immediate considering that 𝓠 comes from the sum of i.i.d. positive random variables each one with Laplace transform given by (2.2) and then divided by W. Therefore 𝓩 = W 𝓠 has pdf given by (2.1) with β =
An immediate corollary of this result is
Corollary 2.1
(Marginal pdf). Consider the pdf in equation(2.5). Then its marginal onQiis given by
As the three-parameter Mittag-Leffler function has a representation as an H-function [15],
for suitable choices of (ai, Ai) and (bi, Bi), the marginal pdf (2.14) can be expressed in terms of an H-function too.
Proposition 2.3
IfQiis the random variable with pdf(2.14), then
Proof
In the integral (2.14) set
Denote with
For δ > 0 and ν ∈ (0, 1] we have
For η ∈ (0, β̄), by using (2.19) and Theorem 2.9 in [12], we have
Set η = β̄–ε in (2.20) with ε ∈ (0, β̄), and use (2.17) to recover
Remark 2.4
By using Properties 2.1, 2.3 and 2.5 of [12], the H-function in (2.16) can be rewritten interchanging βi with β̄–βi and qi with 1–qi, which corresponds to commuting the two Mittag-Leffler functions in (2.14).
According to Theorems 1.3 and 1.4 in [12], since
the H-function in (2.16) has a power series expansion. The following propositions rely on this property.
Proposition 2.4
Forqi < 1/2 andβinot a positive integer
where
Proof
Consider the H-function
The claim follows replacing a1 = 1 – βi, a2 = 1 – β̄ + ε, a3 = ν (ε – βi), b2 = ε – βi and b3 = 1 – νβi in (2.24) and taking the limit as ε ↓ 0.□
Remark 2.5
If ν (βi + k) is not a positive integer for k = 0, 1, 2, … thanks to the reflection formula for the gamma function [12], we might simplify the expansion in (2.22) using
Similarly we get Γ(–νk) Γ(νk) = –π/(kν sin(πνk)) if νk is not a positive integer.
Proposition 2.5
Forqi > 1/2 andβ̄ – βinot a positive integer
where
Proof
Consider again the H-function
The claim follows replacing a1 = 1 – βi, a2 = 1 – β̄ + ε, a3 = ν (ε – βi), b2 = ε – βi and b3 = 1 – νβi in (2.28) and taking the limit as ε ↓ 0.□
By using the reflection formula for the gamma function, also the expansion (2.26) might be simplified similarly to what has been addressed in Remark 2.5.
3 An alternative generalization
We now give an alternative generalization with desirable properties which in addition approximates the fractional Dirichlet distribution with density (2.5) for appropriate values of the parameters.
Let us thus consider a random vector Q = (Q1, …, Qn–1), n ≥ 2, with the following probability density function:
for q1, …, qn–1 ∈ (0, 1), q1 + … + qn–1 < 1, ν > 0, βi > 0, i = 1, …, n, β̄ = β1 + … + βn.
For the sake of clarity we check that fQ(q1, …, qn–1) as given in (3.1) is a genuine probability density function. This will follow by proving that
in the rhs of (3.1) plays the role of a normalization coefficient.
Theorem 3.1
We have
Proof
Observe that the lhs of (3.3) can be rewritten as
where q̄ = q1 + ⋯ + qn–1. Apply the change of variables (1 – q̄)/qi = zi for i = 1, …, n – 1, in multivariate integration. Thus, we have
where 𝓘n–1,k = {1, …, k – 1, k + 1, …, n – 1} for k = 1, …, n – 1, and J is the Jacobian of the transformation. Note that
By putting (3.5) and (3.6) in (3.4) we have
Apply the change of variables
where
Observe that In–2(t1, …, tn–2) in (3.9) can be rewritten as
With the change of variable
Now, replace In–2(t1, …, tn–2) in (3.8) with the closed form (3.11). This leads us to
where
By comparing the integral in (3.13) with that in (3.10), we observe that the former has the same expression of the latter with β̄ replaced by β̄ – βn–1. Thus, by recurring to the same arguments employed to compute In–2(t1, …, tn–2) we recover
Replacing In–3(t1, …, tn–3) in (3.12) with the closed form (3.14) we get
where
which indeed has the same expression of In–2 and In–3 with suitable updates of β̄. The result follows by iterating from i = 4 up to i = n – 1 the computation of
with
We obtain the closed form expression
with
from which the claimed result follows by observing that
□
Remark 3.1
Alternatively, in (3.7) use the transformation
which is in agreement with (3.3).
Remark 3.2
On the n-dimensional simplex Δn the pdf of the random vector (Q1, …, Qn), where ∑nQn = 1 a.s., writes
with
Notice that for ν = 1 the Dirichlet(β) is obtained. In this case the random vector (Q1, …, Qn) is uniformly distributed on Δn for βi = 1, i ∈ ℕ*. If instead we only let βi = 1,
which is symmetric but clearly not uniform. If βi = 1/ν (again symmetric) we obtain
Remark 3.3
The alternative generalized Dirichlet distribution considered in this section (i.e. that with pdf (3.1)) can be derived by the same procedure described in Section 2 with (Zi)ν distributed as Gamma(βi, 1), i = 1, …, n – 1. Note that the random variable X such that Xν, ν > 0, is Gamma(α, 1)-distributed, α > 0, is a special case of the generalized Gamma distribution (see e.g. [11], Section 8.7). In particular, X has pdf
and Laplace transform (from (2.3.23) of [15] and the definition of Wright functions)
Remark 3.4
The generalized Dirichlet pdf (3.1) turns out to be a reasonably good approximation of the fractional Dirichlet pdf (2.5) for βi < 1 (see for example Fig. 3). A partial explanation is that for λ = 1, β < 1, and ν ∈ (0, 1] the fractional Gamma pdf (2.1) has a rather similar shape to the generalized Gamma pdf (3.25), as Fig. 1 shows for β = 0.2 and β = 0.4. For βi > 1, the fractional Dirichlet pdf exhibits a behaviour different from the generalized Dirichlet pdf (see for example Fig. 2). Indeed, Fig. 1 shows a different shape of the fractional Gamma pdf compared to the generalized Gamma pdf for β = 2 and β = 3.
Pdf’s for ν = 0.7, β1 = 2, β2 = 3. The black circles represent the results of a Monte Carlo simulation for the fractional Dirichlet distribution and the blue line is the corresponding theoretical value from equation (2.5). The red circles come from a Monte Carlo simulation of the generalized Dirichlet distribution and the green curve is the plot of the theoretical pdf.
Proposition 3.1
(Conjugate distribution). The generalized Dirichlet distribution GDIR(ν, β) (with pdf(3.22)) is the conjugate prior to a re-parametrized Multinomial distribution with pmf
whereN ∈ ℕ+, xi ∈ {0, …, N}, i = 1, …, n, n ∈ ℕ+, q1 + … + qn = 1, ν > 0. In particular, if the prior is GDIR(ν, β) and the likelihood is as in(3.27), then the posterior becomes GDIR(ν, β + x).
Proof
The proof is a straightforward application of Bayes theorem. The reparametrization in (3.27) is such that
3.1 Representation in terms of Dirichlet random variables
In order to derive a meaningful representation in terms of Dirichlet random variables for the random vector Q, we first recall the definitions of two related classes of random vectors (see [10]).
Definition 3.1
(Liouville distribution of the first kind). Let X = (X1, …, Xn) be an absolutely continuous random vector supported on the n-dimensional positive orthant, i.e. Rn = {(x1, …, xn) : xi > 0 for each i = 1, …, n}. It is said to have Liouville distribution of the first kind if its joint pdf writes
where ai > 0, i = 1, …, n, and f is a positive continuous function satisfying ∫ℝ+ya–1f(y) dy < ∞, with a = a1, …, an. Further, we write
Definition 3.2
(Liouville distribution of the second kind). Let Z = (Z1, …, Zn) be an absolutely continuous random vector supported on Sn = {(z1, …, zn) : zi > 0 for each i = 1, …, n,
where ci > 0, i = 1, …, n, and g is a positive continuous function satisfying ∫ℝ+yc–1g(y) dy < ∞, with c = c1, …, cn. Further, we write Z∼
Remark 3.5
If we let f(t) = (1 + t)–(a+an+1), t > 0, an+1 > 0, in (3.28), then X is distributed as an inverted Dirichlet.
If, in (3.29), we choose g(t) = (1 – t)cn+1–1, 0 < t < 1, an+1 > 0, we have that Z is distributed as a Dirichlet.
Proposition 3.1 of [10] tells us what is the relationship between Liouville distributions of the first and of the second kind (and hence between the Dirichlet and the inverted Dirichlet). Specifically, if Z ∼
then X ∼
Plainly, the converse relation is true as well: inverting (3.31) (letting h = t/(1 + t)) we have
As a simple example, considering f(t) = (1 + t)–(c+cn+1), t > 0 (inverted Dirichlet), we readily obtain g(1 – h)cn+1–1 (Dirichlet).
Now, by exploiting the above definition we prove the following distributional representation for Q.
Proposition 3.2
LetQ = (Q1, …, Qn–1) be distributed with pdf(3.1). Then the random vectorM = (M1, …, Mn–1) such that
is distributed as a Dirichlet(β = (β1, …, βn)).
Conversely, ifM ∼ Dirichlet(β) we have thatQ = (Q1, …, Qn–1) such that
is distributed with pdf(3.1).
Proof
Let us define the random vector Y = (Y1, …, Yn–1) such that
and let β* = β1 + …+ βn–1. Combining the transformations in the proof of Theorem 3.3 and of Remark 3.1 we see that Y has pdf
and hence Y ∼
By using the mentioned transformations between Liouville distributions of first and second kinds we have that
M ∼
and (3.33) follows.
Finally, a rewriting of the components of Q in terms of those of Y, leads easily to (3.34).□
4 Monte Carlo simulations
The simulation of the random variables Q for the fractional Dirichlet distribution is straightforward based on the construction presented in Section 2. First, one needs to generate random variables with density (2.1) and one can use the mixture representation discussed in [4]
where Ui is Gamma(βi, λ)-distributed and Vν is strictly positive-stable distributed with exp(–sν) as the Laplace transform of the pdf. Summing the Xi to get W and dividing Xi by W gives Qi.
Remark 4.1
For β = 1, there is an alternative representation [13, 6]:
where Ξ is Exp(λ)-distributed and Z is Cauchy-distributed.
The behaviour of the fractional Dirichlet distribution in the case N = 2 is shown in Fig. 2 for ν = 0.7, β1 = 2, β2 = 3. In this case, the generalized Dirichlet distribution is not a good approximation.
This is not the case for N = 2, ν = 0.7, β1 = 0.2 and β2 = 0.4 where the generalized Dirichlet distribution is a reasonably good approximation of the fractional Dirichlet distribution. This is represented in Fig. 3.
Pdf’s for ν = 0.7, β1 = 0.2, β2 = 0.4. Circles and curves have the same meaning as in Fig. 2.
For larger values of the parameters βi, one gets a unimodal distribution in both cases as shown in Fig. 4 for N = 2, ν = 0.95, β1 = 10, β2 = 30.
Pdf’s for ν = 0.95, β1 = 10, β2 = 30. Circles and curves have the same meaning as in Fig. 2.
In Fig. 5, the heavy character of the right tail of the generalized Dirichlet distribution is highlighted by the log-log plot.
Double logarithmic plot of the pdf’s for ν = 0.95, β1 = 10, β2 = 30. Circles and curves have the same meaning as in Fig. 2.
Acknowledgements
F. Polito has been partially supported by the project “Memory in Evolving Graphs” (Compagnia di San Paolo/Università degli Studi di Torino).
References
[1] L. Beghin, E. Orsingher, Fractional Poisson processes and related planar random motions. Electron. J. Probab. 14, No 61 (2009), 1790–1827.10.1214/EJP.v14-675Search in Google Scholar
[2] J. Bertoin, Exchangeable Coalescents. ETH Zurich, 2010.Search in Google Scholar
[3] L. Bondesson, A general result on infinite divisibility. Ann. Probab. 7, No 6 (1979), 965–979.10.1214/aop/1176994890Search in Google Scholar
[4] D. O. Cahoy, F. Polito, Renewal processes based on generalized Mittag-Leffler waiting times. Comm. Nonlinear Sci. Numer. Simulat. 18, No 3 (2013), 639–650.10.1016/j.cnsns.2012.08.013Search in Google Scholar
[5] S. Favaro, G. Hadjicharalambous, I. Prünster, On a class of distributions on the simplex. J. Statist. Plann. Inference141, No 9 (2011), 2987–3004.10.1016/j.jspi.2011.03.015Search in Google Scholar
[6] D. Fulger, E. Scalas, G. Germano, Monte Carlo simulation of uncoupled continuous-time random walks yielding a stochastic solution of the space-time fractional diffusion equation. Phys. Rev. E77 (2008), Art. 021122.10.1103/PhysRevE.77.021122Search in Google Scholar
[7] A. Giusti, I. Colombaro, R. Garra, R. Garrappa, F. Polito, M. Popolizio, F. Mainardi, A practical guide to Prabhakar fractional calculus. Fract. Calc. Appl. Anal. 23, No 1 (2020), 9–54; 10.1515/fca-2020-0002; https://www.degruyter.com/view/journals/fca/23/1/fca.23.issue-1.xml.Search in Google Scholar
[8] R. Gorenflo, F. Mainardi, Fractional Calculus: Integral and Differential Equations of Fractional Order. In: Fractals and Fractional Calculus in Continuum Mechanics, A. Carpinteri and F. Mainardi (Eds.), Springer, New York and Wien (1997), 223–276.10.1007/978-3-7091-2664-6_5Search in Google Scholar
[9] I. S. Gradshteyn, I. M. Ryzhik, Table of Integrals, Series, and Products. 8th Edition, Elsevier/Academic Press, Amsterdam (2015).Search in Google Scholar
[10] R. D. Gupta, D. S. P. Richards, Multivariate Liouville distributions. J. Multivariate Anal. 23, No 2 (1987), 233–256.10.1016/0047-259X(87)90155-2Search in Google Scholar
[11] N. L. Johnson, S. Kotz, N. Balakrishnan, Continuous Univariate Distributions. Vol. 1. 2nd Edition, Wiley Ser. in Probability and Mathematical Statistics: Applied Probability and Statistics, John Wiley & Sons, Inc., New York (1994).Search in Google Scholar
[12] A. A. Kilbas, M. Saigo, H-Transforms. Theory and Applications. Ser. Analytical Methods and Special Functions, Vol. 9, Chapman & Hall/CRC, Boca Raton, FL (2004).Search in Google Scholar
[13] T. Kozubowski, Mixture representation of Linnik distributions revisited. Stat. Probab. Lett. 38 (1998), 157–160.10.1016/S0167-7152(97)00167-3Search in Google Scholar
[14] F. Mainardi, R. Gorenflo, E. Scalas, A fractional generalization of the Poisson processes. Vietnam J. Math. 32 (Special Issue) (2004), 53–64.Search in Google Scholar
[15] A. M. Mathai, H. J. Haubold, Special Functions for Applied Scientists. Springer, New York (2008).10.1007/978-0-387-75894-7Search in Google Scholar
[16] T. M. Michelitsch, A. P. Riascos, Generalized fractional Poisson process and related stochastic dynamics. Fract. Calc. Appl. Anal. 23, No 3 (2020), 656–693; 10.1515/fca-2020-0034; https://www.degruyter.com/view/journals/fca/23/3/fca.23.issue-3.xml.Search in Google Scholar
[17] T. M. Michelitsch, A. P. Riascos, Continuous time random walk and diffusion with generalized fractional Poisson process. Phys. A545 (2020), Art. 123294.10.1016/j.physa.2019.123294Search in Google Scholar
[18] T. M. Michelitsch, F. Polito, A. P. Riascos, Biased continuous-time random walks with Mittag-Leffler jumps. Fractal and Fractional4, No 4 (2020), Art. 51.10.3390/fractalfract4040051Search in Google Scholar
[19] R. N. Pillai, On Mittag-Leffler functions and related distributions. Ann. Inst. Statist. Math. 42, No 1 (1990), 157–161.10.1007/BF00050786Search in Google Scholar
© 2021 Diogenes Co., Sofia
