A New Class of Non-Central Dirichlet Distributions

Carlo Orsi¹¹1e-mail: \mail[email protected]
Ph.D. in Statistics at University of Milano-Bicocca (Milan, Italy),
Casella Postale 13, Ufficio Postale Fermo, 63900 Fermo (Italy)

Abstract

In the present paper new light is shed on the non-central extensions of the Dirichlet distribution. Due to several probabilistic and inferential properties and to the easiness of parameter interpretation, the Dirichlet distribution proves the most well-known and widespread model on the unitary simplex. However, despite its many good features, such distribution is inadequate for modeling the data portions next to the vertices of the support due to the strictness of the limiting values of its joint density. To replace this gap, a new class of distributions, called Conditional Non-central Dirichlet, is presented herein. This new model stands out for being a more easily tractable version of the existing Non-central Dirichlet distribution which maintains the ability of this latter to capture the tails of the data by allowing its own density to have arbitrary positive and finite limits.

Keywords: Unitary simplex; finite limits of the density; mixture representation; perturbation representation; simulation; computational efficiency.
MSC 2010 subject classifications: 62E15

1 Introduction

Thanks to its simplicity and good mathematical properties, the Dirichlet distribution has historically represented the first tool for modeling data on the unitary simplex, this latter being the subset of $\mathbb{R}^{\,D}$ , $D\geq 1$ , given by

\mathcal{S}^{\,D}=\left\{\underline{x}=\left(x_{1},\,\ldots\,,x_{D}\right)\,:\;0<x_{i}<1\,,\;i=1,\,\ldots\,,D\;,\;\sum_{i=1}^{D}x_{i}<1\right\}

(1)

and corresponding to the Real interval $(0,1)$ for $D=1$ . Such data typically consist of components which represent parts of a whole and naturally arise in many areas of application such as medical, biological and social sciences. In this regard, a comprehensive review of the great variety of disciplines in which the above distribution has been employed is provided for example by [Ng et al. (2011)]. Despite its several statistical properties fundamental for the analysis of this type of data, one of the weakest points of the Dirichlet family lies in the poorness of its parametrization. Amongst others, this fact implies low flexibility of the density at the vertices of Eq. (1). More precisely, only the $i$ -th component of the parameter vector is devoted to modeling the limiting value of the Dirichlet density at the $i$ -th vertex of the simplex, this latter being the Real vector whose elements are all equal to zero except for the $i$ -th element which is one. This characteristic feature necessarily restricts the applicative potential of such distribution. In fact, the Dirichlet density tipically performs the identification of the tails of the data in an unsatisfactory way by showing limits equal to 0 or $+\infty$ in the cases that positive finite limits would be more appropriate instead.

Several generalizations of the Dirichlet distribution have been proposed in the literature. However, only a few of them enable to overcome the aforementioned limitation of the Dirichlet by allowing their densities to take on arbitrary positive and finite limits which are potentially apt to properly model the data portions lying in the neighborhood of the vertices of the support. Specifically, the bivariate extension of the Kummer-Beta distribution [Bran-Cardona et al. (2011)] and the Non-central Dirichlet distribution [Sánchez et al. (2006)] are right worthy of consideration in the above sense due to the major flexibility achieved by the limiting values of their densities by means of richer parametrizations. As a matter of fact, due to its further parameter, the density of the former model shows a more flexible behavior than the Dirichlet at the vertices of the unit simplex, whereas the density of the latter, despite its analytical hardness and poor interpretability, enables to fully overcome the above mentioned drawback of the Dirichlet thanks to its additional vector of non-centrality parameters.

That said, the present paper sheds new light on the variety of the probabilistic models on the unitary simplex which accommodate for a more sophisticated structure for the limits of the density by proposing a new non-central generalization of the Dirichlet distribution. This novel model, called Conditional Non-central Dirichlet, can be considered as a more tractable analogue of the existing Non-central Dirichlet model. In fact, the former preserves the ability of the latter to potentially capture the tails of the data by letting its density have arbitrary positive and finite limits and, at the same time, offers the advantage of having a density function which is more straightforward and easily handeable than the existing non-central one.

More precisely, our work is organized as follows. Section 2 provides the mathematical tool kit, which essentially hinges on the notions of ascending factorial and generalized hypergeometric function. A glance at the unit simplex distributions involved in the definition and the analysis of the new model is given as well; particular attention is focused on the limiting values of the densities of such models on varying the parameter vector. An in-depth study of the novel model is carried out in Section 3. Specifically, this latter is properly specified in Subsection 3.1, where definition, distribution and theoretical properties are discussed in detail. In this regard, the Conditional Non-central Dirichlet density is proved to exhibit the same complex mixture type structure as the standard non-central one; however, from its perturbation representation, the greater tractability and interpretability of this new density surprisingly emerge. Then, an investigation of the mixing distribution of the mixture type form of the new Non-central Dirichlet density is reported in Subsection 3.2. Evidences of the similarities between the new and the standard non-central models are provided in Subsections 3.3 and 3.4, where some representations of the new model are presented and significant plots of its density are displayed, respectively. In particular, a special focus is given to the case of unitary shape parameters, under which the density of the new model shows the attractive feature of taking on arbitrary positive and finite limits at the vertices of the unit simplex. Subsection 3.5 derives a closed-form expression for the mixed raw moments, whereas the identifiability of the model is tackled in Subsection 3.6. In Section 4 an application to a real data set highlights the potential of the new model with respect to the aforementioned alternative ones. In conclusion, Section 5 contains some final remarks and Section 6 provides an appendix which bears an implementation in the programming enviromnment R of the multiple infinite sum involved in the perturbation representation of the standard Non-central Dirichlet density.

2 Preliminaries

This section is intended to briefly review the probabilistic models involved in the present study. However, first of all we shall provide the mathematical tool kit that paves the way for a thorough analysis of the new distribution we are going to define. In this regard, it is essential remembering that

\left(a\right)_{l}=\frac{\Gamma\left(a+l\right)}{\Gamma\left(a\right)}=\left\{\begin{array}[]{ll}1&\mbox{ if }l=0\\ \\ a\left(a+1\right)\,\ldots\,\left(a+l-1\right)&\mbox{ if }l\in\mathbb{N}\end{array}\right.

(2)

is the $l$ -th ascending factorial or Pochhammer’s symbol of $a>0$ [Johnson et al. (2005)], where $l\in\mathbb{N}$ indicates the number of terms and $\Gamma\left(\cdot\right)$ denotes the gamma function. Some interesting properties of the above notion are listed in the following. Firstly, by additively decomposing $l$ into $l_{1},\,l_{2}\in\mathbb{N}$ , in light of Eq. (2) one has that

\left(a\right)_{l_{1}+l_{2}}=\frac{\Gamma\left(a+l_{1}+l_{2}\right)}{\Gamma\left(a\right)}=\left\{\begin{array}[]{l}\frac{\Gamma\left(a+l_{1}\right)}{\Gamma\left(a\right)}\,\frac{\Gamma\left(a+l_{1}+l_{2}\right)}{\Gamma\left(a+l_{1}\right)}=\left(a\right)_{l_{1}}\,\left(a+l_{1}\right)_{l_{2}}\\ \\ \frac{\Gamma\left(a+l_{2}\right)}{\Gamma\left(a\right)}\,\frac{\Gamma\left(a+l_{2}+l_{1}\right)}{\Gamma\left(a+l_{2}\right)}=\left(a\right)_{l_{2}}\,\left(a+l_{2}\right)_{l_{1}}\end{array}\right.\,.

(3)

Then, the ratio of two ascending factorials of $a$ can be expressed as

\frac{\left(a\right)_{l_{1}}}{\left(a\right)_{l_{2}}}=\left\{\begin{array}[]{ll}\left(a+l_{2}\right)_{l_{1}-l_{2}}&\mbox{ if }l_{1}\geq l_{2}\\ \\ \left[\left(a+l_{1}\right)_{l_{2}-l_{1}}\right]^{-1}&\mbox{ if }l_{1}<l_{2}\end{array}\right.\,,

(4)

whereas the following expansion holds true for the Pochhammer’s symbol of a binomial:

\left(a+b\right)_{l}=\sum_{j=0}^{l}{l\choose j}\left(a\right)_{l-j}\left(b\right)_{j}\,.

(5)

By virtue of the foregoing concept, the generalized hypergeometric function with $p$ upper parameters and $q$ lower parameters, $p,\,q\in\mathbb{N}\cup\{0\}$ , can be accordingly defined as

{}_{p}^{\,}F_{q}\left(a_{1},\,\ldots\,,a_{p};b_{1},\,\ldots\,,b_{q};x\right)=\sum_{i=0}^{+\infty}\frac{(a_{1})_{i}\,\ldots\,(a_{p})_{i}}{(b_{1})_{i}\,\ldots\,(b_{q})_{i}}\frac{x^{i}}{i\,!}\,,\quad x\in\mathbb{R}\,;

(6)

for more details on the convergence of the hypergeometric series in Eq. (6) as well as for further results and properties of ${}_{p}^{\,}F_{q}\,$ , the reader is recommended referring to [Srivastava and Karlsson (1985)]. In this regard, the special case of Eq. (6) corresponding to $p=0$ and $q=1$ , namely

{}_{0}^{\,}F_{1}\left(b;x\right)=\sum_{i=0}^{+\infty}\frac{1}{(b)_{i}}\frac{x^{i}}{i\,!}\,,\quad x\in\mathbb{R}\,,

(7)

plays a prominent role in the analysis of the new model at study. Hence, some features of this latter that make it an extremely regular function are emphasised herein. Specifically, the function in Eq. (7) is increasing in $x\geq 0$ with values in $\left[1,+\infty\right)$ for any given $b>0$ and decreasing in $b$ for any given $x\geq 0$ as well as its first derivative with respect to $x$ . Finally, another important special case of Eq. (6) is

{}_{1}^{\,}F_{1}\left(a;b;x\right)=\sum_{i=0}^{+\infty}\frac{(a)_{i}}{(b)_{i}}\frac{x^{i}}{i\,!}\,,\quad x\in\mathbb{R}\,,

(8)

which is otherwise called Kummer’s confluent hypergeometric function; this latter satisfies the following transformation:

{}_{1}^{\,}F_{1}\left(a;b;x\right)=e^{x}\,_{1}{}^{\,}F_{1}\left(b-a;b;-x\right)\,,

(9)

which is known as Kummer’s First Theorem.

That said, denoting independence by $\bot$ , we remember that the $D$ -dimensional Dirichlet model with vector of shape parameters $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D+1}\right)$ , $\alpha_{i}>0$ , $i=1,\,\ldots\,,D+1$ , denoted by $\mbox{Dir}^{\,D}\left(\underline{\alpha}\right)$ , is defined as follows:

	$\displaystyle\left\{\begin{array}[]{l}Y_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\,2}_{2\alpha_{i}}\quad i=1,\,\ldots\,,D+1\\ \\ Y^{+}=\sum_{i=1}^{D+1}Y_{i}\end{array}\right.\quad\Rightarrow$				(13)
		$\displaystyle\Rightarrow$	$\displaystyle\quad\underline{X}=\left(X_{1},\,\ldots\,,X_{D}\right)=\left(\frac{Y_{1}}{Y^{+}},\,\ldots\,,\frac{Y_{D}}{Y^{+}}\right)\sim\mbox{Dir}^{\,D}\left(\alpha_{1},\,\ldots\,,\alpha_{D+1}\right)$		(14)

[Kotz et al. (2000)] and its joint density function is given by

\mbox{Dir}^{\,D}\left(\underline{x};\underline{\alpha}\right)=\frac{\Gamma\left(\alpha^{+}\right)}{\prod_{i=1}^{D+1}\Gamma\left(\alpha_{i}\right)}\left[\prod_{i=1}^{D}x_{i}^{\alpha_{i}-1}\right]\left(1-\sum_{i=1}^{D}x_{i}\right)^{\alpha_{D+1}-1}\,,\qquad\underline{x}\in\mathcal{S}^{\,D}\,,

(15)

where $\alpha^{+}=\sum_{i=1}^{D+1}\alpha_{i}$ and $\mathcal{S}^{\,D}$ is specified in Eq. (1). For the ends of this paper it is useful bearing in mind that, in the notation of Eq. (14), the Dirichlet distribution can be also obtained as conditional distribution of $\underline{X}$ given $Y^{+}$ , the former being independent of the latter by virtue of a characterizing property of independent Gamma random variables [Lukacs (1955)]. The relevance of this property is remarkable; therefore, it is made explicit in the following.

Property 2.1 ([Lukacs (1955)] Characterizing property of independent Gamma random variables)

Let $Y_{1}$ , $Y_{2}$ be two nondegenerate and positive random variables and suppose that they are independently distributed. The random variables $Y^{+}=Y_{1}+Y_{2}$ and $X_{i}=Y_{i}\,/\,Y^{+}$ , for every $i=1,2$ , are independently distributed if and only if both $Y_{1}$ and $Y_{2}$ has Gamma distributions with the same scale parameter.

Clearly, Property 2.1 is valid also in the case that the number of the involved independent Gamma random variables is greater than two. In particular, by setting the common scale parameter of these latter equal to $1/2$ , such property holds true also for any finite number of independent Chi-Squared random variables.

Now let $D$ equal 2. The Dirichlet density, despite the great variety of shapes reachable by it on varying the shape parameters, shows poor flexibility at the vertices of the unit simplex. In this regard, the limiting values of the density of $\left(X_{1},X_{2}\right)\sim\mbox{{Dir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3}\right)$ are as follows:

$\displaystyle\mbox{{Dir}}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3}\right)\quad\rightarrow$			(21)
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(1,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}=\alpha_{3}=1\,:&\alpha_{1}\left(\alpha_{1}+1\right)\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}>2\,:&0\end{array}\right.$	(21)
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,1\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{3}=1\,:&\alpha_{2}\left(\alpha_{2}+1\right)\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}>2\,:&0\end{array}\right.$	(27)
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{2}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{2}=1\,:&\alpha_{3}\left(\alpha_{3}+1\right)\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}>2\,:&0\end{array}\right.\,,$	(33)

where the symbol $\not\exists\,$ indicates that the corresponding limit does not exist because depends on the direction followed to reach the accumulation point under consideration. From Eq. (33) it is noticeable that the limiting value at the $i$ -th vertex of the simplex, when different from 0 or $+\infty$ , depends only on the $i$ -th shape parameter under unitary values for the remaining ones. Such state gets worse if $\alpha_{i}=1$ for every $i=1,2,3$ ; in fact, in this latter case the bivariate Dirichlet density reduces to the Uniform on $\mathcal{S}^{2}$ and all its limiting values equal 2. Finally, by Eq. (2), the mixed raw moment of order $(r_{1},r_{2})$ of the $\mbox{{Dir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3}\right)$ distribution can be stated as

\mathbb{E}\left(X_{1}^{\,r_{1}}\,X_{2}^{\,r_{2}}\right)=\frac{\left(\alpha_{1}\right)_{r_{1}}\,\left(\alpha_{2}\right)_{r_{2}}}{\left(\alpha^{+}\right)_{r^{+}}}\,,\qquad\left\{\begin{array}[]{l}r_{1},\,r_{2}\in\mathbb{N}\\ \\ r^{+}=r_{1}+r_{2}\end{array}\right.\,.

(34)

The bivariate Kummer-Beta distribution is the first generalization of the Dirichlet taken into account herein. This model is briefly recalled in the following for what is of interest in the present analysis; however, for further results and properties, the reader can be referred to [Bran-Cardona et al. (2011)]. Specifically, a bidimensional random vector is said to follow a bivariate Kummer-Beta distribution with shape parameters $\alpha_{1}$ , $\alpha_{2}$ , $\alpha_{3}$ and additional parameter $\delta\in\mathbb{R}$ , denoted by $\mbox{KB}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)$ , if its joint density can be expressed as the following infinite mixture of Dirichlet densities weighted by the probabilities obtained by normalizing the terms of the infinite sum in Eq. (8) where $a=\alpha_{3}$ , $b=\alpha^{+}$ , $x=\delta$ :

	$\displaystyle\mbox{KB}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)=\quad$				(35)
		$\displaystyle=$	$\displaystyle\sum_{j=0}^{+\infty}\frac{\frac{\left(\alpha_{3}\right)_{j}}{\left(\alpha^{+}\right)_{j}}\frac{\delta^{j}}{j!}}{{}_{1}F_{1}\left(\alpha_{3};\alpha^{+};\delta\right)}\;\mbox{Dir}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3}+j\right)\,,\qquad\left(x_{1},x_{2}\right)\in\mathcal{S}^{\,2}\,.$		(35)

By simple computations, Eq. (35) can be equivalently stated in terms of the following perturbation of the bivariate Dirichlet density:

\mbox{KB}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)=\mbox{Dir}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3}\right)\,\frac{e^{\,\left(1-x_{1}-x_{2}\right)\,\delta}}{{}_{1}F_{1}\left(\alpha_{3};\alpha^{+};\delta\right)}\,,\quad\left(x_{1},x_{2}\right)\in\mathcal{S}^{\,2}\,,

which, in turn, by Eq. (9), can be rewritten in the final form reported by the above cited reference:

	$\displaystyle\mbox{KB}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)=\quad$				(36)
		$\displaystyle=$	$\displaystyle\mbox{Dir}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3}\right)\,\frac{e^{\,-\left(x_{1}+x_{2}\right)\,\delta}}{{}_{1}F_{1}\left(\alpha_{1}+\alpha_{2};\alpha^{+};-\delta\right)}\,,\quad\left(x_{1},x_{2}\right)\in\mathcal{S}^{\,2}\,.$		(36)

Under unitary shape parameters the behavior of the KB² density at the vertices of the unit simplex turns out to be more flexible than that of the bivariate Dirichlet thanks to its additional parameter $\delta$ . More precisely, the limits of the KB² density take on the following values:

$\displaystyle\mbox{{KB}}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)\quad\rightarrow$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(1,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}=\alpha_{3}=1\,:&\frac{\alpha_{1}\left(\alpha_{1}+1\right)\,e^{-\delta}}{{}_{1}F_{1}\left(\alpha_{1}+1;\alpha_{1}+2;-\delta\right)}\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}>2\,:&0\end{array}\right.$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,1\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{3}=1\,:&\frac{\alpha_{2}\left(\alpha_{2}+1\right)\,e^{-\delta}}{{}_{1}F_{1}\left(1+\alpha_{2};2+\alpha_{2};-\delta\right)}\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}>2\,:&0\end{array}\right.$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{2}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{2}=1\,:&\frac{\alpha_{3}\left(\alpha_{3}+1\right)}{{}_{1}F_{1}\left(2;2+\alpha_{3};-\delta\right)}\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}>2\,:&0\end{array}\right.\,,$

where the symbol $\not\exists\,$ indicates that the corresponding limit does not exist because depends on the direction followed to reach the accumulation point under consideration.

The non-central extension of the Chi-Squared model [Johnson et al. (1995)] represents the main ingredient for the definition and the analysis of the existing Non-central Dirichlet distribution, this latter being the second generalization of the Dirichlet considered in the present paper. Specifically, a Non-central Chi-Squared random variable $Y^{\prime}$ with $g>0$ degrees of freedom and non-centrality parameter $\lambda\geq 0$ , denoted by $\chi^{\prime\,2}_{g}\left(\lambda\right)$ , can be characterized by means of the following mixture representation:

Y^{\prime}\sim\chi^{\prime\,2}_{g}\left(\lambda\right)\qquad\Leftrightarrow\qquad Y^{\prime}\,|\,M\;\sim\;\chi^{\,2}_{g+2M}\,,\qquad\mbox{where}\;\;M\sim\mbox{Poisson}\left(\lambda/2\right)\,,

(43)

the case $\lambda=0$ corresponding to the $\chi^{\,2}_{g}$ distribution. Moreover, such a random variable can be additively decomposed into two independent parts, a central one with $g$ degrees of freedom and a purely non-central one with non-centrality parameter $\lambda$ , namely

Y^{\prime}=Y+\sum_{j=1}^{M}F_{j}\,,\quad\mbox{\normalsize{where}}\;\;\;Y\sim\chi^{\,2}_{g}\quad\bot\quad M\sim\mbox{{Poisson}}\left(\lambda/2\right)\quad\bot\quad\{F_{j}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\,2}_{2}\}\,.

(44)

By virtue of Eq. (44), the random variable $Y^{\prime}_{pnc}=\sum_{j=1}^{M}F_{j}$ is said to have a Purely Non-central Chi-Squared distribution with non-centrality parameter $\lambda$ . Indeed, we shall denote it by $\chi^{\prime\,2}_{0}\left(\lambda\right)$ , its number of degrees of freedom being equal to zero [Siegel (1979)]. By Eq. (43), the density function $f_{Y^{\prime}}$ of $Y^{\prime}\sim\chi^{\prime\,2}_{g}\left(\lambda\right)$ can be expressed as

f_{Y^{\prime}}\left(y;g,\lambda\right)=\sum_{i=0}^{+\infty}\frac{e^{-\frac{\lambda}{2}}\left(\frac{\lambda}{2}\right)^{i}}{i!}\frac{y^{\frac{g+2i}{2}-1}\,e^{-\frac{y}{2}}}{\Gamma\left(\frac{g+2i}{2}\right)2^{\frac{g+2i}{2}}},\quad y>0\,,

(45)

i.e. as the infinite series of the $\chi^{2}_{g+2i}$ densities, $i\in\mathbb{N}\cup\{0\}$ , weighted by the probabilities of $M\sim\mbox{{Poisson}}\left(\lambda/2\right)$ . In this regard, the case $g=2$ is of prominent interest in the present setting; in fact, in this case the density in Eq. (45) exhibits a flexible limiting behavior on varying the non-centrality parameter by allowing its limit at $0$ to take on the following expression:

\lim_{y\rightarrow 0^{+}}f_{Y^{\prime}}\left(y;2,\lambda\right)=\frac{1}{2}\,e^{-\frac{\lambda}{2}}\,.

(46)

Finally, the Non-central Chi-Squared distribution is reproductive with respect to both the number of degrees of freedom and the non-centrality parameter; specifically:

Y^{\prime}_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\prime\,2}_{g_{i}}(\lambda_{i})\quad i=1,\,\ldots\,,m\qquad\Rightarrow\qquad\left\{\begin{array}[]{l}Y^{\prime+}=\sum_{i=1}^{m}Y^{\prime}_{i}\sim\chi^{\prime\,2}_{g^{+}}(\lambda^{+})\\ \\ g^{+}=\sum_{i=1}^{m}g_{i}\,,\;\lambda^{+}=\sum_{i=1}^{m}\lambda_{i}\end{array}.\right.

(47)

That said, the $D$ -dimensional Non-central Dirichlet model with vector of shape parameters $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D+1}\right)$ and vector of non-centrality parameters $\underline{\lambda}=\left(\lambda_{1},\,\ldots\,,\lambda_{D+1}\right)$ , $\lambda_{i}\geq 0$ , $i=1,\,\ldots\,,D+1$ , denoted by $\mbox{NcDir}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ , can be easily defined by replacing the $Y_{i}\,$ ’s by $Y^{\prime}_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\prime\,2}_{2\alpha_{i}}\left(\lambda_{i}\right)$ in Eq. (14) as follows:

	$\displaystyle\left\{\begin{array}[]{l}Y^{\prime}_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\prime\,2}_{2\alpha_{i}}\left(\lambda_{i}\right)\quad i=1,\,\ldots\,,D+1\\ \\ Y^{\prime+}=\sum_{i=1}^{D+1}Y^{\prime}_{i}\end{array}\right.\Rightarrow$				(51)
		$\displaystyle\Rightarrow$	$\displaystyle\quad\underline{X}^{\prime}=\left(X^{\prime}_{1},\,\ldots\,,X^{\prime}_{D}\right)=\left(\frac{Y^{\prime}_{1}}{Y^{\prime+}},\,\ldots\,,\frac{Y^{\prime}_{D}}{Y^{\prime+}}\right)\sim\mbox{NcDir}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$		(52)

[Sánchez et al. (2006)]. The density function of the $\mbox{NcDir}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ distribution can be simply derived by following the next arguments; in this regard, let:

	$\displaystyle\underline{M}=\left(M_{1},\,\ldots\,,M_{D+1}\right)\sim\mbox{{Multi-Poisson}}^{\,D+1}(\underline{\lambda}\,/\,2)\qquad\Leftrightarrow$				(53)
		$\displaystyle\qquad\Leftrightarrow\qquad$	$\displaystyle M_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\mbox{{Poisson}}(\lambda_{i}\,/\,2)\,,\;\quad i=1,\,\ldots\,,D+1\;.$		(53)

Then, in the notation of Eq. (52), by Eq. (43) the conditional distribution of $Y^{\prime}_{i}$ given $\underline{M}$ is of the type $\chi^{\,2}_{2\alpha_{i}+2M_{i}}$ , $i=1,\,\ldots\,,D+1$ ; therefore, by Eq. (14):

\underline{X}^{\prime}\sim\mbox{NcDir}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)\qquad\Leftrightarrow\qquad\left\{\begin{array}[]{l}\left.\underline{X}^{\prime}\,\right|\,\underline{M}\sim\mbox{Dir}^{\,D}\left(\underline{\alpha}+\underline{M}\right)\\ \\ \underline{M}\sim\mbox{{Multi-Poisson}}^{\,D+1}(\underline{\lambda}\,/\,2)\end{array}\right.\,,

(54)

which is the mixture representation of the NcDir distribution. Hence, the joint density of $\underline{X}^{\prime}\sim\mbox{NcDir}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ can be stated as

\mbox{NcDir}^{\,D}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)=\sum_{\underline{j}\in\mathbb{N}_{0}^{\,D+1}}\left[\Pr\left(\underline{M}=\underline{j}\right)\cdot\mbox{Dir}^{\,D}\left(\underline{x};\underline{\alpha}+\underline{j}\right)\right]\,,\qquad\underline{x}\in\mathcal{S}^{\,D}\,,

(55)

	$\displaystyle\mbox{NcDir}^{\,D}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)=\mbox{Dir}^{\,D}\left(\underline{x};\underline{\alpha}\right)\quad\cdot$				(56)
		$\displaystyle\cdot$	$\displaystyle e^{-\frac{\lambda^{+}}{2}}\,\Psi_{2}^{\,(D+1)}\left[\alpha^{+};\underline{\alpha};\frac{\lambda_{1}}{2}x_{1},\,\ldots\,,\frac{\lambda_{D}}{2}x_{D},\frac{\lambda_{D+1}}{2}\left(1-\sum_{i=1}^{D}x_{i}\right)\right]\,,\qquad\underline{x}\in\mathcal{S}^{\,D}\qquad$		(56)

where

\Psi_{2}^{\,(m)}\left[a;b_{1},\,\ldots\,,b_{m};x_{1},\,\ldots\,,x_{m}\right]=\sum_{j_{1},\,\ldots\,,\,j_{m}=\,0}^{+\infty}\frac{(a)_{j_{1}+\,\ldots\,+j_{m}}}{(b_{1})_{j_{1}}\,\ldots\,(b_{m})_{j_{m}}}\,\frac{x_{1}^{\,j_{1}}}{j_{1}!}\,\ldots\,\frac{x_{m}^{\,j_{m}}}{j_{m}!}

(57)

is the $m$ -dimensional ( $m>2$ ) generalization of the Humbert’s confluent hypergeometric function $\Psi_{2}$ [Srivastava and Karlsson (1985)]. Unfortunately, the perturbation representation of the Non-central Dirichlet density in Eq. (56) highlights the uneasy tractability and interpretability of this latter. Indeed, regardless of the constant term, the Dirichlet density is perturbed by a function in $D+1$ variables given by the sum of the multiple power series in Eq. (57), that, as far as we know, cannot be reduced to a more easily tractable analytical form. In this regard, the code of a routine implemented in the programming environment R for the computation of the function in Eq. (57) with $m=3$ is provided in the appendix contained in Section 6.
Now let $D$ equal 2. By simple computations it is easy to see that, despite its poor tractability from a mathematical standpoint, the limiting values of the NcDir² density interestingly take on the following simple expressions:

$\displaystyle\mbox{{NcDir}}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)\quad\rightarrow$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(1,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}=\alpha_{3}=1\,:&e^{\,-\frac{\lambda_{2}+\lambda_{3}}{2}}\cdot\\ &\left[\left(\frac{\lambda_{1}}{2}\right)^{2}+\left(\alpha_{1}+1\right)\left(\alpha_{1}+\lambda_{1}\right)\right]\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}>2\,:&0\end{array}\right.$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,1\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{3}=1\,:&e^{\,-\frac{\lambda_{1}+\lambda_{3}}{2}}\cdot\\ &\left[\left(\frac{\lambda_{2}}{2}\right)^{2}+\left(\alpha_{2}+1\right)\left(\alpha_{2}+\lambda_{2}\right)\right]\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}>2\,:&0\end{array}\right.$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{2}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{2}=1\,:&e^{\,-\frac{\lambda_{1}+\lambda_{2}}{2}}\cdot\\ &\left[\left(\frac{\lambda_{3}}{2}\right)^{2}+\left(\alpha_{3}+1\right)\left(\alpha_{3}+\lambda_{3}\right)\right]\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}>2\,:&0\end{array}\right.$

Hence, thanks to its richer parametrization based on three non-centrality parameters, the density of the bivariate Non-central Dirichlet distribution enables to overcome the aforementioned limitations of the Dirichlet by allowing its own limits at the vertices of the unit simplex to have arbitrary positive and finite values under unitary shape parameters. This feature, which results from the flexibility of the limit at 0 of the Non-central Chi-Squared density previously noticed in Eq. (46), makes the NcDir model potentially apt to properly capture the data portions having values next to the vertices of the support.

3 A new Non-central Dirichlet model

3.1 Definition and distribution

The above recalled characterizing property of independent Gamma random variables is no longer valid in the non-central setting. This fact covers a key role in the present paper and, more specifically, is assumed as the starting point in the derivation of the new family of non-central generalizations of the Dirichlet distribution we are interested in. Indeed, in the notation of Eq. (52), a new Non-central Dirichlet model can be achieved by conditioning the existing one $\underline{X}^{\prime}$ on the sum $Y^{\prime+}$ of the $D+1$ independent Non-central Chi-Squared random variables involved in its definition. Such a model is thus distributed according to the conditional distribution of $\underline{X}^{\prime}$ given $Y^{\prime+}$ ; moreover, as $\underline{X}^{\prime}$ is not independent of $Y^{\prime+}$ , this distribution must be different from the $\mbox{NcDir}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ except for $\underline{\lambda}=\underline{0}$ , this latter case corresponding to the Dirichlet. In this regard, the conditional density of $\underline{X}^{\prime}$ given $Y^{\prime+}=y$ can be obtained for any fixed $y>0$ by mixing the conditional density of $\underline{X}^{\prime}$ given $(\underline{M},Y^{\prime+}=y)$ with respect to the conditional probability mass function of $\underline{M}$ given $Y^{\prime+}=y$ . As $\underline{X}^{\prime}$ is conditionally independent of $Y^{\prime+}$ given $\underline{M}$ by Property 2.1, in light of Eq. (54) we have that $\underline{X}^{\prime}\left|\right.(\underline{M},Y^{\prime+}=y)\sim\mbox{Dir}^{\,D}\left(\underline{\alpha}+\underline{M}\right)$ , while a direct application of Bayes’ Theorem shows that:

$\displaystyle\Pr\left(\left.\underline{M}=\underline{j}\;\right\|\;Y^{\prime+}=y\right)=$			(66)
	$\displaystyle=$	$\displaystyle\Pr\left(\left.M_{1}=j_{1}\,,\,\ldots\,,\,M_{D+1}=j_{D+1}\;\right\|\;Y^{\prime+}=y\right)\qquad\left\{\begin{array}[]{l}j_{i}\in\mathbb{N}_{0}\,,\quad i=1,\,\ldots\,,D+1\\ j^{+}=\sum_{i=1}^{D+1}j_{i}\end{array}\right.$	(66)
	$\displaystyle=$	$\displaystyle\frac{1}{{}_{0}F_{1}\left(\alpha^{+};\frac{y\,\lambda^{+}}{4}\right)}\,\frac{1}{\left(\alpha^{+}\right)_{j^{+}}}\prod_{i=1}^{D+1}\frac{\left(\frac{y\,\lambda_{i}}{4}\right)^{j_{i}}}{j_{i}\,!}\;.$	(67)

Therefore, by incorporating the conditioning value $y$ in the $\lambda_{i}$ ’s, without loss of generality we are led to the following definitions.

Definition 3.1 ( $D$ -variate Conditional Non-central Dirichlet distribution)

The $D$ -variate Conditional Non-central Dirichlet distribution with vector of shape parameters $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D+1}\right)$ and vector of non-centrality parameters $\underline{\lambda}=\left(\lambda_{1},\,\ldots\,,\lambda_{D+1}\right)$ , denoted by $\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ , is the distribution of the $D$ -dimensional random vector $\underline{X}^{\prime\prime}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{D}\right)$ defined as follows:

	$\displaystyle\left\{\begin{array}[]{l}Y^{\prime}_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\prime\,2}_{2\alpha_{i}}\left(\lambda_{i}\right)\quad i=1,\,\ldots\,,D+1\\ \\ Y^{\prime+}=\sum_{i=1}^{D+1}Y^{\prime}_{i}\end{array}\right.\quad\Rightarrow$				(71)
		$\displaystyle\Rightarrow$	$\displaystyle\quad\underline{X}^{\prime\prime}=\left.\left(\frac{Y^{\prime}_{1}}{Y^{\prime+}},\,\ldots\,,\frac{Y^{\prime}_{D}}{Y^{\prime+}}\right)\,\right\|Y^{\prime+}=1\quad\sim\;\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)\,.$		(72)

Definition 3.2 ( $(D+1)$ -variate Mixture Weight distribution)

The $(D+1)$ -variate Mixture Weight distribution with shape parameter $\alpha^{+}$ and vector of non-centrality parameters $\underline{\lambda}=\left(\lambda_{1},\,\ldots\,,\lambda_{D+1}\right)$ , denoted by $\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ , is the distribution of the $(D+1)$ -dimensional random vector $\underline{N}=\left(N_{1},\,\ldots\,,N_{D+1}\right)$ obtained as follows:

	$\displaystyle\left\{\begin{array}[]{l}\underline{M}\sim\mbox{{Multi-Poisson}}^{\,D+1}(\underline{\lambda}\,/\,2)\\ \\ Y^{\prime}_{i}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\prime\,2}_{2\alpha_{i}}\left(\lambda_{i}\right)\,,\;\;i=1,\,\ldots\,,D+1\\ \\ Y^{\prime+}=\sum_{i=1}^{D+1}Y^{\prime}_{i}\end{array}\right.\quad\Rightarrow$				(78)
		$\displaystyle\qquad\qquad\qquad\Rightarrow\qquad$	$\displaystyle\underline{N}=\left.\underline{M}\;\right\|\,Y^{\prime+}=1\quad\sim\;\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)\;.$		(79)

The special case of Eq. (79) where $\underline{\lambda}$ is set equal to $\underline{0}$ leads to a random vector degenerate at $\underline{0}$ . By Eq. (67) and in view of the foregoing arguments, the joint probability mass function of the $\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ distribution takes on the following form:

	$\displaystyle\Pr\left(\underline{N}=\underline{j}\right)\qquad=$				(83)
		$\displaystyle=$	$\displaystyle\frac{1}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,\frac{1}{\left(\alpha^{+}\right)_{j^{+}}}\prod_{i=1}^{D+1}\frac{\left(\frac{\lambda_{i}}{4}\right)^{j_{i}}}{j_{i}\,!}\,,\qquad\left\{\begin{array}[]{l}\underline{j}=\left(j_{1},\,\ldots\,,j_{D+1}\right)\in\mathbb{N}_{0}^{\,D+1}\\ \\ j^{+}=\sum_{i=1}^{D+1}j_{i}\end{array}\right..$		(83)

Therefore, the density function of $\underline{X}^{\prime\prime}\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ can be accordingly expressed as follows:

\mbox{CNcDir}^{\,D}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)=\sum_{\underline{j}\in\mathbb{N}_{0}^{\,D+1}}\left[\Pr\left(\underline{N}=\underline{j}\right)\cdot\mbox{Dir}^{\,D}\left(\underline{x};\underline{\alpha}+\underline{j}\right)\right]\,,\qquad\underline{x}\in\mathcal{S}^{\,D}\,,

(84)

i.e. as the multiple infinite series of the $\mbox{Dir}^{\,D}\left(\underline{\alpha}+\underline{j}\right)$ densities, $\underline{j}\in\mathbb{N}_{0}^{\,D+1}$ , weighted by the joint probabilities of the random vector $\underline{N}$ in Eq. (83). Hence, the CNcDir density shows the same mixture type form as the NcDir one; the only difference lies in the mixing distribution, given by the MW for the former and the Multi-Poisson for the latter. In this regard, an overview of the properties of the multivariate Mixture Weight distribution is deferred to the next subsection. More interestingly, by simple computations it is easy to see that the CNcDir density turns out to be equivalently stated in terms of the following perturbation of the corresponding central case:

	$\displaystyle\mbox{CNcDir}^{\,D}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)\quad=\quad\mbox{Dir}^{\,D}\left(\underline{x};\underline{\alpha}\right)\cdot$				(85)
		$\displaystyle\cdot$	$\displaystyle\frac{\left[\prod_{i=1}^{D}\,{}_{0}F_{1}\left(\alpha_{i};\frac{\lambda_{i}}{4}x_{i}\right)\right]\,_{0}F_{1}\left[\alpha_{D+1};\frac{\lambda_{D+1}}{4}\left(1-\sum_{i=1}^{D}x_{i}\right)\right]}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,,\qquad\underline{x}\in\mathcal{S}^{D}\,.$		(85)

From Eq. (85) the major tractability and interpretability of the CNcDir density over the NcDir one prove evident. Contrarily to the standard non-central case where the perturbing factor of the Dirichlet density shows an inner structure that is too complex to be handled analytically, with reference to the density of the conditional model the perturbing effect can be clearly noticed for each value of the parameter vector. More precisely, in Eq. (85), regardless of the constant term, the Dirichlet density is perturbed by the product of $D+1$ generalized hypergeometric functions ${}_{0}F_{1}$ which show perfectly symmetric behaviors. Specifically, as $\lambda_{i}$ gets higher, $i=1,\,\ldots\,,D+1$ , the corresponding function ${}_{0}F_{1}$ gives more weight to the tail of the Dirichlet density relative to the $i$ -th vertex of the unit simplex; moreover, the less is the corresponding shape parameter $\alpha_{i}$ , the larger is the extent of this phenomenon.

The Conditional Non-central Dirichlet model shares with the Dirichlet and the existing Non-central Dirichlet some important properties that are discussed below.

The first one is the closure under permutation, an apparently trivial property which expresses the fact that the CNcDir distribution treats all its components in a completely symmetric way.

Proposition 3.1 (Closure under permutation)

Let $\underline{X}^{\prime\prime}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{D}\right)\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ where $\underline{\alpha}=(\alpha_{1},\,\ldots\,,\alpha_{D},\alpha_{D+1})$ , $\underline{\lambda}=(\lambda_{1},\ldots,$ $\lambda_{D},\lambda_{D+1})$ and $(\underline{k})=(k_{(1)},\,\ldots\,,k_{(D)},k_{(D+1)})$ be a permutation of $(1,\,\ldots\,,D,D+1)$ . Then, $\underline{X}^{\prime\prime}_{(\underline{k})}=(X^{\prime\prime}_{k_{(1)}},\,\ldots\,,X^{\prime\prime}_{k_{(D)}})\sim\mbox{{CNcDir}}^{\,D}(\underline{\alpha}_{(\underline{k})},\underline{\lambda}_{(\underline{k})})$ where $\underline{\alpha}_{(\underline{k})}=(\alpha_{k_{(1)}},\,\ldots\,,\alpha_{k_{(D)}},$ $\alpha_{k_{(D+1)}})$ and $\underline{\lambda}_{(\underline{k})}=(\lambda_{k_{(1)}},\,\ldots\,,\lambda_{k_{(D)}},\lambda_{k_{(D+1)}})$ .

Proof.

The proof directly follows from the definition of the model at study in Eq. (72) by applying the permutation $(\underline{k})$ to the $Y^{\prime}_{i}$ ’s. ∎

Then, like the Dirichlet and the Non-central Dirichlet, the present model satisfies the aggregation property, which is crucial to the derivation of the results established in the following.

Proposition 3.2 (Aggregation property)

Let $\underline{X}^{\prime\prime}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{D}\right)\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ where $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D},\alpha_{D+1}\right)$ , $\underline{\lambda}=\left(\lambda_{1},\ldots,\right.$ $\left.\lambda_{D},\lambda_{D+1}\right)$ . Let $\mathcal{P}=\left\{A_{1},\,\ldots\,,A_{m},A_{m+1}\right\}$ be a partition of $\left\{1,\,\ldots\,,D,D+1\right\}$ , $1\leq m\leq D$ . Then:

\left(\sum_{i\in A_{1}}X^{\prime\prime}_{i},\,\ldots\,,\sum_{i\in A_{m}}X^{\prime\prime}_{i}\right)\sim\mbox{{CNcDir}}^{\,m}\left(\sum_{i\in A_{1}}\alpha_{i},\,\ldots\,,\sum_{i\in A_{m+1}}\alpha_{i},\sum_{i\in A_{1}}\lambda_{i},\,\ldots\,,\sum_{i\in A_{m+1}}\lambda_{i}\right)

(86)

Proof.

The proof directly ensues from the definition of the model at study in Eq. (72) and from the reproductive property of the Non-central Chi-Squared distribution in Eq. (47) by applying the partition $\mathcal{P}$ to the $Y^{\prime}_{i}$ ’s. ∎

By suitably resorting to Proposition 3.2, the $k$ -dimensional marginals of the $D$ -dimensional Conditional Non-central Dirichlet model can be simply achieved for every $k=1,\ldots,D-1$ . Specifically, like the Dirichlet and the existing Non-central Dirichlet models, the Conditional Non-central Dirichlet is closed under marginalization and simple relationships hold between the parameters of the joint and the marginal distributions. In particular, the univariate marginals are of Conditional Doubly Non-central Beta (CDNcB) type, this latter being the analogue of the CNcDir on the Real interval $(0,1)$ . In this regard, a first analysis and an in-depth study of the CDNcB distribution are given in [Ongaro and Orsi (2015)] and [Orsi (2021)], respectively. In the following, only the univariate and bivariate marginals are made explicit but as far as the investigation of the $k$ -dimensional marginals, $k>2$ , is concerned, similar results apply.

Proposition 3.3 (Univariate and bivariate marginals)

Let $\underline{X}^{\prime\prime}\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ . Then, for every $p,\,q=1,\,\ldots\,,D$ , $p\neq q$ :

X^{\prime\prime}_{p}\sim\mbox{{CDNcB}}\left(\alpha_{p},\alpha^{+}-\alpha_{p},\lambda_{p},\lambda^{+}-\lambda_{p}\right)\,,

(87)

\left(X^{\prime\prime}_{p},X^{\prime\prime}_{q}\right)\sim\mbox{{CNcDir}}^{\,2}\left(\alpha_{p},\alpha_{q},\alpha^{+}-(\alpha_{p}+\alpha_{q}),\lambda_{p},\lambda_{q},\lambda^{+}-(\lambda_{p}+\lambda_{q})\right)\,.

(88)

Proof.

The proof of Eq. (87) ensues from Eq. (86) by setting $m=1$ and taking $A_{1}=\left\{p\right\}$ , $A_{2}=\left\{1,\,\ldots\,,p-1,p+1,\,\ldots\,,D,D+1\right\}$ for every $p=1,\,\ldots\,,D$ , whereas the proof of Eq. (88) follows from Eq. (86) by setting $m=2$ and taking $A_{1}=\left\{p\right\}$ , $A_{2}=\left\{q\right\}$ , $A_{3}=\left\{1,\,\ldots\,,D,D+1\right\}-\left\{p,q\right\}$ , for every $p,q=1,\,\ldots\,,D$ , $p\neq q$ . ∎

Another consequence of the aggregation property is about the distribution of the sum of the components of a CNcDir random vector.

Proposition 3.4 (Sum of the components)

Let $\underline{X}^{\prime\prime}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{D}\right)\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ where $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D},\alpha_{D+1}\right)$ , $\underline{\lambda}=\left(\lambda_{1},\ldots,\right.$ $\left.\lambda_{D},\lambda_{D+1}\right)$ . Then, $X^{\prime\prime+}=\sum_{j=1}^{D}X^{\prime\prime}_{j}\sim\mbox{{CDNcB}}\left(\sum_{j=1}^{D}\alpha_{j},\alpha_{D+1},\sum_{j=1}^{D}\lambda_{j},\lambda_{D+1}\right)$ .

Proof.

The proof ensues from Proposition 3.2 by setting $m=1$ and taking $A_{1}=\left\{1,\,\ldots\,,D\right\}$ , $A_{2}=\left\{D+1\right\}$ . ∎

We end this subsection pointing out that, unlike the standard Non-central Dirichlet case, the Conditional Non-central Dirichlet distribution is closed under conditioning after normalization.

Proposition 3.5 (Closure under normalized conditioning)

Let $\underline{X}^{\prime\prime}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{D}\right)\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ where $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D+1}\right)$ , $\underline{\lambda}=\left(\lambda_{1},\,\ldots\,,\right.$ $\left.\lambda_{D+1}\right)$ . Let $\underline{X}^{\prime\prime}_{1}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{k}\right)$ and $\underline{X}^{\prime\prime}_{2}=\left(X^{\prime\prime}_{k+1},\,\ldots\,,X^{\prime\prime}_{D}\right)$ , $1\leq k\leq D-1$ . Then, the normalized conditional distribution of $\underline{X}^{\prime\prime}_{2}$ given $\underline{X}^{\prime\prime}_{1}=\underline{x}_{1}$ is:

	$\displaystyle\left.\frac{\underline{X}^{\prime\prime}_{2}}{1-x_{1}^{+}}\right\|\,\underline{X}^{\prime\prime}_{1}=\underline{x}_{1}=\left(x_{1},\,\ldots\,,x_{k}\right),\;x_{1}^{+}=\sum_{j=1}^{k}x_{j}\quad\sim$				(92)
		$\displaystyle\sim\quad$	$\displaystyle\mbox{{CNcDir}}^{\,D-k}\left(\underline{\alpha}_{2},\underline{\lambda}_{2}\left(1-x_{1}^{+}\right)\right)\,,\qquad\left\{\begin{array}[]{l}\underline{\alpha}_{2}=\left(\alpha_{k+1},\,\ldots\,,\alpha_{D+1}\right)\\ \\ \underline{\lambda}_{2}=\left(\lambda_{k+1},\,\ldots\,,\lambda_{D+1}\right)\end{array}\right.\,.$		(92)

Proof.

By setting $m=k$ and taking $A_{j}=\left\{j\right\}$ for every $j=1,\,\ldots\,,k$ , $A_{k+1}=\left\{k+1,\,\ldots\,,\right.$ $\left.D+1\right\}$ in Eq. (86), one has $\underline{X}^{\prime\prime}_{1}\sim\mbox{{CNcDir}}^{\,k}\left(\underline{\alpha}_{1},\alpha_{2}^{+},\underline{\lambda}_{1},\lambda_{2}^{+}\right)$ , where $\underline{\alpha}_{1}=\left(\alpha_{1},\,\ldots\,,\alpha_{k}\right)$ , $\underline{\lambda}_{1}=\left(\lambda_{1},\,\ldots\,,\lambda_{k}\right)$ , $\alpha_{2}^{+}=\sum_{i=k+1}^{D+1}\alpha_{i}$ and $\lambda_{2}^{+}=\sum_{i=k+1}^{D+1}\lambda_{i}$ . For every $\underline{x}_{2}=\left(x_{k+1},\,\ldots\,,\right.$ $\left.x_{D}\right)$ , the conditional density of $\underline{X}^{\prime\prime}_{2}$ given $\underline{X}^{\prime\prime}_{1}=\underline{x}_{1}$ admits the following expression:

$\displaystyle f_{\underline{X}^{\prime\prime}_{2}\,\|\,\underline{X}^{\prime\prime}_{1}=\,\underline{x}_{1}}\left(\underline{x}_{2};\underline{\alpha}_{2},\underline{\lambda}_{2}\right)\quad=$
	$\displaystyle=$	$\displaystyle\frac{\Gamma\left(\alpha_{2}^{+}\right)}{\prod_{i=k+1}^{D+1}\Gamma\left(\alpha_{i}\right)}\left[\prod_{i=k+1}^{D}\left(\frac{x_{i}}{1-x_{1}^{+}}\right)^{\alpha_{i}-1}\right]\left(1-\frac{x_{2}^{+}}{1-x_{1}^{+}}\right)^{\alpha_{D+1}-1}$
	$\displaystyle\cdot$	$\displaystyle\frac{\left[\prod_{i=k+1}^{D}\,{}_{0}F_{1}\left(\alpha_{i};\frac{\lambda_{i}}{4}x_{i}\right)\right]\cdot\,_{0}F_{1}\left[\alpha_{D+1};\frac{\lambda_{D+1}}{4}\left(1-x^{+}\right)\right]}{\left(1-x_{1}^{+}\right)^{D-k}\,_{0}F_{1}\left[\alpha_{2}^{+};\frac{\lambda_{2}^{+}}{4}\left(1-x_{1}^{+}\right)\right]}\;.$

Now consider the $(D-k)$ -dimensional random vector $\underline{Y}$ defined by the following linear transformation:

\underline{Y}=h\left(\underline{X}^{\prime\prime}_{2}\right)=\frac{\underline{X}^{\prime\prime}_{2}}{1-x_{1}^{+}}\quad\Leftrightarrow\quad\left\{\begin{array}[]{l}\underline{X}^{\prime\prime}_{2}=h^{-1}\left(\underline{Y}\right)=\left(1-x_{1}^{+}\right)\,\underline{Y}\\ \\ |J|=\mbox{detJac}\,h^{-1}\left(\underline{Y}\right)=(1-x_{1}^{+})^{D-k}\end{array}\right.\,;

note that the conditional density of this latter given $\underline{X}^{\prime\prime}_{1}=\underline{x}_{1}$ takes on the form of

$\displaystyle f_{\underline{Y}\,\|\,\underline{X}^{\prime\prime}_{1}=\,\underline{x}_{1}}\left(\underline{y};\underline{\alpha}_{2},\underline{\lambda}_{2}\left(1-x_{1}^{+}\right)\right)=$			(93)
	$\displaystyle=$	$\displaystyle f_{\underline{X}^{\prime\prime}_{2}\,\|\,\underline{X}^{\prime\prime}_{1}=\,\underline{x}_{1}}\left(\left(1-x_{1}^{+}\right)\underline{y};\underline{\alpha}_{2},\underline{\lambda}_{2}\right)\cdot\|J\|=\quad\mbox{{Dir}}^{\,D-k}\left(\underline{y};\underline{\alpha}_{2}\right)\cdot$
	$\displaystyle\cdot$	$\displaystyle\frac{\prod_{i=k+1}^{D}\,{}_{0}F_{1}\left[\alpha_{i};\frac{\lambda_{i}}{4}\left(1-x_{1}^{+}\right)y_{i}\right]\cdot\,_{0}F_{1}\left[\alpha_{D+1};\frac{\lambda_{D+1}}{4}\left(1-x_{1}^{+}\right)\left(1-y^{+}\right)\right]}{{}_{0}F_{1}\left[\alpha_{2}^{+};\frac{\lambda_{2}^{+}}{4}\left(1-x_{1}^{+}\right)\right]}\,,\quad$

where $y^{+}=\sum_{i=k+1}^{D}y_{i}$ . Clearly, Eq. (93) corresponds to the density of the distribution in Eq. (92). ∎

3.2 Mixing distribution

The $(D+1)$ -variate Mixture Weight distribution specified in Eq. (83) plays the role of mixing distribution in the mixture type form of the CNcDir ${}^{\,D}$ density. Such model represents the generalization to $D+1$ dimensions $(D>1)$ of the $\mbox{{MW}}^{\,2}$ , this latter being the mixing distribution in the mixture type form of the CDNcB density. Hence, the properties of the $\mbox{{MW}}^{\,D+1}$ distribution illustrated in the present subsection are extensions of the ones of the bivariate case, an in-depth study of which is provided in [Orsi (2021)].

First of all, the $\mbox{{MW}}^{\,D+1}$ distribution is closed under permutation and the non-centrality parameters of the permuted random vector are obtained by applying the same permutation to the original ones.

Proposition 3.6 (Closure under permutation)

Let $\underline{N}=\left(N_{1},\,\ldots\,,N_{D+1}\right)\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ where $\underline{\lambda}=\left(\lambda_{1},\,\ldots\,,\lambda_{D+1}\right)$ and $(\underline{k})=\left(k_{(1)},\,\ldots\,,k_{(D+1)}\right)$ be a permutation of $\left\{1,\,\ldots\,,D+1\right\}$ . Then, $\underline{N}_{(\underline{k})}=\left(N_{k_{(1)}},\,\ldots\,,\right.$ $\left.N_{k_{(D+1)}}\right)\sim\mbox{{MW}}^{\,D+1}(\alpha^{+},\underline{\lambda}_{(\underline{k})})$ where $\underline{\lambda}_{(\underline{k})}=(\lambda_{k_{(1)}},\,\ldots\,,\lambda_{k_{(D+1)}})$ .

Proof.

The proof directly follows from the definition of the model at study in Eq. (79) by applying the permutation $(\underline{k})$ to the $Y^{\prime}_{i}$ ’s. ∎

Proposition 3.6 fully expresses the symmetric nature of the $\mbox{{MW}}^{\,D+1}$ distribution, implying that the particular choice of its components does not affect any of the following results.

The distribution under consideration is not closed under marginalization, but it is closed under conditioning.

Proposition 3.7 (Marginals and conditionals)

Let $\underline{N}=\left(N_{1},\,\ldots\,,N_{D+1}\right)\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ . Then, for every $m=1,\,\ldots\,,D$ , the marginal probability mass function of $\left(N_{1},\,\ldots\,,N_{m}\right)$ is

	$\displaystyle\Pr\left(\,\left(N_{1},\,\ldots\,,N_{m}\right)=\left(j_{1},\,\ldots\,,j_{m}\right)\,\right)\quad=$
		$\displaystyle=$	$\displaystyle\frac{1}{\left(\alpha^{+}\right)_{j^{+}}}\left[\prod_{i=1}^{m}\frac{\left(\frac{\lambda_{j_{i}}}{4}\right)^{j_{i}}}{j_{i}\,!}\right]\,\frac{{}_{0}F_{1}\left(\alpha^{+}+j^{+};\frac{\lambda^{+}-\sum_{i=1}^{m}\lambda_{j_{i}}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,,\quad\left\{\begin{array}[]{l}\left(j_{1},\,\ldots\,,j_{m}\right)\in\mathbb{N}_{0}^{\,m}\\ \\ j^{+}=\sum_{i=1}^{m}j_{i}\end{array}\right.$

and the conditional distribution of $\left(N_{1},\,\ldots\,,N_{m}\right)$ given $\left(N_{m+1},\,\ldots\,,N_{D+1}\right)$ is

\left.\left(N_{1},\,\ldots\,,N_{m}\right)\,\right|\,\left(N_{m+1},\,\ldots\,,N_{D+1}\right)\;\sim\;\mbox{{MW}}^{\,m}\left(\alpha^{+}+\sum_{i=m+1}^{D+1}N_{i}\,,\lambda_{1},\,\ldots\,,\lambda_{m}\right)\,.

Proof.

The proof is straightforward and is omitted. ∎

The sum of the components of $\underline{N}\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ belongs to the univariate Mixture Weight family.

Proposition 3.8 (Sum of the components)

Let $\underline{N}=\left(N_{1},\,\ldots\,,N_{D+1}\right)\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ . Then:

N^{+}=\sum_{i=1}^{D+1}N_{i}\;\sim\;\mbox{{MW}}^{\,1}\left(\alpha^{+},\lambda^{+}\right)\,.

(95)

Proof.

By employing the probability mass function of the sum of the components of the $\mbox{{Multi-Poisson}}^{\,D+1}$ distribution in Eq. (53), for every $s\in\mathbb{N}_{0}$ the probabilities of $N^{+}$ take on the following form:

$\displaystyle\Pr\left(N^{+}=s\right)=$
	$\displaystyle=$	$\displaystyle\sum_{\sum_{i=1}^{D}s_{i}\leq s}\Pr\left(N_{1}=s_{1},\,\ldots\,,N_{D}=s_{D},N_{D+1}=s-\sum_{i=1}^{D}s_{i}\right)=$
	$\displaystyle=$	$\displaystyle\frac{1}{\left(\alpha^{+}\right)_{s}\,{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\sum_{\sum_{i=1}^{D}s_{i}\leq s}\left\{\left[\prod_{i=1}^{D}\frac{\left(\frac{\lambda_{i}}{4}\right)^{s_{i}}}{s_{i}\,!}\right]\,\frac{\left(\frac{\lambda_{D+1}}{4}\right)^{s-\sum_{i=1}^{D}s_{i}}}{\left(s-\sum_{i=1}^{D}s_{i}\right)!}\right\}=$
	$\displaystyle=$	$\displaystyle\frac{e^{\frac{\lambda^{+}}{4}}}{\left(\alpha^{+}\right)_{s}\,{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\sum_{\sum_{i=1}^{D}s_{i}\leq s}\left\{\left[\prod_{i=1}^{D}e^{-\frac{\lambda_{i}}{4}}\,\frac{\left(\frac{\lambda_{i}}{4}\right)^{s_{i}}}{s_{i}\,!}\right]\frac{e^{-\frac{\lambda_{D+1}}{4}}\left(\frac{\lambda_{D+1}}{4}\right)^{s-\sum_{i=1}^{D}s_{i}}}{\left(s-\sum_{i=1}^{D}s_{i}\right)!}\right\}$
	$\displaystyle=$	$\displaystyle\frac{1}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\frac{1}{\left(\alpha^{+}\right)_{s}}\frac{\left(\frac{\lambda^{+}}{4}\right)^{s}}{s\,!}\qquad s\in\mathbb{N}_{0}\,,$

which corresponds to the special case of Eq. (83) where $D+1$ is set equal to 1. ∎

Despite the differences existing between their joint probability mass functions, the $\mbox{{MW}}^{\,D+1}$ and the $\mbox{{Multi-Poisson}}^{\,D+1}$ share the same conditional distribution given the sum of their components.

Proposition 3.9 (Conditional distribution given the sum of the components)

\left.\left(N_{1},\,\ldots\,,N_{D}\right)\,\right|\,N^{+}\sim\mbox{{Multinomial}}^{\,D}\left(N^{+},\frac{\lambda_{1}}{\lambda^{+}},\,\ldots\,,\frac{\lambda_{D}}{\lambda^{+}}\right)\,.

(96)

Proof.

Let $\left(j_{1},\,\ldots\,,j_{D}\right)\in\mathbb{N}_{0}^{\,D}$ , $s\in\mathbb{N}_{0}$ and $\sum_{i=1}^{D}j_{i}=j^{+}\leq s$ ; then, the proof ensues from Eqs. (83) and (95) by carrying out the following computations:

$\displaystyle\Pr\left(\,\left.\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}=s\right)=$
	$\displaystyle=$	$\displaystyle\frac{\Pr\left(N_{1}=j_{1},\,\ldots\,,N_{D}=j_{D},N_{D+1}=s-j^{+}\right)}{\Pr\left(N^{+}=s\right)}=$
	$\displaystyle=$	$\displaystyle\frac{1}{\left(\alpha^{+}\right)_{s}\,{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\left[\prod_{i=1}^{D}\frac{\left(\frac{\lambda_{i}}{4}\right)^{j_{i}}}{j_{i}\,!}\right]\,\frac{\left(\frac{\lambda_{D+1}}{4}\right)^{s-j^{+}}}{\left(s-j^{+}\right)!}\left[\frac{1}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\frac{1}{\left(\alpha^{+}\right)_{s}}\frac{\left(\frac{\lambda^{+}}{4}\right)^{s}}{s\,!}\right]^{-1}$
	$\displaystyle=$	$\displaystyle{s\choose j_{1},\,\ldots\,,j_{D}}\left[\prod_{i=1}^{D}\left(\frac{\lambda_{i}}{\lambda^{+}}\right)^{j_{i}}\right]\left(1-\sum_{i=1}^{D}\frac{\lambda_{i}}{\lambda^{+}}\right)^{s-j^{+}}\,,$

this latter being the probability mass function of the distribution in Eq. (96). ∎

Finally, the $\mbox{{MW}}^{\,D+1}$ distribution can be simulated by first generating the random variable $N^{+}$ according to Eq. (95) and then the $D$ -variate Multinomial distribution in Eq. (96). In particular, the former can be easily simulated by using, for example, the Inverse-Transform method [Law and Kelton (2000)]. Specifically, after generating a realization $u$ from a Uniform random variable on the Real interval $(0,1)$ , an exact realization $s\in\mathbb{N}_{0}$ from the random variable $N^{+}$ is obtainable by taking the generalized inverse of its distribution function $F$ in $u$ , i.e. $F^{-}\left(u\right)=\inf\{x\,|\,F\left(x\right)\geq u\}$ .

3.3 Representations

In the present subsection a number of fundamental representations of the Conditional Non-central Dirichlet model are provided, disclosing its remarkable tractability.

The first one is the mixture representation, which summarizes the arguments leading to Definitions 3.1 and 3.2.

Proposition 3.10 (Mixture representation)

Let $\underline{X}^{\prime\prime}\sim$ $\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ and $\underline{N}\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ . Then:

\left.\underline{X}^{\prime\prime}\,\right|\,\underline{N}\sim\mbox{{Dir}}^{\,D}\left(\underline{\alpha}+\underline{N}\right)\;.

(97)

Proof.

The proof is straightforward and is omitted. ∎

Proposition 3.10 proves significant in several respects. Indeed, the above representation constitutes a valid practical method for generating from a Conditional Non-central Dirichlet random vector. Specifically, the algorithm to obtain a realization from the $\mbox{{CNcDir}}^{\,D}$ model using such a representation requires to first generate from the $\mbox{MW}^{\,D+1}$ distribution and then from the conditional distribution in Eq. (97). In particular, the former operation can be accomplished by following the lines made explicit in the previous subsection. In this regard, an interesting restatement of the combination of the two aforementioned generation steps is provided by the following result, which shows that, conditionally on $N^{+}$ , the density of the model of interest is given by a mixture of Dirichlet densities weighted by the Multinomial probabilities in Eq. (96).

Proposition 3.11 (Conditional density given $N^{+}$ )

Let $\underline{X}^{\prime\prime}\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ , $\underline{N}\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ and $N^{+}=\sum_{i=1}^{D+1}N_{i}$ . Then, the conditional density of $\underline{X}^{\prime\prime}$ given $N^{+}$ is

	$\displaystyle f_{\left.\underline{X}^{\prime\prime}\,\right\|\,N^{+}}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)=\sum_{j^{+}\,\leq\,N^{+}}\left[\mbox{{Multinomial}}^{\,D}\left(j_{1},\,\ldots\,,j_{D};\,N^{+},\frac{\lambda_{1}}{\lambda^{+}},\,\ldots\,,\frac{\lambda_{D}}{\lambda^{+}}\right)\cdot\right.$				(101)
		$\displaystyle\cdot$	$\displaystyle\left.\mbox{{Dir}}^{\,D}\left(\underline{x};\alpha_{1}+j_{1},\,\ldots\,,\alpha_{D}+j_{D},\alpha_{D+1}+N^{+}-j^{+}\right)\right]\qquad\left\{\begin{array}[]{l}\left(j_{1},\,\ldots\,,j_{D}\right)\in\mathbb{N}_{0}^{\,D}\\ \\ j^{+}=\sum_{i=1}^{D}j_{i}\end{array}\right.$		(101)

Proof.

The proof ensues from the statement

	$\displaystyle f_{\left.\underline{X}^{\prime\prime}\,\right\|\,N^{+}}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)$	$\displaystyle=$	$\displaystyle\sum_{j^{+}\,\leq\,N^{+}}\left[\,\Pr\left(\left.\,\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}\right)\right.\cdot$
		$\displaystyle\cdot$	$\displaystyle\left.f_{\left.\underline{X}^{\prime\prime}\,\right\|\,\left(\,N_{1},\,\ldots\,,N_{D},N^{+}\right)}\left(\underline{x};\underline{\alpha},\underline{\lambda},\left(j_{1},\,\ldots\,,j_{D}\right)\right)\,\right]\,,$

where $\left(j_{1},\,\ldots\,,j_{D}\right)\in\mathbb{N}_{0}^{\,D}$ , $j^{+}=\sum_{i=1}^{D}j_{i}$ , by using Eq. (96) and noting that Eq. (97) can be rearranged as follows:

\left.\underline{X}^{\prime\prime}\,\right|\,\left(\,N_{1},\,\ldots\,,N_{D},N^{+}\right)\sim\mbox{{Dir}}^{\,D}\left(\alpha_{1}+N_{1},\,\ldots\,,\alpha_{D}+N_{D},\alpha_{D+1}+N^{+}-\sum_{i=1}^{D}N_{i}\right)\;.

∎

Finally, the mixture representation in Eq. (97) leads to the following representation of the CNcDir model in terms of unconditional composition.

Proposition 3.12 (Unconditional composition type representation)

Let $\underline{N}\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ and $\left.Z^{\prime}_{i}\,\right|\,\underline{N}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\,2}_{2\alpha_{i}+2N_{i}}$ , $i=1,\ldots,D+1$ with $Z^{\prime+}=\sum_{i=1}^{D+1}Z^{\prime}_{i}$ . Then:

i)

$\left(\frac{Z^{\prime}_{1}}{Z^{\prime+}},\,\ldots\,,\frac{Z^{\prime}_{D}}{Z^{\prime+}}\right)\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ ,
ii)

$Z^{\prime}_{i}=Z_{i}+\sum_{j=1}^{N_{i}}F_{j}$ , $i=1,\,\ldots\,,D+1$ , where $Z_{i}\sim\chi^{\,2}_{2\alpha_{i}}$ , $\underline{N}$ , $\{F_{j}\stackrel{{\scriptstyle\bot}}{{\sim}}\chi^{\,2}_{2}\}$ are mutually independent.

Proof.

The proof is straightforward and is omitted. ∎

The issue of generating from the Conditional Non-central Dirichlet model can be also addressed by resorting to i) in Proposition 3.12. Specifically, the simulation methods based on the mixture representation and on the above unconditional composition type representation differ only in the way that the Dirichlet distribution is generated: directly by the former and via Chi-Squareds by the latter. Moreover, from ii) in Proposition 3.12 it is noticeable that the $Z^{\prime}_{i}$ ’s play the same role as the Non-central Chi-Squared random variables involved in the definition of the standard Non-central Dirichlet model. Hence, the present representation guides us to view the CNcDir model as an analogue of the NcDir one. Indeed, like the $Y^{\prime}_{i}$ ’s in Eq. (44), the $Z^{\prime}_{i}$ ’s can be additively decomposed into two independent parts, a central one and a purely non-central one. Moreover, by resorting to the $Z^{\prime}_{i}$ ’s, the following generalization of the independence relationship stated in Property 2.1 holds true.

Proposition 3.13 (Conditional independence)

Let $\underline{X}^{\prime\prime}\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ and $\underline{N}=\left(N_{1},\,\ldots\,,N_{D+1}\right)\sim\mbox{{MW}}^{\,D+1}\left(\alpha^{+},\underline{\lambda}\right)$ with $N^{+}=\sum_{i=1}^{D+1}N_{i}$ . Furthermore, let $\left(Z^{\prime}_{1},\,\ldots\,,Z^{\prime}_{D+1}\right)$ be defined as in Proposition 3.12 with $Z^{\prime+}=\sum_{i=1}^{D+1}Z^{\prime}_{i}$ . Then, $\underline{X}^{\prime\prime}$ and $Z^{\prime+}$ are conditionally independent given $N^{+}$ .

Proof.

The ingredients of the representation of the CNcDir ${}^{\,D}$ model in Proposition 3.12 can be analogously defined as follows:

Z^{\prime}_{i}\left.\,\right|\,\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)\;\stackrel{{\scriptstyle\bot}}{{\sim}}\;\left\{\begin{array}[]{ll}\chi^{2}_{2\alpha_{i}+2N_{i}}&,\;i=1,\,\ldots\,,D\\ \\ \chi^{2}_{2\alpha_{D+1}+2\left(N^{+}-\sum_{r=1}^{D}N_{r}\right)}&,\;i=D+1\end{array}\right.\;;

hence, by Eq. (47):

\left.Z^{\prime+}\,\right|\,\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)\stackrel{{\scriptstyle d}}{{=}}\left.Z^{\prime+}\,\right|\,N^{+}\sim\chi^{2}_{2\alpha^{+}+2N^{+}}\,,

(102)

where $\stackrel{{\scriptstyle d}}{{=}}$ denotes equality in distribution. By Property 2.1, $\underline{X}^{\prime\prime}$ and $Z^{\prime+}$ are independent conditionally on $\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)$ and therefore, given this latter, the joint distribution of $\left(\underline{X}^{\prime\prime},Z^{\prime+}\right)$ factores into the marginal distributions of $\underline{X}^{\prime\prime}$ and $Z^{\prime+}$ . That said, the proof follows by noting that the joint density of $\left.\left(\underline{X}^{\prime\prime},Z^{\prime+}\right)\,\right|\,N^{+}$ factores into the marginal densities of $\left.\underline{X}^{\prime\prime}\,\right|\,N^{+}$ and $\left.Z^{\prime+}\,\right|\,N^{+}$ ; indeed, under Eq. (102):

$\displaystyle f_{\left.\left(\underline{X}^{\prime\prime},\,Z^{\prime+}\right)\,\right\|\,N^{+}}\left(\underline{x},z\right)$	$\displaystyle=$	$\displaystyle\sum_{j^{+}\leq N^{+}}\left[\,\Pr\left(\left.\,\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}\right)\cdot\right.$
	$\displaystyle\cdot$	$\displaystyle\left.f_{\left.\left(\underline{X}^{\prime\prime},\,Z^{\prime+}\right)\,\right\|\,\left(\,N_{1},\,\ldots\,,N_{D},N^{+}\,\right)}\left(\underline{x},z;\,\left(j_{1},\,\ldots\,,j_{D}\right)\,\right)\,\right]=$
	$\displaystyle=$	$\displaystyle\sum_{j^{+}\leq N^{+}}\left[\,\Pr\left(\left.\,\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}\right)\cdot\right.$
	$\displaystyle\cdot$	$\displaystyle\left.f_{\left.\underline{X}^{\prime\prime}\,\right\|\,\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)}\left(\underline{x},\left(j_{1},\,\ldots\,,j_{D}\right)\,\right)\,\right]\cdot f_{\left.Z^{\prime+}\,\right\|\,\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)}\left(z\right)$
	$\displaystyle=$	$\displaystyle f_{\left.\underline{X}^{\prime\prime}\,\right\|\,N^{+}}\left(\underline{x}\right)\cdot f_{\left.Z^{\prime+}\,\right\|\,N^{+}}\left(z\right)\,,$

where $\left(j_{1},\,\ldots\,,j_{D}\right)\in\mathbb{N}_{0}^{\,D}$ , $j^{+}=\sum_{i=1}^{D}j_{i}$ and the density $f_{\left.X^{\prime\prime}\,\right|\,N^{+}}$ is given in Eq. (101). ∎

3.4 Density plots

In this subsection we shall focus on the case $D=2$ in order to make it possible the graphical depiction of the joint density of the model of interest.

A key advantage of the $\mbox{{CNcDir}}^{\,2}$ distribution over the bivariate Dirichlet is the much larger variety of shapes reachable by the density of the former on varying the non-centrality parameters. In this regard, some significant plots of the bivariate CNcDir density are displayed in Figure 1 for selected values of the parameter vector.

Refer to caption — Figure 1: Plots of the density of $\left(X^{\prime\prime}_{1},X^{\prime\prime}_{2}\right)\sim\mbox{{CNcDir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ for selected values of the parameter vector.

Specifically, the main attraction of the CNcDir² distribution lies in the following feature. When at least two shape parameters are set equal to one, the three non-centrality parameters control the height of the density at the vertices of the unit simplex by allowing the limits of this latter to have arbitrary positive and finite values. That said, the limiting values of the bivariate Conditional Non-central Dirichlet density are established herein.

Proposition 3.14 (Limiting values of the bivariate density)

Let $\left(X^{\prime\prime}_{1},X^{\prime\prime}_{2}\right)\sim\mbox{{CNcDir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ . Then:

$\displaystyle\mbox{{CNcDir}}^{\,2}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)\quad\rightarrow$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(1,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{2}=\alpha_{3}=1\,:&\alpha_{1}\left(\alpha_{1}+1\right)\frac{{}_{0}F_{1}\left(\alpha_{1};\frac{\lambda_{1}}{4}\right)}{{}_{0}F_{1}\left(\alpha_{1}+2;\frac{\lambda^{+}}{4}\right)}\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{2}+\alpha_{3}>2\,:&0\end{array}\right.$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,1\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{3}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{3}=1\,:&\alpha_{2}\left(\alpha_{2}+1\right)\frac{{}_{0}F_{1}\left(\alpha_{2};\frac{\lambda_{2}}{4}\right)}{{}_{0}F_{1}\left(\alpha_{2}+2;\frac{\lambda^{+}}{4}\right)}\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{3}>2\,:&0\end{array}\right.$
	$\displaystyle\stackrel{{\scriptstyle\left(x_{1},\,x_{2}\right)\,\rightarrow\,\left(0,0\right)}}{{\rightarrow}}$	$\displaystyle\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}+\alpha_{2}<2\,:&+\,\infty\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}=2\,:&\left\{\begin{array}[]{ll}\mbox{{if} }\,\alpha_{1}=\alpha_{2}=1\,:&\alpha_{3}\left(\alpha_{3}+1\right)\frac{{}_{0}F_{1}\left(\alpha_{3};\frac{\lambda_{3}}{4}\right)}{{}_{0}F_{1}\left(\alpha_{3}+2;\frac{\lambda^{+}}{4}\right)}\\ \mbox{{otherwise:} }&\not\exists\end{array}\right.\\ \mbox{{if} }\,\alpha_{1}+\alpha_{2}>2\,:&0\end{array}\right.$

where the symbol $\not\exists\,$ indicates that the corresponding limit does not exist because depends on the direction followed to reach the accumulation point under consideration.

Proof.

The proof ensues from Eqs. (33) and (85) by observing that, in the notation of Eq. (7), ${}_{0}F_{1}\left(b;0\right)=1$ . ∎

The surface representing the density function of the CNcDir² distribution is depicted in Figure 2 for shape parameters equal one and selected values of the non-centrality parameters. The above mentioned feature of the density at study can be clearly detected from these plots.

All the elements that enable us to affirm that the CNcDir distribution preserves the applicative potential of the NcDir one are finally at hand. Indeed, there exist deep similarities between the varieties of shapes reachable by the two densities in the case that, under the same values for the shape parameters, the non-centrality parameters of the Conditional model are assigned much higher values than those used for the existing one. In light of the mixture representations of the two models, in fact, the differences between such distributions can be entirely explained by the differences existing between their respective mixing distributions. In light of Eq. (96), these latter are completely due to the differences between the sums of their components $N^{+}$ and $M^{+}$ , the probability mass functions of which depend on $\underline{\lambda}$ by means of $\lambda^{+}$ . By virtue of the arguments put forward in [Orsi (2021)] (Section 4) with reference to the comparison between the Conditional and the standard Doubly Non-central Beta distributions, the link between the $\lambda^{+}$ values of $N^{+}$ and $M^{+}$ can be understood by comparing their modes as functions of $\lambda^{+}$ . Specifically, for a given value of the mode of $M^{+}$ , one can find the value $\lambda^{+}$ relative to $N^{+}$ which produces similar mode; as shown in the above cited reference, this $\lambda^{+}$ value is much larger than the one relative to $M^{+}$ and gets even larger when $\alpha^{+}$ increases.

3.5 Moments

For the sake of simplicity, in all the present subsection we shall focus on the case $D=2$ ; however, by virtue of Eq. (88), the results derived herein also apply for any two-dimensional marginal of the CNcDir ${}^{\,D}$ model with $D>2$ .

By analogy with the form of the density in Eq. (84) and by Eq. (34), for every $r_{1},\,r_{2}\in\mathbb{N}$ , the mixed raw moment of order $(r_{1},r_{2})$ of $\left(X^{\prime\prime}_{1},X^{\prime\prime}_{2}\right)\sim\mbox{{CNcDir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ can be stated as

	$\displaystyle\mathbb{E}\left[\left(X_{1}^{\prime\prime}\right)^{\,r_{1}}\left(X_{2}^{\prime\prime}\right)^{\,r_{2}}\,\right]=\qquad\qquad\qquad\qquad\qquad\qquad\left\{\begin{array}[]{l}r^{+}=r_{1}+r_{2}\\ j^{+}=j_{1}+j_{2}+j_{3}\end{array}\right.$				(111)
		$\displaystyle=$	$\displaystyle\sum_{j_{1},\,j_{2},\,j_{3}\,=\,0}^{+\infty}\left\{\,\Pr\left[\;\left(N_{1},N_{2},N_{3}\right)=\left(j_{1},j_{2},j_{3}\right)\;\right]\cdot\frac{\left(\alpha_{1}+j_{1}\right)_{r_{1}}\,\left(\alpha_{2}+j_{2}\right)_{r_{2}}}{\left(\alpha^{+}+j^{+}\right)_{r^{+}}}\,\right\}\,,$		(112)

i.e. as the multiple infinite series of the mixed raw moments of order $(r_{1},r_{2})$ of the $\mbox{Dir}^{\,2}\left(\alpha_{1}+j_{1},\alpha_{2}+j_{2},\alpha_{3}+j_{3}\right)$ distributions weighted by the $\mbox{MW}^{\,3}\left(\alpha^{+},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ probabilities. Moreover, by Eq. (3), Eq. (112) can be equivalently expressed as the following doubly infinite sum of generalized hypergeometric functions with one numerator parameter and two denominator parameters:

	$\displaystyle\mathbb{E}\left[\left(X_{1}^{\prime\prime}\right)^{\,r_{1}}\left(X_{2}^{\prime\prime}\right)^{\,r_{2}}\,\right]=\frac{\left(\alpha_{1}\right)_{r_{1}}\left(\alpha_{2}\right)_{r_{2}}}{\left(\alpha^{+}\right)_{r^{+}}}\cdot\sum_{j_{2},\,j_{3}\,=\,0}^{+\infty}\left[\,\frac{\left(\alpha_{2}+r_{2}\right)_{j_{2}}}{\left(\alpha_{2}\right)_{j_{2}}\,\left(\alpha^{+}+r^{+}\right)_{j_{2}+j_{3}}}\right.\cdot$				(113)
		$\displaystyle\cdot$	$\displaystyle\left.\frac{\left(\frac{\lambda_{2}}{4}\right)^{j_{2}}}{j_{2}!}\frac{\left(\frac{\lambda_{3}}{4}\right)^{j_{3}}}{j_{3}!}\cdot\,\frac{{}_{1}^{\,}F^{\,}_{2}\left(\alpha_{1}+r_{1};\alpha_{1},\alpha^{+}+r^{+}+j_{2}+j_{3};\frac{\lambda_{1}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,\right]\,.$		(113)

Unfortunately, the above formulas are computationally cumbersome. Hence, in the following we shall provide an interesting general expression for the quantity at study that allows the computation of this latter to be reduced from the doubly infinite sum in Eq. (113) to a surprisingly simple form given by a doubly finite sum.

Proposition 3.15 (Mixed raw moments of the bivariate distribution)

The mixed raw moment of order $(r_{1},r_{2})$ of $\left(X^{\prime\prime}_{1},X^{\prime\prime}_{2}\right)\sim\mbox{{CNcDir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ admits the following expression:

	$\displaystyle\mathbb{E}\left[\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right]\quad=\quad\frac{\left(\alpha_{1}\right)_{r_{1}}\left(\alpha_{2}\right)_{r_{2}}}{\left(\alpha^{+}\right)_{r^{+}}}\;\cdot\qquad\quad\left\{\begin{array}[]{l}r_{1},r_{2}\in\mathbb{N},\,r^{+}=r_{1}+r_{2}\\ j^{+}=j_{1}+j_{2}\end{array}\right.$				(116)
		$\displaystyle\cdot$	$\displaystyle\sum_{j_{1}=0}^{r_{1}}\sum_{j_{2}=0}^{r_{2}}\frac{{r_{1}\choose j_{1}}{r_{2}\choose j_{2}}\left(\frac{\lambda_{1}}{4}\right)^{j_{1}}\left(\frac{\lambda_{2}}{4}\right)^{j_{2}}}{\left(\alpha_{1}\right)_{j_{1}}\left(\alpha_{2}\right)_{j_{2}}\left(\alpha^{+}+r^{+}\right)_{j^{+}}}\frac{{}_{0}F_{1}\left(\alpha^{+}+r^{+}+j^{+};\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,.$		(117)

Proof.

Let $\left(L_{1},L_{2}\right)$ have a $\mbox{Multinomial}^{\,2}\left(N^{+},\theta_{1},\theta_{2}\right)$ distribution conditionally on $N^{+}$ with $\theta_{i}=\lambda_{i}\,/\,\lambda^{+}$ , $i=1,2$ . By Eq. (101), one has:

$\displaystyle\mathbb{E}\left[\left.\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right\|\,N^{+}\right]=$			(118)
	$\displaystyle=$	$\displaystyle\int_{\left(x_{1},\,x_{2}\right)\,\in\,\mathcal{S}^{2}}x_{1}^{\,r_{1}}x_{2}^{\,r_{2}}\cdot f_{\left.\left(X^{\prime\prime}_{1},\,X^{\prime\prime}_{2}\right)\,\right\|\,N^{+}}\left(x_{1},x_{2};\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)\,dx_{1}\,dx_{2}=$
	$\displaystyle=$	$\displaystyle\sum_{l_{1}=0}^{N^{+}}\sum_{l_{2}=0}^{N^{+}-\,l_{1}}\frac{\left(\alpha_{1}+l_{1}\right)_{r_{1}}\left(\alpha_{2}+l_{2}\right)_{r_{2}}}{\left(\alpha^{+}+N^{+}\right)_{r^{+}}}\,{N^{+}\choose l_{1}\;l_{2}}\;\theta_{1}^{\,l_{1}}\,\theta_{2}^{\,l_{2}}\left(1-\theta_{1}-\theta_{2}\right)^{N^{+}-\,l_{1}-l_{2}}=$
	$\displaystyle=$	$\displaystyle\frac{\mathbb{E}\left[\left.\left(\alpha_{1}+L_{1}\right)_{r_{1}}\left(\alpha_{2}+L_{2}\right)_{r_{2}}\right\|\,N^{+}\right]}{\left(\alpha^{+}+N^{+}\right)_{r^{+}}}\,;$

in light of Eq. (5):

\left(\alpha_{i}+L_{i}\right)_{r_{i}}=\left[\left(\alpha_{i}-1\right)+\left(L_{i}+1\right)\right]_{r_{i}}=\sum_{j_{i}=0}^{r_{i}}{r_{i}\choose j_{i}}\left(\alpha_{i}-1\right)_{r_{i}-j_{i}}\left(L_{i}+1\right)_{j_{i}}\,,\quad\;i=1,2\,,

so that one obtains:

	$\displaystyle\mathbb{E}\left[\left.\left(\alpha_{1}+L_{1}\right)_{r_{1}}\left(\alpha_{2}+L_{2}\right)_{r_{2}}\right\|\,N^{+}\right]=$
		$\displaystyle=$	$\displaystyle\sum_{j_{1}=0}^{r_{1}}\sum_{j_{2}=0}^{r_{2}}{r_{1}\choose j_{1}}{r_{2}\choose j_{2}}\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\mathbb{E}\left[\left.\left(L_{1}+1\right)_{j_{1}}\left(L_{2}+1\right)_{j_{2}}\right\|\,N^{+}\right]\,,$

where:

	$\displaystyle\mathbb{E}\left[\left.\left(L_{1}+1\right)_{j_{1}}\left(L_{2}+1\right)_{j_{2}}\right\|\,N^{+}\right]=$				(120)
		$\displaystyle=$	$\displaystyle\sum_{l_{1}=0}^{N^{+}}\sum_{l_{2}=0}^{N^{+}-\,l_{1}}\left(l_{1}+1\right)_{j_{1}}\left(l_{2}+1\right)_{j_{2}}\,{N^{+}\choose l_{1}\;l_{2}}\;\theta_{1}^{\,l_{1}}\,\theta_{2}^{\,l_{2}}\left(1-\theta_{1}-\theta_{2}\right)^{N^{+}-\,l_{1}-\,l_{2}}\,.$		(120)

By Eq. (2), for every $j_{i}=0,\,\ldots\,,r_{i}$ , $i=1,2$ :

\left(l_{i}+1\right)_{j_{i}}=\frac{\Gamma\left(l_{i}+j_{i}+1\right)}{\Gamma\left(l_{i}+1\right)}=\frac{\left(l_{i}+j_{i}\right)!}{l_{i}!}={l_{i}+j_{i}\choose l_{i}}\,j_{i}!\;;

(121)

under Eq. (121), Eq. (120) can be thus rewritten as follows:

$\displaystyle\mathbb{E}\left[\left.\left(L_{1}+1\right)_{j_{1}}\left(L_{2}+1\right)_{j_{2}}\right\|\,N^{+}\right]=$			(122)
	$\displaystyle=$	$\displaystyle j_{1}!\;j_{2}!\,\sum_{l_{1}=0}^{N^{+}}{l_{1}+j_{1}\choose l_{1}}{N^{+}\choose l_{1}}\left(1-\theta_{1}\right)^{N^{+}-\,l_{1}}\theta_{1}^{\,l_{1}}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{l_{2}=0}^{N^{+}-\,l_{1}}{l_{2}+j_{2}\choose l_{2}}{N^{+}-\,l_{1}\choose l_{2}}\left(1-\frac{\theta_{2}}{1-\theta_{1}}\right)^{\left(N^{+}-\,l_{1}\right)-\,l_{2}}\left(\frac{\theta_{2}}{1-\theta_{1}}\right)^{\,l_{2}}\,.$

In carrying out the prove, reference must be made to Ljunggren’s Identity, namely

\sum_{k=0}^{n}{\alpha+k\choose k}{n\choose k}\left(x-y\right)^{n-k}y^{k}=\sum_{k=0}^{n}{\alpha\choose k}{n\choose k}\,x^{n-k}\,y^{k}\,,

(123)

which is (3.18) in [Gould (1972)]. By the special case of Eq. (123) where $k=l_{2}$ , $n=N^{+}-\,l_{1}$ , $\alpha=j_{2}$ , $x=1$ , $y=\theta_{2}/(1-\theta_{1})$ , the final sum in Eq. (122) can be equivalently expressed as

\sum_{l_{2}=0}^{N^{+}-\,l_{1}}{j_{2}\choose l_{2}}\,{N^{+}-\,l_{1}\choose l_{2}}\left(\frac{\theta_{2}}{1-\theta_{1}}\right)^{l_{2}}\,,

so that, for every $j_{i}=0,\,\ldots\,,r_{i}$ , $i=1,2$ , Eq. (122) can be restated as

$\displaystyle\mathbb{E}\left[\left.\left(L_{1}+1\right)_{j_{1}}\left(L_{2}+1\right)_{j_{2}}\right\|\,N^{+}\right]=$			(124)
	$\displaystyle=$	$\displaystyle j_{1}!\,j_{2}!\,\sum_{l_{1}=0}^{N^{+}}\sum_{l_{2}=0}^{N^{+}-\,l_{1}}{j_{2}\choose l_{2}}{l_{1}+j_{1}\choose l_{1}}{N^{+}\choose l_{1}\;l_{2}}\,\theta_{1}^{\,l_{1}}\,\theta_{2}^{\,l_{2}}\left(1-\theta_{1}\right)^{N^{+}-\,l_{1}-\,l_{2}}=$
	$\displaystyle=$	$\displaystyle j_{1}!\,j_{2}!\,\sum_{l_{2}=0}^{N^{+}}{j_{2}\choose l_{2}}\,{N^{+}\choose l_{2}}\,\theta_{2}^{\,l_{2}}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{l_{1}=0}^{N^{+}-\,l_{2}}{l_{1}+j_{1}\choose l_{1}}\,{N^{+}-l_{2}\choose l_{1}}\left(1-\theta_{1}\right)^{(N^{+}-\,l_{2})-\,l_{1}}\,\theta_{1}^{\,l_{1}}\,.$

Under the special case of Eq. (123) where $k=l_{1}$ , $n=N^{+}-\,l_{2}$ , $\alpha=j_{1}$ , $x=1$ , $y=\theta_{1}$ , Eq. (124) can be rewritten in the form of

	$\displaystyle\mathbb{E}\left[\left.\left(L_{1}+1\right)_{j_{1}}\left(L_{2}+1\right)_{j_{2}}\right\|\,N^{+}\right]=$				(125)
		$\displaystyle=$	$\displaystyle j_{1}!\,j_{2}!\,\sum_{l_{2}=0}^{N^{+}}{j_{2}\choose l_{2}}{N^{+}\choose l_{2}}\,\theta_{2}^{\,l_{2}}\sum_{l_{1}=0}^{N^{+}-\,l_{2}}{j_{1}\choose l_{1}}{N^{+}-\,l_{2}\choose l_{1}}\,\theta_{1}^{\,l_{1}}\,,$		(125)

so that, by Eqs. (118), (LABEL:eq:mom.dim1b) and (125), conditionally on $N^{+}$ , the mixed raw moment of interest takes on the following expression:

$\displaystyle\mathbb{E}\left[\left.\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right\|\,N^{+}\right]=$			(126)
	$\displaystyle=$	$\displaystyle\frac{1}{\left(\alpha^{+}+N^{+}\right)_{r^{+}}}\sum_{j_{1}=0}^{r_{1}}\sum_{j_{2}=0}^{r_{2}}{r_{1}\choose j_{1}}j_{1}!\,{r_{2}\choose j_{2}}j_{2}!\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{l_{2}=0}^{N^{+}}{j_{2}\choose l_{2}}{N^{+}\choose l_{2}}\,\theta_{2}^{\,l_{2}}\sum_{l_{1}=0}^{N^{+}-\,l_{2}}{j_{1}\choose l_{1}}{N^{+}-\,l_{2}\choose l_{1}}\,\theta_{1}^{\,l_{1}}\,;$

applying the law of iterated expectations to Eq. (126) finally leads to:

$\displaystyle\mathbb{E}\left[\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right]=$			(127)
	$\displaystyle=$	$\displaystyle\frac{\Gamma\left(\alpha^{+}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\sum_{j_{1}=0}^{r_{1}}\sum_{j_{2}=0}^{r_{2}}{r_{1}\choose j_{1}}j_{1}!\,{r_{2}\choose j_{2}}j_{2}!\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{n=0}^{+\infty}\frac{\left(\frac{\lambda^{+}}{4}\right)^{n}}{n!\,\Gamma\left(\alpha^{+}+r^{+}+n\right)}\sum_{l_{2}=0}^{n}{j_{2}\choose l_{2}}{n\choose l_{2}}\,\theta_{2}^{\,l_{2}}\sum_{l_{1}=0}^{n-\,l_{2}}{j_{1}\choose l_{1}}{n-\,l_{2}\choose l_{1}}\,\theta_{1}^{\,l_{1}}\,.$

Since ${j_{i}\choose l_{i}}=0$ for $l_{i}>j_{i}$ , $i=1,2$ , one has:

l_{2}\leq j_{2}\leq n\;\Rightarrow\;\left\{\begin{array}[]{l}l_{2}=0,\,\ldots\,,j_{2}\\ \\ n=l_{2},\,\ldots\,,+\infty\end{array}\right.\,\quad l_{1}\leq j_{1}\leq n-l_{2}\;\Rightarrow\;\left\{\begin{array}[]{l}l_{1}=0,\,\ldots\,,j_{1}\\ \\ n=l_{1}+l_{2},\,\ldots\,,+\infty\end{array}\right.

so that Eq. (127) is equivalent to:

$\displaystyle\mathbb{E}\left[\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right]=\frac{\Gamma\left(\alpha^{+}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{j_{1}=0}^{r_{1}}{r_{1}\choose j_{1}}j_{1}!\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\,\sum_{l_{1}=0}^{j_{1}}\frac{\theta_{1}^{\,l_{1}}}{l_{1}!}\,{j_{1}\choose l_{1}}\,\sum_{j_{2}=0}^{r_{2}}{r_{2}\choose j_{2}}j_{2}!\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{l_{2}=0}^{j_{2}}\frac{\theta_{2}^{\,l_{2}}}{l_{2}!}\,{j_{2}\choose l_{2}}\sum_{n=l_{1}+l_{2}}^{+\infty}\frac{\left(\frac{\lambda^{+}}{4}\right)^{n}}{\left(n-\,l_{1}-\,l_{2}\right)!\;\Gamma\left(\alpha^{+}+r^{+}+n\right)}\,;$

then, by Eqs. (2) and (7), setting $k=n-(l_{1}+l_{2})\Leftrightarrow n=k+(l_{1}+l_{2})$ yields:

		$\displaystyle=$	$\displaystyle\Gamma\left(\alpha^{+}\right)\,\sum_{j_{1}=0}^{r_{1}}{r_{1}\choose j_{1}}j_{1}!\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\,\sum_{l_{1}=0}^{j_{1}}\frac{\left(\frac{\lambda_{1}}{4}\right)^{l_{1}}{j_{1}\choose l_{1}}}{l_{1}!}\,\sum_{j_{2}=0}^{r_{2}}{r_{2}\choose j_{2}}j_{2}!\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\cdot$
		$\displaystyle\cdot$	$\displaystyle\sum_{l_{2}=0}^{j_{2}}\frac{\left(\frac{\lambda_{2}}{4}\right)^{l_{2}}\,{j_{2}\choose l_{2}}}{l_{2}!\;\Gamma\left(\alpha^{+}+r^{+}+l_{1}+l_{2}\right)}\,\sum_{k=0}^{+\infty}\frac{\left(\frac{\lambda^{+}}{4}\right)^{k}}{k!\,\left(\alpha^{+}+r^{+}+l_{1}+l_{2}\right)_{k}}\frac{1}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}=$
		$\displaystyle=$	$\displaystyle\Gamma\left(\alpha^{+}\right)\sum_{j_{1}=0}^{r_{1}}{r_{1}\choose j_{1}}j_{1}!\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\,\sum_{l_{1}=0}^{j_{1}}\frac{\left(\frac{\lambda_{1}}{4}\right)^{l_{1}}{j_{1}\choose l_{1}}}{l_{1}!}\,\sum_{j_{2}=0}^{r_{2}}{r_{2}\choose j_{2}}j_{2}!\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\cdot$
		$\displaystyle\cdot$	$\displaystyle\sum_{l_{2}=0}^{j_{2}}\frac{\left(\frac{\lambda_{2}}{4}\right)^{l_{2}}\,{j_{2}\choose l_{2}}}{l_{2}!\;\Gamma\left(\alpha^{+}+r^{+}+l_{1}+l_{2}\right)}\,\frac{{}_{0}F_{1}\left(\alpha^{+}+r^{+}+l_{1}+l_{2};\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}$

and, by observing that:

\left\{\begin{array}[]{l}j_{i}=0,\,\ldots\,,r_{i}\\ \\ l_{i}=0,\,\ldots\,,j_{i}\end{array}\right.\quad\Rightarrow\quad 0\,\leq\,l_{i}\,\leq\,j_{i}\,\leq\,r_{i}\quad\Rightarrow\quad\left\{\begin{array}[]{l}l_{i}=0,\,\ldots\,,r_{i}\\ \\ j_{i}=l_{i},\,\ldots\,,r_{i}\end{array}\right.\,,

one can obtain:

$\displaystyle\mathbb{E}\left[\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right]=$			(128)
	$\displaystyle=$	$\displaystyle\Gamma\left(\alpha^{+}\right)\sum_{l_{1}=0}^{r_{1}}\frac{\left(\frac{\lambda_{1}}{4}\right)^{l_{1}}}{\left(l_{1}!\right)^{2}}\,\sum_{l_{2}=0}^{r_{2}}\frac{\left(\frac{\lambda_{2}}{4}\right)^{l_{2}}\,_{0}F_{1}\left(\alpha^{+}+r^{+}+l_{1}+l_{2};\frac{\lambda^{+}}{4}\right)}{\left(l_{2}!\right)^{2}\,\Gamma\left(\alpha^{+}+r^{+}+l_{1}+l_{2}\right)\,_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\cdot$
	$\displaystyle\cdot$	$\displaystyle\sum_{j_{1}=l_{1}}^{r_{1}}{r_{1}\choose j_{1}}\left(\alpha_{1}-1\right)_{r_{1}-j_{1}}\frac{\left(j_{1}!\right)^{2}}{\left(j_{1}-l_{1}\right)!}\,\sum_{j_{2}=l_{2}}^{r_{2}}{r_{2}\choose j_{2}}\left(\alpha_{2}-1\right)_{r_{2}-j_{2}}\frac{\left(j_{2}!\right)^{2}}{\left(j_{2}-l_{2}\right)!}\,.$

By setting $p_{i}=j_{i}-l_{i}\Leftrightarrow j_{i}=p_{i}+l_{i}$ , $i=1,2$ and by Eq. (5), each of the two final sums in Eq. (128) turns out to be tantamount to:

$\displaystyle\sum_{j_{i}=l_{i}}^{r_{i}}{r_{i}\choose j_{i}}\left(\alpha_{i}-1\right)_{r_{i}-j_{i}}\frac{\left(j_{i}!\right)^{2}}{\left(j_{i}-l_{i}\right)!}=$			(129)
	$\displaystyle=$	$\displaystyle\sum_{p_{i}=0}^{r_{i}-l_{i}}{r_{i}\choose p_{i}+l_{i}}\left(\alpha_{i}-1\right)_{r_{i}-p_{i}-l_{i}}\frac{\left[\left(p_{i}+l_{i}\right)!\right]^{2}}{p_{i}!}=$
	$\displaystyle=$	$\displaystyle\frac{r_{i}!}{\left(r_{i}-l_{i}\right)!}\sum_{p_{i}=0}^{r_{i}-l_{i}}{r_{i}-l_{i}\choose p_{i}}\left(\alpha_{i}-1\right)_{r_{i}-p_{i}-l_{i}}\,\left(l_{i}+p_{i}\right)!=$
	$\displaystyle=$	$\displaystyle\frac{r_{i}!\,l_{i}!}{\left(r_{i}-l_{i}\right)!}\sum_{p_{i}=0}^{r_{i}-l_{i}}{r_{i}-l_{i}\choose p_{i}}\left(\alpha_{i}-1\right)_{r_{i}-p_{i}-l_{i}}\,\left(l_{i}+1\right)_{p_{i}}=$
	$\displaystyle=$	$\displaystyle\frac{r_{i}!\,l_{i}!}{\left(r_{i}-l_{i}\right)!}\left[\left(\alpha_{i}-1\right)+\left(l_{i}+1\right)\right]_{r_{i}-l_{i}}=\frac{r_{i}!\,l_{i}!}{\left(r_{i}-l_{i}\right)!}\left(\alpha_{i}+l_{i}\right)_{r_{i}-l_{i}}\,;$

hence, by Eq. (129), Eq. (128) can be stated in the following form:

$\displaystyle\mathbb{E}\left[\left(X^{\prime\prime}_{1}\right)^{r_{1}}\left(X^{\prime\prime}_{2}\right)^{r_{2}}\right]=$			(130)
	$\displaystyle=$	$\displaystyle\Gamma\left(\alpha^{+}\right)\cdot\sum_{l_{1}=0}^{r_{1}}\sum_{l_{2}=0}^{r_{2}}\left[\frac{{r_{1}\choose l_{1}}\,{r_{2}\choose l_{2}}\left(\frac{\lambda_{1}}{4}\right)^{l_{1}}\left(\frac{\lambda_{2}}{4}\right)^{l_{2}}\,\left(\alpha_{1}+l_{1}\right)_{r_{1}-l_{1}}\,\left(\alpha_{2}+l_{2}\right)_{r_{2}-l_{2}}}{\Gamma\left(\alpha^{+}+r^{+}+l_{1}+l_{2}\right)}\cdot\right.$
	$\displaystyle\cdot$	$\displaystyle\left.\frac{{}_{0}F_{1}\left(\alpha^{+}+r^{+}+l_{1}+l_{2};\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\right]\,.$

Eq. (130) can be finally exhibited in the form reported in Eq. (117) by noting that:

	$\displaystyle\Gamma\left(\alpha^{+}+r^{+}+l_{1}+l_{2}\right)$	$\displaystyle=$	$\displaystyle\Gamma\left(\alpha^{+}\right)\,\left(\alpha^{+}\right)_{r^{+}}\,\left(\alpha^{+}+r^{+}\right)_{l_{1}+l_{2}}\,,$
	$\displaystyle\left(\alpha_{i}+l_{i}\right)_{r_{i}-l_{i}}$	$\displaystyle=$	$\displaystyle\frac{\left(\alpha_{i}\right)_{r_{i}}}{\left(\alpha_{i}\right)_{l_{i}}}\,\qquad i=1,2$

in light of Eqs. (2), (3) and (4). ∎

Hence, the formula of the mixed raw moment of order $(1,1)$ of the $\mbox{{CNcDir}}^{\,2}$ distribution can be obtained by taking $r_{1}=r_{2}=1$ in Eq. (117) as follows:

	$\displaystyle\mathbb{E}\left(X^{\prime\prime}_{1}\,X^{\prime\prime}_{2}\right)$	$\displaystyle=$	$\displaystyle\frac{\alpha_{1}\,\alpha_{2}}{\left(\alpha^{+}\right)_{2}}\frac{{}_{0}F_{1}\left(\alpha^{+}+2;\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}+\frac{\alpha_{1}\,\frac{\lambda_{2}}{4}+\alpha_{2}\,\frac{\lambda_{1}}{4}}{\left(\alpha^{+}\right)_{3}}\,\frac{{}_{0}F_{1}\left(\alpha^{+}+3;\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}+$		(131)
		$\displaystyle+$	$\displaystyle\frac{\frac{\lambda_{1}}{4}\,\frac{\lambda_{2}}{4}}{\left(\alpha^{+}\right)_{4}}\frac{{}_{0}F_{1}\left(\alpha^{+}+4;\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,.$		(131)

This latter can be algebraically manipulated with the aim to reduce the number of distinct functions $\,{}_{0}F_{1}$ appearing in it. More precisely, in carrying out the computations, reference must be made to the recurrence identity

{}_{0}F_{1}\left(\alpha^{+}+2;\frac{\lambda^{+}}{4}\right)=\,_{0}F_{1}\left(\alpha^{+}+3;\frac{\lambda^{+}}{4}\right)+\frac{\frac{\lambda^{+}}{4}}{\left(\alpha^{+}+2\right)\left(\alpha^{+}+3\right)}\,_{0}F_{1}\left(\alpha^{+}+4;\frac{\lambda^{+}}{4}\right)\,,

which holds for consecutive values of the denominator parameter of $\,{}_{0}F_{1}$ and is easy to verify upon noting that, by Eq. (3), $\left(\alpha^{+}+2\right)_{i}\left(\alpha^{+}+2+i\right)=\left(\alpha^{+}+2\right)\left(\alpha^{+}+3\right)_{i}$ for every $i\in\mathbb{N}\cup\{0\}$ . Specifically, the aforementioned improvement of Eq. (131) depends on three functions ${}_{0}F_{1}$ instead of four and is given by

	$\displaystyle\mathbb{E}\left(X^{\prime\prime}_{1}\,X^{\prime\prime}_{2}\right)$	$\displaystyle=$	$\displaystyle\frac{\alpha_{1}\,\alpha_{2}+\frac{\lambda_{1}\,\lambda_{2}}{4\,\lambda^{+}}}{\left(\alpha^{+}\right)_{2}}\frac{{}_{0}F_{1}\left(\alpha^{+}+2;\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}+$
		$\displaystyle+$	$\displaystyle\frac{\alpha_{1}\,\frac{\lambda_{2}}{4}+\alpha_{2}\,\frac{\lambda_{1}}{4}-\frac{\lambda_{1}\,\lambda_{2}}{4\,\lambda^{+}}\left(\alpha^{+}+2\right)}{\left(\alpha^{+}\right)_{3}}\,\frac{{}_{0}F_{1}\left(\alpha^{+}+3;\frac{\lambda^{+}}{4}\right)}{{}_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)}\,.$

3.6 Identifiability

Before facing the applicative issues included in the next section, in the present subsection we briefly discuss the identifiability of the CNcDir model.

Proposition 3.16 (Identifiability)

Let $\underline{X}^{\prime\prime}=\left(X^{\prime\prime}_{1},\,\ldots\,,X^{\prime\prime}_{D}\right)\sim\mbox{{CNcDir}}^{\,D}\left(\underline{\alpha},\underline{\lambda}\right)$ and $\underline{\tilde{X}}^{\prime\prime}=(\tilde{X}^{\prime\prime}_{1},\,\ldots\,,\tilde{X}^{\prime\prime}_{D})\sim\mbox{{CNcDir}}^{\,D}(\underline{\tilde{\alpha}},$ $\underline{\tilde{\lambda}})$ where $\underline{\alpha}=\left(\alpha_{1},\,\ldots\,,\alpha_{D+1}\right)$ , $\underline{\tilde{\alpha}}=\left(\tilde{\alpha}_{1},\,\ldots\,,\tilde{\alpha}_{D+1}\right)$ , $\underline{\lambda}=\left(\lambda_{1},\,\ldots\,,\lambda_{D+1}\right)$ and $\underline{\tilde{\lambda}}=(\tilde{\lambda}_{1},\,\ldots\,,\tilde{\lambda}_{D+1})$ . Then, $\underline{X}^{\prime\prime}\sim\underline{\tilde{X}}^{\prime\prime}$ if and only if $\underline{\alpha}=\underline{\tilde{\alpha}}$ and $\underline{\lambda}=\underline{\tilde{\lambda}}$ .

Proof.

We need to show that if $\underline{X}^{\prime\prime}\sim\underline{\tilde{X}}^{\prime\prime}$ then $\underline{\alpha}=\underline{\tilde{\alpha}}$ and $\underline{\lambda}=\underline{\tilde{\lambda}}$ , the converse being obvious. Clearly, if $\underline{X}^{\prime\prime}\sim\underline{\tilde{X}}^{\prime\prime}$ , then $X^{\prime\prime}_{i}\sim\tilde{X}^{\prime\prime}_{i}$ for every $i=1,\,\ldots\,,D$ . By Eq. (87), $X^{\prime\prime}_{i}\sim\mbox{{CDNcB}}(\alpha_{i},\alpha^{+}-\alpha_{i},\lambda_{i},\lambda^{+}-\lambda_{i})$ and $\tilde{X}^{\prime\prime}_{i}\sim\mbox{{CDNcB}}(\tilde{\alpha}_{i},\tilde{\alpha}^{+}-\tilde{\alpha}_{i},\tilde{\lambda}_{i},\tilde{\lambda}^{+}-\tilde{\lambda}_{i})$ ; therefore, the proof follows from the identifiability of the CDNcB model ([Orsi (2021)], Section 7). ∎

4 Applications

The potential of the CNcDir model is now illustrated by means of an application to a real data set. To this end we turned our attention to the Cartesian coordinates of the locations of 584 longleaf pine trees (Pinus palustris) placed in a $200\times 200$ metre region within Southern Georgia (USA). This data set was originally collected and analyzed by [Platt et al. (1988)], was later quoted in [Rathbun et al. (1994)] and finally was made available by [Baddeley et al. (2020)] as part of the package “spatstat.data” in the programming environment R. More precisely, in a comparative perspective, the four distributions on the unit simplex considered in the present paper were fitted to the data set obtained as follows. The above region was rescaled to the unit square and cut by the diagonal joining the vertices $(1,0)$ and $(0,1)$ into two congruent triangles; hence, the aforementioned analysis was carried out on the subset of $n=346$ locations lying in the interior of the upper triangle after undergoing the transformation $(x_{1},x_{2})\mapsto(x^{\prime}_{1},x^{\prime}_{2})=(1-x_{2},1-x_{1})$ that mapped this latter into the unit simplex in $\mathbb{R}^{2}$ (see Figure 3).

What is of interest to us is to check if the evidence suggests that the data portions having values next to the vertices of $\mathcal{S}^{2}$ are large enough to let the densities of such models take on positive and finite limits as $(x^{\prime}_{1},x^{\prime}_{2})$ tends to $(1,0)$ , $(0,1)$ and $(0,0)$ . Hence, the problem at study can be formalized by testing the hypotheses that all three shape parameters of each model or two of them, at the very least, are unitary.

That said, the method of Maximum-Likelihood (ML) was applied in order to derive the estimates for the parameter vectors of all the above models. In this regard, let $\left\{\underline{X}_{i}=\left(X_{i1},X_{i2}\right)\right\}_{i=1}^{n}$ be a sequence of independent random vectors with identical distribution depending on an unknown parameter vector $\underline{\theta}$ . First, suppose that $\left(X_{i1},X_{i2}\right)\stackrel{{\scriptstyle\bot}}{{\sim}}\mbox{{Dir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3}\right)$ , $i=1,\,\ldots\,,n$ ; by Eq. (15), given the observed sample $\left\{\underline{x}_{i}=\left(x_{i1},x_{i2}\right)\right\}_{i=1}^{n}$ , the log-likelihood function for the parameter vector $\left(\alpha_{1},\alpha_{2},\alpha_{3}\right)$ of the $\mbox{{Dir}}^{\,2}$ model is

$\displaystyle l_{\,\mbox{{Dir}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3};\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)=$			(132)
	$\displaystyle=$	$\displaystyle n\left[\,\log\Gamma\left(\alpha^{+}\right)-\log\Gamma\left(\alpha_{1}\right)-\log\Gamma\left(\alpha_{2}\right)-\log\Gamma\left(\alpha_{3}\right)\,\right]+\left(\alpha_{1}-1\right)\sum_{i=1}^{n}\log x_{i1}+$
	$\displaystyle+$	$\displaystyle\left(\alpha_{2}-1\right)\sum_{i=1}^{n}\log x_{i2}+\left(\alpha_{3}-1\right)\sum_{i=1}^{n}\log\left(1-x_{i1}-x_{i2}\right)\,.$

By taking $\left(X_{i1},X_{i2}\right)\stackrel{{\scriptstyle\bot}}{{\sim}}\mbox{{KB}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)$ , $i=1,\,\ldots\,,n$ , in light of Eqs. (36) and (132), given the observed sample $\left\{\underline{x}_{i}=\left(x_{i1},x_{i2}\right)\right\}_{i=1}^{n}$ , the log-likelihood function for the parameter vector $\left(\alpha_{1},\alpha_{2},\alpha_{3},\delta\right)$ of the $\mbox{{KB}}^{\,2}$ model is

	$\displaystyle l_{\,\mbox{{KB}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3},\delta;\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)\quad=$
		$\displaystyle=$	$\displaystyle l_{\,\mbox{{Dir}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3};\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)-\left[\sum_{i=1}^{n}\left(x_{i1}+x_{i2}\right)\right]\delta-n\log\,_{1}F_{1}\left(\alpha_{1}+\alpha_{2};\alpha^{+};-\delta\right)\,.$

Now let $\left(X_{i1},X_{i2}\right)\stackrel{{\scriptstyle\bot}}{{\sim}}\mbox{{NcDir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ , $i=1,\,\ldots\,,n$ ; by Eqs. (56) and (132), given the observed sample $\left\{\underline{x}_{i}=\left(x_{i1},x_{i2}\right)\right\}_{i=1}^{n}$ , the log-likelihood function for the parameter vector $\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ of the $\mbox{{NcDir}}^{\,2}$ model is

	$\displaystyle l_{\,\mbox{{NcDir}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3};\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)=l_{\,\mbox{{Dir}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3};\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)+$
		$\displaystyle-$	$\displaystyle\frac{n\,\lambda^{+}}{2}+\sum_{i=1}^{n}\log\Psi_{2}^{(3)}\left[\alpha^{+};\alpha_{1},\alpha_{2},\alpha_{3};\frac{\lambda_{1}}{2}x_{i1},\frac{\lambda_{2}}{2}x_{i2},\frac{\lambda_{3}}{2}\left(1-x_{i1}-x_{i2}\right)\right]\,,$

where, in light of Eqs. (3) and (57), by means of simple computations, the above specified function $\Psi_{2}^{\,(3)}$ can be expressed as

$\displaystyle\Psi_{2}^{\,(3)}\left[\alpha^{+};\alpha_{1},\alpha_{2},\alpha_{3};\frac{\lambda_{1}}{2}x_{1},\frac{\lambda_{2}}{2}x_{2},\frac{\lambda_{3}}{2}\left(1-x_{1}-x_{2}\right)\right]=$			(133)
	$\displaystyle=$	$\displaystyle\sum_{j_{1}=0}^{+\infty}\frac{\left(\alpha^{+}\right)_{j_{1}}}{\left(\alpha_{1}\right)_{j_{1}}}\frac{\left(\frac{\lambda_{1}}{2}x_{1}\right)^{j_{1}}}{j_{1}!}\,\Psi_{2}\left[\alpha^{+}+j_{1};\alpha_{2},\alpha_{3};\frac{\lambda_{2}}{2}x_{2},\frac{\lambda_{3}}{2}\left(1-x_{1}-x_{2}\right)\right]=$
	$\displaystyle=$	$\displaystyle\sum_{j_{1}=0}^{+\infty}\frac{\left(\alpha^{+}\right)_{j_{1}}}{\left(\alpha_{1}\right)_{j_{1}}}\frac{\left(\frac{\lambda_{1}}{2}x_{1}\right)^{j_{1}}}{j_{1}!}\sum_{j_{2}=0}^{+\infty}\frac{\left(\alpha^{+}+j_{1}\right)_{j_{2}}}{\left(\alpha_{2}\right)_{j_{2}}}\frac{\left(\frac{\lambda_{2}}{2}x_{2}\right)^{j_{2}}}{j_{2}!}\cdot$
	$\displaystyle\cdot$	${}_{1}F_{1}\left[\alpha^{+}+j_{1}+j_{2};\alpha_{3};\frac{\lambda_{3}}{2}\left(1-x_{1}-x_{2}\right)\right]\,,$

i.e. as a doubly infinite sum of weighted Kummer’s confluent hypergeometric functions. This formula can be usefully adopted as a natural basis for implementing such a function in any statistical package where the generalized hypergeometric function ${}_{1}F_{1}$ is already implemented. In this regard, the code of a possible implementation in R of the algorithm in Eq. (133) is proposed in Section 6. Similar reasonings clearly apply as far as the computation of the function $\Psi_{2}^{\,(m)}$ with $m>3$ is concerned.

Finally, let $\left(X_{i1},X_{i2}\right)\stackrel{{\scriptstyle\bot}}{{\sim}}\mbox{{CNcDir}}^{\,2}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ , $i=1,\,\ldots\,,n$ ; by Eqs. (85) and (132), given the observed sample $\left\{\underline{x}_{i}=\left(x_{i1},x_{i2}\right)\right\}_{i=1}^{n}$ , the log-likelihood function for the parameter vector $\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3}\right)$ of the $\mbox{{CNcDir}}^{\,2}$ model is

$\displaystyle l_{\,\mbox{{CNcDir}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3},\lambda_{1},\lambda_{2},\lambda_{3};\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)=$
	$\displaystyle+$	$\displaystyle l_{\,\mbox{{Dir}}^{\,2}}\left(\alpha_{1},\alpha_{2},\alpha_{3};\underline{x}_{1},\,\ldots\,,\underline{x}_{n}\right)+\sum_{i=1}^{n}\log\,_{0}F_{1}\left(\alpha_{1};\frac{\lambda_{1}}{4}x_{i1}\right)+$
	$\displaystyle+$	$\displaystyle\sum_{i=1}^{n}\log\,_{0}F_{1}\left(\alpha_{2};\frac{\lambda_{2}}{4}x_{i2}\right)+\sum_{i=1}^{n}\log\,_{0}F_{1}\left[\alpha_{3};\frac{\lambda_{3}}{4}\left(1-x_{i1}-x_{i2}\right)\right]+$
	$\displaystyle-$	$\displaystyle n\log\,_{0}F_{1}\left(\alpha^{+};\frac{\lambda^{+}}{4}\right)\,.$

Under well known regularity conditions, the ML estimator $\underline{\hat{\Theta}}$ for the parameter vector $\underline{\theta}$ has an asymptotic multivariate Normal distribution with mean vector $\underline{\theta}$ and covariance matrix which can be approximated by the inverse of the observed information matrix $I(\underline{\hat{\theta}};\,\underline{x})=\{-\partial^{\,2}l(\underline{\theta};\,\underline{x})\,/\,\partial\underline{\theta}\,\partial\underline{\theta}^{T}\}_{\underline{\theta}=\underline{\hat{\theta}}}$ , $\underline{\hat{\theta}}$ being the ML estimate for $\underline{\theta}$ . The maximization of the above log-likelihoods was numerically accomplished by using the built-in-function optim available in the statistical software R. The global maximum of the log-likelihoods is ensured to be actually achieved by using a wide range of starting values for the estimates of the parameters in the optimization algorithm.

The null hypotheses ( $H_{0}$ ) $\alpha_{1}=\alpha_{2}=\alpha_{3}=1$ , $\alpha_{1}=\alpha_{2}=1$ , $\alpha_{1}=\alpha_{3}=1$ and $\alpha_{2}=\alpha_{3}=1$ are issued towards all the models at study with the aim to verify if the tails of the data corresponding to all three vertices of the unit simplex or to a few couples of them are valuable enough to be captured by the densities of the aforementioned distributions. If all the above hypotheses are rejected, the model with all the parameters to be estimated is chosen; otherwise, in case of non-rejection of at least one of them, the model corresponding to the highest $p$ -value is selected amongst the non-rejected ones with fewest parameters. Such hypotheses can be easily verified by means of a suitable Likelihood Ratio (LR) test with asymptotic Chi-Squared distribution. The construction of such LR test requires the unconstrained maximization of the log-likelihoods as well as the constrained maximization of these latter under the hypotheses of interest. Specifically, the LR statistic is given by $w=-2(\hat{l}_{0}-\hat{l})$ , where $\hat{l}_{0}$ and $\hat{l}$ are, respectively, the maximum of the log-likelihood under $H_{0}$ and the maximum of the unrestricted log-likelihood. This statistic can be used to check if the fit of the models with two or three shape parameters set equal to 1 is statistically superior to the one of the model containing all the parameters. According to Wilks’ theorem [Wilks (1938)], on varying the sample, $w$ describes the random variable $W$ , the distribution of which has been approximated by a Chi-Squared with number of degrees of freedom equal to the number of parameters fixed by $H_{0}$ .

Details of the LR tests carried out herein are listed in Table 1, which indicates that the CNcDir and the NcDir models are the only for which the null hypothesis $\alpha_{1}=\alpha_{2}=\alpha_{3}=1$ is not rejected with the highest $p-$ value; therefore, their densities allow for the tails of the data corresponding to all three vertices of $\mathcal{S}^{2}$ to be captured. On the contrary, none of the four hypotheses at study is significantly supported by the data with reference to the remaining models; hence, all three shape parameters of these models are estimated. In this regard, the ML estimates for the parameters of the chosen models and the asymptotic standard errors of the corresponding estimators are reported in Table 2.

Table 1: Likelihood Ratio (LR) tests to evaluate the null hypotheses (

H_{0}

) that all three or two of the shape parameters

\alpha_{1},\,\alpha_{2},\,\alpha_{3}

are equal to 1 for the bivariate CNcDir, NcDir, Dir and KB models: maximum value of the log-likelihood under

H_{0}

(

\hat{l}_{0}

), maximum value of the unrestricted log-likelihood (

\hat{l}

), LR statistic (

w

), degrees of freedom (df) of the corresponding approximated

\chi^{2}

distribution and

p

-value. The

p

-values corresponding to the chosen models are written in bold (

n=346

Model	$H_{0}$	$\hat{l}_{0}$	$\hat{l}$	$w$	df	$p$ -value
$\mbox{{CNcDir}}^{\,2}$	$\alpha_{1}=\alpha_{2}=\alpha_{3}=1$	262.0291	262.7820	1.5058	3	0.6809
	$\alpha_{1}=\alpha_{2}=1$	262.3874		0.7892	2	0.6739
	$\alpha_{1}=\alpha_{3}=1$	262.0826		1.3988		0.4969
	$\alpha_{2}=\alpha_{3}=1$	262.1861		1.1918		0.5511
$\mbox{{NcDir}}^{\,2}$	$\alpha_{1}=\alpha_{2}=\alpha_{3}=1$	262.7941	263.8433	2.0984	3	0.5522
	$\alpha_{1}=\alpha_{2}=1$	263.0525		1.5816	2	0.4535
	$\alpha_{1}=\alpha_{3}=1$	262.9253		1.8360		0.3993
	$\alpha_{2}=\alpha_{3}=1$	263.1430		1.4006		0.4964
$\mbox{{Dir}}^{\,2}$	$\alpha_{1}=\alpha_{2}=\alpha_{3}=1$	239.8289	255.2603	30.8628	3	$<.0001$
	$\alpha_{1}=\alpha_{2}=1$	240.7546		29.0114	2
	$\alpha_{1}=\alpha_{3}=1$	243.9976		22.5254
	$\alpha_{2}=\alpha_{3}=1$	240.3968		29.7270
$\mbox{{KB}}^{\,2}$	$\alpha_{1}=\alpha_{2}=\alpha_{3}=1$	239.8813	255.3020	30.8414	3	$<.0001$
	$\alpha_{1}=\alpha_{2}=1$	245.2593		20.0854	2
	$\alpha_{1}=\alpha_{3}=1$	244.5826		21.4388
	$\alpha_{2}=\alpha_{3}=1$	240.4078		29.7884

Table 2: Maximum-Likelihood (ML) estimates for the parameters of the bivariate CNcDir, NcDir, Dir and KB models selected on the basis of the LR tests carried out in Table 1 and asymptotic standard errors (SE) of the corresponding estimators between round brackets (

n=346

Model		ML estimate (SE)	Model		ML estimate (SE)
$\mbox{{CNcDir}}^{\,2}$	$\alpha_{1}$	1	$\mbox{{NcDir}}^{\,2}$	$\alpha_{1}$	1
	$\alpha_{2}$	1		$\alpha_{2}$	1
	$\alpha_{3}$	1		$\alpha_{3}$	1
	$\hat{\lambda_{1}}$	42.7802 (8.5087)		$\hat{\lambda_{1}}$	3.0478 (0.4093)
	$\hat{\lambda_{2}}$	48.7569 (9.3385)		$\hat{\lambda_{2}}$	3.4644 (0.4382)
	$\hat{\lambda_{3}}$	44.1538 (8.6944)		$\hat{\lambda_{3}}$	3.1084 (0.4121)
$\mbox{{Dir}}^{\,2}$	$\hat{\alpha_{1}}$	1.2671 (0.0717)	$\mbox{{KB}}^{\,2}$	$\hat{\alpha_{1}}$	1.2810 (0.0870)
	$\hat{\alpha_{2}}$	1.3594 (0.0772)		$\hat{\alpha_{2}}$	1.3748 (0.0944)
	$\hat{\alpha_{3}}$	1.2818 (0.0726)		$\hat{\alpha_{3}}$	1.2543 (0.1187)
				$\hat{\delta}$	0.1884 (0.6519)

Moreover, the graphical depiction of the fitted densities of the four models under consideration is displayed in Figure 4, from which it is remarkable how much similar are the fitting of the two-dimensional CNcDir and NcDir densities on the one hand (top left and top right panels) and of the two-dimensional Dir and KB densities from the other (bottom left and bottom right panels). In particular, it is noticeable how much comparable are the abilities of the two non-central densities to identify the presence of the data portions next to the vertices of the unit simplex on one side and the behaviors of the two remaining densities on the other, which show all zero limits despite the four-parameter structure of the two-dimensional KB.

The main lesson to be learnt from the results reported herein is an essential equivalence between the CNcDir and the NcDir models. However, the fitting of the former is computationally less demanding than the one of the latter. Indeed, the amount of execution time that the CNcDir model enables to save with regard to the maximization of the log-likelihood is remarkable. This raises the question of further investigating the computational performances of the two models on the basis of efficiency arguments regarding the time to accomplish the Maximum-Likelihood fitting. In this respect, $n=30$ series of random values were generated from each of the two above models for selected combinations of the series size and of the parameter vector in the case $D=2$ . Then, the machine-time spent to fit the two models was measured for each replication and means and standard deviations were computed for this latter within each stratum. Hence, the null hypothesis of non-inferiority of the true mean of the CNcDir fitting time with respect to the one of the NcDir is checked by using the one-tailed $Z$ test for large samples. The results of the above described simulation study are listed in Table 3, from which we conclude that the Conditional model is beyond doubt to be preferred to the existing one for the major computational efficiency this permits in carrying out inferences based on the Maximum-Likelihood theory. More precisely, the obtained founds suggest that on average the NcDir fitting is approximately of the order of 90 times slower than the one of the CNcDir model.

Table 3: Means and standard deviations (between round brackets) of the machine-times (in seconds) to execute the Maximum-Likelihood fitting of the

\mbox{CNcDir}^{\,2}

and the

\mbox{NcDir}^{\,2}

densities to

n=30

series of

N

random draws and

p

-values of the one-tailed

Z

test to assess the null hypothesis

H_{0}:\,\mu_{CNcDir}-\mu_{NcDir}\,\geq\,0

for selected values of

N

and of the shape parameters and the non-centrality parameters vectors.

$\underline{\alpha}$	$\underline{\lambda}$		Time (^′′)
$\underline{\alpha}$	$\underline{\lambda}$		$N=25$	$N=50$	$N=100$
1.8	0.7	$\mbox{CNcDir}^{\,2}$	0.9053 (0.5852)	1.3360 (1.0834)	1.5253 (1.2979)
1.2	1.0	$\mbox{NcDir}^{\,2}$	54.5803 (41.8910)	72.9673 (76.5182)	107.7123 (126.4231)
1.5	0.9	$p$ -value	$<.0001$	$<.0001$	$<.0001$
0.7	4.6	$\mbox{CNcDir}^{\,2}$	1.0907 (0.8587)	1.7207 (1.0176)	1.3610 (0.9568)
1.3	0.5	$\mbox{NcDir}^{\,2}$	72.3550 (93.5070)	110.6507 (138.8951)	105.8680 (89.4877)
1.9	3.8	$p$ -value	$<.0001$	$<.0001$	$<.0001$
2.1	0.8	$\mbox{CNcDir}^{\,2}$	1.1537 (0.6813)	1.2613 (0.9024)	1.4517 (0.7845)
0.2	2.9	$\mbox{NcDir}^{\,2}$	91.0197 (88.6285)	142.1647 (101.2835)	172.0060 (98.3391)
0.6	4.2	$p$ -value	$<.0001$	$<.0001$	$<.0001$
0.8	2.4	$\mbox{CNcDir}^{\,2}$	0.7337 (0.3396)	1.0263 (0.5822)	1.2010 (0.6801)
0.9	1.7	$\mbox{NcDir}^{\,2}$	82.7857 (86.6345)	121.7610 (71.8105)	158.0427 (78.6970)
0.4	2.8	$p$ -value	$<.0001$	$<.0001$	$<.0001$

5 Conclusions

In this paper a new non-central extension of the Dirichlet model was obtained by conditioning the existing Non-central Dirichlet on the sum of the Non-central Chi-Squared random variables involved in its definition. The resulting Conditional Non-central Dirichlet distribution represents the multivariate generalization of the Conditional Doubly Non-central Beta one, namely the new non-central version of the Beta distribution which has been recently studied by [Orsi (2021)]. An in-depth analysis of this new multidimensional model was carried out herein in order to highlight similarities and differences with the existing one. Specifically, despite the densities of the new and the standard Non-central Dirichlet distributions share the same complex mixture type form, from the comparison between their perturbation representations the surprisingly greater tractability and interpretability of the former over the latter emerge. Hence, the Conditional Non-central Dirichlet distribution enables to overcome the strong analytical and mathematical limitations affecting the existing Non-central Dirichlet by exhibiting a substantially simpler density function than that of the latter. At the same time, the new non-central generalization of the Dirichlet model preserves the applicative potential of the existing one by allowing its density to take on arbitrary positive and finite limits which enables to properly capture the data portions next to the vertices of the support. These are the reasons why we hope this model may attract a wide range of applications in statistics.

6 Appendix

In the present section we shall propose a possible translation into R code of the algorithm specified in Eq. (133) to compute the generalization to three dimensions of the Humbert’s confluent hypergeometric function. In light of Eq. (56), this issue is instrumental in assessing the bivariate Non-central Dirichlet density. Specifically, the present implementation consists of two routines. The first routine, named psi2.3, implements the external infinite sum in Eq. (133). This routine internally calls the other one, named intpsi2.3, which computes the Humbert’s confluent hypergeometric function as infinite sum of Kummer’s confluent hypergeometric functions. The two routines share the same following arguments:

•

x1, x2: vectors of the first and of the second components of the data on $\mathcal{S}^{\,2}$ ;
•

shape1, shape2, shape3: the three shape parameters $\alpha_{1}$ , $\alpha_{2}$ , $\alpha_{3}$ of the bivariate Non-central Dirichlet density;
•

ncp1, ncp2, ncp3: the three non-centrality parameters $\lambda_{1}$ , $\lambda_{2}$ , $\lambda_{3}$ of the bivariate Non-central Dirichlet density;
•

tol: tolerance, where zero means keeping on the iterations until the additional terms do not change the partial sum of the infinite series under consideration;
•

maxiter: maximum number of iterations to perform;
•

debug: Boolean, where TRUE means returning debugging information whereas FALSE means returning only the final evaluation.

The code of the first routine is as follows:

psi2.3<-

function(x1,x2,shape1,shape2,shape3,ncp1,ncp2,ncp3,tol,maxiter,debug) {

L=c(shape1,shape2,shape3)

S=sum(L)

y1=(ncp1/2)*x1

coef<-1

temp<-intpsi2.3(x1=x1,x2=x2,shape1=shape1,shape2=shape2,shape3=shape3,

ncp1=ncp1,ncp2=ncp2,ncp3=ncp3,tol=tol,maxiter=maxiter,debug=debug)

out<-NULL

for(n in seq_len(maxiter)) {

coef<-coef*((S/L[1])*y1/n)

fac<-coef*

intpsi2.3(x1=x1,x2=x2,shape1=shape1+n,shape2=shape2,shape3=shape3,

ncp1=ncp1,ncp2=ncp2,ncp3=ncp3,tol=tol,maxiter=maxiter,debug=debug)

series<-temp+fac

if(debug) {

out<-c(out,fac)

}

if(hypergeo::isgood(series-temp,tol)) {

if(debug) {

return(list(series,out))

}

else {

return(series)

}

temp<-series

S<-S+1

L[1]<-L[1]+1

}

if(debug) {

return(list(series,out))

}

The code of the second routine is as follows:

intpsi2.3<-

function(x1,x2,shape1,shape2,shape3,ncp1,ncp2,ncp3,tol,maxiter,debug) {

L=c(shape1,shape2,shape3)

S=sum(L)

y2=(ncp2/2)*x2

y3=(ncp3/2)*(1-x1-x2)

coef<-1

temp<-hypergeo::genhypergeo(U=S,L=L[3],z=y3,tol=tol,maxiter=maxiter,

check_mod=TRUE,polynomial=FALSE,debug=FALSE)

out<-NULL

for(m in seq_len(maxiter)) {

coef<-coef*((S/L[2])*y2/m)

fac<-coef*hypergeo::genhypergeo(U=S+1,L=L[3],z=y3,tol=tol,

maxiter=maxiter,check_mod=TRUE,polynomial=FALSE,debug=FALSE)

series<-temp+fac

if(debug) {

out<-c(out,fac)

}

if(hypergeo::isgood(series-temp,tol)) {

if(debug) {

return(list(series,out))

}

else {

return(series)

}

temp<-series

S<-S+1

L[2]<-L[2]+1

}

if(debug) {

return(list(series,out))

}

References

[Baddeley et al. (2020)] Baddeley, A., Turner, R. and Rubak, E. 2020. spatstat.data: datasets for ’spatstat’. R package version 1.4-3, https://cran.r-project.org/web/packages/spatstat.data/index.html.
[Bran-Cardona et al. (2011)] Bran-Cardona, P. A., Orozco-Casta $\widetilde{\mbox{{n}}}$ eda, J. M. and Nagar, D. K. 2011. Bivariate Generalization of the Kummer-Beta Distribution. Revista Colombiana de Estadística 34:(3) 497-512.
[Gould (1972)] Gould, H. W. 1972. Combinatorial Identities. Morgantown (WV): Morgantown Printing and Binding Company.
[Johnson et al. (1995)] Johnson, N. L., Kotz, S. and Balakrishnan, N. 1995. Vol. 2 of Continuous Univariate Distributions, 2nd edn. New York: John Wiley & Sons.
[Johnson et al. (2005)] Johnson, N. L., Kemp, A. W. and Kotz, S. 2005. Univariate Discrete Distributions, 3rd edn. Hoboken (NJ): John Wiley & Sons.
[Kotz et al. (2000)] Kotz, S., Balakrishnan, N. and Johnson, N. L. 2000. Vol. 1 of Continuous Multivariate Distributions: Models and Applications, 2nd edn. New York: John Wiley & Sons.
[Law and Kelton (2000)] Law, A. M. and Kelton, W. D. 2000. Simulation Modeling and Analysis, 3rd edn. Boston: McGraw Hill.
[Lukacs (1955)] Lukacs, E. 1955. A Characterization of the Gamma Distribution. Annals of Mathematical Statistics 26:(2) 319-324.
[Ng et al. (2011)] Ng, K. W., Tian, G.-L., Tang, M.-L. 2011. Dirichlet and Related Distributions: Theory, Methods and Applications. Wiley Series in Probability and Statistics.
[Ongaro and Orsi (2015)] Ongaro, A. and Orsi, C. 2015. Some Results on Non-central Beta Distributions. Statistica 75:(1) 85-100.
[Orsi (2021)] Orsi, C. 2021. On the Conditional Non-central Beta Distribution. Accepted for publication at Statistica Neerlandica on 01 May 2021, doi:10.1111/stan.12249.
[Platt et al. (1988)] Platt, W. J., Evans, G. W. and Rathbun, S. L. 1988. The Population Dynamics of a Long-Lived Conifer (Pinus palustris). The American Naturalist, 131:(4) 491-525
[Rathbun et al. (1994)] Rathbun, S. L. and Cressie, N. 1994. A Space-Time Survival Point Process for a Longleaf Pine Forest in Southern Georgia. Journal of the American Statistical Association, 89:(428) 1164-1173
[Sánchez et al. (2006)] Sánchez, L. E., Nagar, D. K., Gupta, A. K. 2006. Properties of Noncentral Dirichlet Distributions. Computers and Mathematics with Applications, 52:(12) 1671-1682
[Siegel (1979)] Siegel, A. F. 1979. The Noncentral Chi-Squared Distribution with Zero Degrees of Freedom and Testing for Uniformity. Biometrika, 66:(2) 381-386
[Srivastava and Karlsson (1985)] Srivastava, H. M. and Karlsson, W. 1985. Multiple Gaussian Hypergeometric Series. Chichester (UK): Ellis Horwood.
[Wilks (1938)] Wilks, S. S. 1938. The Large-Sample Distribution of the Likelihood Ratio for Testing Composite Hypotheses. Annals of Mathematical Statistics, 9:(1) 60

	$\displaystyle f_{\left.\underline{X}^{\prime\prime}\,\right\|\,N^{+}}\left(\underline{x};\underline{\alpha},\underline{\lambda}\right)$	$\displaystyle=$	$\displaystyle\sum_{j^{+}\,\leq\,N^{+}}\left[\,\Pr\left(\left.\,\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}\right)\right.\cdot$
		$\displaystyle\cdot$	$\displaystyle\left.f_{\left.\underline{X}^{\prime\prime}\,\right\|\,\left(\,N_{1},\,\ldots\,,N_{D},N^{+}\right)}\left(\underline{x};\underline{\alpha},\underline{\lambda},\left(j_{1},\,\ldots\,,j_{D}\right)\right)\,\right]\,,$

$\displaystyle f_{\left.\left(\underline{X}^{\prime\prime},\,Z^{\prime+}\right)\,\right\|\,N^{+}}\left(\underline{x},z\right)$	$\displaystyle=$	$\displaystyle\sum_{j^{+}\leq N^{+}}\left[\,\Pr\left(\left.\,\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}\right)\cdot\right.$
	$\displaystyle\cdot$	$\displaystyle\left.f_{\left.\left(\underline{X}^{\prime\prime},\,Z^{\prime+}\right)\,\right\|\,\left(\,N_{1},\,\ldots\,,N_{D},N^{+}\,\right)}\left(\underline{x},z;\,\left(j_{1},\,\ldots\,,j_{D}\right)\,\right)\,\right]=$
	$\displaystyle=$	$\displaystyle\sum_{j^{+}\leq N^{+}}\left[\,\Pr\left(\left.\,\left(N_{1},\,\ldots\,,N_{D}\right)=\left(j_{1},\,\ldots\,,j_{D}\right)\,\right\|\,N^{+}\right)\cdot\right.$
	$\displaystyle\cdot$	$\displaystyle\left.f_{\left.\underline{X}^{\prime\prime}\,\right\|\,\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)}\left(\underline{x},\left(j_{1},\,\ldots\,,j_{D}\right)\,\right)\,\right]\cdot f_{\left.Z^{\prime+}\,\right\|\,\left(N_{1},\,\ldots\,,N_{D},N^{+}\right)}\left(z\right)$
	$\displaystyle=$	$\displaystyle f_{\left.\underline{X}^{\prime\prime}\,\right\|\,N^{+}}\left(\underline{x}\right)\cdot f_{\left.Z^{\prime+}\,\right\|\,N^{+}}\left(z\right)\,,$

A New Class of Non-Central Dirichlet Distributions

Abstract

1 Introduction

2 Preliminaries

Property 2.1 ([Lukacs (1955)] Characterizing property of independent Gamma random variables)

3 A new Non-central Dirichlet model

3.1 Definition and distribution

Definition 3.1 (DD-variate Conditional Non-central Dirichlet distribution)

Definition 3.2 ((D+1)(D+1)-variate Mixture Weight distribution)

Proposition 3.1 (Closure under permutation)

Proof.

Proposition 3.2 (Aggregation property)

Proof.

Proposition 3.3 (Univariate and bivariate marginals)

Proof.

Proposition 3.4 (Sum of the components)

Proof.

Proposition 3.5 (Closure under normalized conditioning)

Proof.

3.2 Mixing distribution

Proposition 3.6 (Closure under permutation)

Proof.

Proposition 3.7 (Marginals and conditionals)

Proof.

Proposition 3.8 (Sum of the components)

Proof.

Proposition 3.9 (Conditional distribution given the sum of the components)

Proof.

3.3 Representations

Proposition 3.10 (Mixture representation)

Proof.

Proposition 3.11 (Conditional density given N+N^{+})

Proof.

Proposition 3.12 (Unconditional composition type representation)

Proof.

Proposition 3.13 (Conditional independence)

Proof.

3.4 Density plots

Proposition 3.14 (Limiting values of the bivariate density)

Proof.

3.5 Moments

Proposition 3.15 (Mixed raw moments of the bivariate distribution)

Proof.

3.6 Identifiability

Proposition 3.16 (Identifiability)

Proof.

4 Applications

5 Conclusions

6 Appendix

References

Definition 3.1 ( $D$ -variate Conditional Non-central Dirichlet distribution)

Definition 3.2 ( $(D+1)$ -variate Mixture Weight distribution)

Proposition 3.11 (Conditional density given $N^{+}$ )