2. Bivariate Poisson generalized gamma process and its basic properties

In this section, the bivariate Poisson generalized gamma process is defined and its basic properties are derived. To define the new counting process model, we first introduce the generalized gamma distribution proposed in Agarwal and Kalla [Reference Agarwal and Kalla1] and Ghitany [Reference Ghitany15]. A random variable Φ is said to follow the generalized gamma distribution (GGD) with parameters $(\nu, k, \alpha, l)$, where $\nu\geq0, k, \alpha, l \gt 0$, if its probability density function (pdf) is given by

(2.1)\begin{align} f(\phi)=\frac{\alpha^{k-\nu}}{\Gamma_{\nu} (k, \alpha l)}\frac{\phi^{k-1}\exp\{-\alpha\phi\}}{(\phi+l)^{\nu}}, \phi \gt 0, \end{align}

where

\begin{equation*} \Gamma_{\nu} (k, \beta)=\int_{0}^{\infty} \frac{y^{k-1}\exp\{-y\}} {(y+\beta)^{\nu}}dy, \end{equation*}

for all β > 0, with

(2.2)\begin{align} \Gamma_{\nu} (k, \alpha l)=\int_{0}^{\infty} \frac{ y^{k-1}\exp\{-y\}}{(y+\alpha l)^{\nu}}dy=\int_{0}^{\infty} \frac{\alpha ^{k-\nu} y^{k-1}\exp\{-\alpha y\}}{(y+l)^{\nu}}dy. \end{align}

The function $\Gamma_{\nu} (k, \beta)$ is called the generalized gamma function (see Kobayashi [Reference Kobayashi20]) and if ν = 0, then

\begin{align*} \Gamma_{0} (k, \beta)=\int_{0}^{\infty} y^{k-1}\exp\{-y\}dy=\Gamma(k), ~\forall k \gt 0. \end{align*}

Thus, when ν = 0, it can be seen that the pdf in (2.1) becomes that of a gamma distribution with parameter $(k,\alpha)$. Hence, the GGD includes the gamma distribution as a special case.

Coming back to the general case, one can note from Gupta and Ong [Reference Gupta and Ong17] that

\begin{equation*} \Gamma_{\nu}(k,\beta)=\frac{\Gamma\left( k\right) }{\beta^{\nu-k}} \varphi\left( k,k-\nu+1;\beta\right), \end{equation*}

where

\begin{equation*} \varphi\left( a,c;x\right) =\frac{1}{\Gamma\left( a\right) }\int _{0}^{\infty}\frac{e^{-xt}t^{a-1}}{\left( 1+t\right) ^{a-c+1}}dt, \end{equation*}

is the confluent hypergeometric function of the second kind. This allows an easy computation of $\Gamma_{\nu}(k,\beta)$, as $\varphi\left( a,c;x\right) $ and $\Gamma\left( k\right) $ are implemented in most statistical or mathematical software. From Ghitany [Reference Ghitany15], the moment generating function of Φ is given by

(2.3)\begin{equation} M_\Phi(s)=E(e^{s\Phi})=\left(1-\frac s\alpha\right)^{\nu-k}\frac{\Gamma_\nu\left(k,(\alpha-s)l\right)}{\Gamma_\nu\left(k,\alpha l\right)},\:s \lt \alpha, \end{equation}

and its r-th moment about the origin is

(2.4)\begin{equation}E(\Phi^r)=\alpha^{-r}\frac{\Gamma_\nu\left(k+r,\alpha l\right)}{\Gamma_\nu\left(k,\alpha l\right)},\:r\in\mathbb{N}^\ast.\end{equation}

Now, we will define the bivariate Poisson generalized gamma process using the GGD. Let $\{{\mathbf{N}}(t),~t\geq0\}$, where ${\mathbf{N}}(t)=(N_{1} (t),N_{2}(t))$, be a bivariate process. Then, the corresponding “pooled” point process $\{M(t),~t\geq0\}$, where $M(t)=N_{1}(t) + N_{2}(t)$, can be defined. In this paper, we will consider regular (also known as orderly) multivariate point processes. In a univariate point process $\{N(t),~t\geq 0\}$, regularity is intuitively the nonoccurrence of multiple events in a small interval (see e.g., Cox and Lewis [Reference Cox and Lewis10], Finkelstein [Reference Finkelstein13, Reference Finkelstein14]; see also Cha and Giorgio [Reference Cha and Giorgio8]). Note that there are two types of regularity in multivariate point processes: (i) marginal regularity and (ii) regularity. For a multivariate point process, we say the process is marginally regular if its marginal processes, considered as univariate point processes, are all regular. The multivariate process is said to be regular if the pooled process is regular. Throughout this paper, we will assume that the multivariate process $\{{\mathbf{N}}(t),~t\geq0\}$ of our interest is a regular process. In the following, we shall use the notation $\Phi \sim{\mathcal{GG}} (\nu, k, \alpha, l)$ to represent that the continuous random variable Φ follows the GGD with parameters $(\nu, k, \alpha, l)$ and $\{N(t),~t\geq0\} \sim{\mathcal{NHPP}} (\eta(t)) $ to indicate that the counting process $\{N(t),~t\geq0\}$ follows the nonhomogeneous Poisson process (NHPP) with intensity function $\eta(t)$.

Throughout this paper, the BPGGP with the set of parameters $(\lambda_{1}(t),\lambda_{2}(t), \nu, k, \alpha, l)$ will be denoted by BPGGP$(\lambda_{1}(t),\lambda_{2}(t),\nu, k, \alpha, l)$.

Based on Definition 2.1, we now derive some basic properties of BPGGP and with that aim, let us introduce $\Lambda_{i}(t)\equiv\int_{0}^{t}\lambda_{i} (s)ds$,  $i=1,2$, $t \geq0$.

The proof is similar to that of Proposition 1 from Cha and Mercier [Reference Cha and Mercier9] in the univariate case and it is omitted. As in that paper, the second point of Proposition 2.2 shows that the BPGGP model as given in Definition 2.1 is not identifiable. An additional constraint should hence be added such as $l \equiv1$ for instance, wherever statistical procedures are studied (which is not the case in the present paper).

We now study the distributions for the numbers of events, which are of major interest for any counting process model.

Theorem 2.3. For $0\leq u_{i1} \lt u_{i2} \lt \cdots \lt u_{im}$, $i=1,2$,

\begin{align*} & P(N_{i}(u_{i2})-N_{i}(u_{i1}) =n_{i},~i=1,2)\\ & =\Bigg[\prod_{i=1}^{2}\frac{(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1} ))^{n_{i}}}{{n_{i}}!}\Bigg]\frac{\alpha^{k-\nu}}{(\alpha+\sum_{i=1} ^{2}(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1})))^{k+n_{1}+n_{2}-\nu}}\\ & ~~~~~~~~~~~~~~~~~~~~~~~~~~~\times\frac{\Gamma_{\nu}(k+n_{1}+n_{2} ,(\alpha+\sum_{i=1}^{2}(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1})))l)} {\Gamma_{\nu}(k,\alpha l)}, \end{align*}

and

\begin{align*} & P(N_{i}(u_{ij})-N_{i}(u_{ij-1}) =n_{ij},~i=1,2,j=1,2,\cdots,m)\\ & =\Bigg[\prod_{i=1}^{2}\prod_{j=1}^{m}\frac{(\Lambda_{i}(u_{ij})-\Lambda _{i}(u_{ij-1}))^{n_{ij}}}{{n_{ij}}!}\Bigg]\frac{\alpha^{k-\nu}}{(\alpha +\sum_{i=1}^{2}\sum_{j=1}^{m}(\Lambda_{i}(u_{ij})-\Lambda_{i}(u_{ij-1} )))^{k+\sum_{i=1}^{2}\sum_{j=1}^{m}n_{ij}-\nu}}\\ & ~~~~~~~~~~~~~~~~~~~~~~~~~~~\times\frac{\Gamma_{\nu}(k+\sum_{i=1}^{2} \sum_{j=1}^{m}n_{ij},(\alpha+\sum_{i=1}^{2}\sum_{j=1}^{m}(\Lambda_{i} (u_{ij})-\Lambda_{i}(u_{ij-1})))l)}{\Gamma_{\nu}(k,\alpha l)}. \end{align*}

Proof. From the definition of BPGGP$(\lambda _{1}(t),\lambda_{2}(t),\nu,k,\alpha,l)$,

\begin{align*} & P(N_{i}(u_{i2})-N_{i}(u_{i1})=n_{i},~i=1,2)\\ & =\int_{0}^{\infty}\Bigg[\prod_{i=1}^{2}\frac{(\phi(\Lambda_{i} (u_{i2})-\Lambda_{i}(u_{i1})))^{n_{i}}\exp\{-\phi(\Lambda_{i}(u_{i2} )-\Lambda_{i}(u_{i1}))\}}{n_{i}!}\Bigg]\frac{\alpha^{k-\nu}}{\Gamma_{\nu }(k,\alpha l)}\frac{\phi^{k-1}\exp\{-\alpha\phi\}}{(\phi+l)^{\nu}}d\phi\\ & =\Bigg[\prod_{i=1}^{2}\frac{(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1} ))^{n_{i}}}{{n_{i}}!}\Bigg]\frac{\alpha^{k-\nu}}{\Gamma_{\nu}(k,\alpha l)}\\ & ~~~~~~~~~~~~~~~~\times\int_{0}^{\infty}\frac{\phi^{k+n_{1}+n_{2}-1} \exp\{-\phi(\alpha+\sum_{i=1}^{2}(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1} )))\}}{(\phi+l)^{\nu}}d\phi\\ & =\Bigg[\prod_{i=1}^{2}\frac{(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1} ))^{n_{i}}}{{n_{i}}!}\Bigg]\frac{\alpha^{k-\nu}}{(\alpha+\sum_{i=1} ^{2}(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1})))^{k+n_{1}+n_{2}-\nu}}\\ & ~~~~~~~~~~~~~~~~~~~~~~~~~~~\times\frac{\Gamma_{\nu}(k+n_{1}+n_{2} ,(\alpha+\sum_{i=1}^{2}(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1})))l)} {\Gamma_{\nu}(k,\alpha l)}. \end{align*}

In a similar way,

\begin{align*} & P(N_{i}(u_{ij})-N_{i}(u_{ij-1}) =n_{ij},~i=1,2,j=1,2,\cdots,m)\\ & =\int_{0}^{\infty}\Bigg[\prod_{i=1}^{2}\prod_{j=1}^{m}\frac{(\phi (\Lambda_{i}(u_{ij})-\Lambda_{i}(u_{ij-1})))^{n_{ij}}\exp\{-\phi(\Lambda _{i}(u_{ij})-\Lambda_{i}(u_{ij-1}))\}}{n_{ij}!}\Bigg]\\ & ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\times \frac{\alpha^{k-\nu}}{\Gamma_{\nu}(k,\alpha l)}\frac{\phi^{k-1}\exp \{-\alpha\phi\}}{(\phi+l)^{\nu}}d\phi\\ & =\Bigg[\prod_{i=1}^{2}\prod_{j=1}^{m}\frac{(\Lambda_{i}(u_{ij})-\Lambda _{i}(u_{ij-1}))^{n_{ij}}}{{n_{ij}}!}\Bigg]\frac{\alpha^{k-\nu}}{(\alpha +\sum_{i=1}^{2}\sum_{j=1}^{m}(\Lambda_{i}(u_{ij})-\Lambda_{i}(u_{ij-1} )))^{k+\sum_{i=1}^{2}\sum_{j=1}^{m}n_{ij}-\nu}}\\ & ~~~~~~~~~~~~~~~~~~~~~~~~~~~\times\frac{\Gamma_{\nu}(k+\sum_{i=1}^{2} \sum_{j=1}^{m}n_{ij},(\alpha+\sum_{i=1}^{2}\sum_{j=1}^{m}(\Lambda_{i} (u_{ij})-\Lambda_{i}(u_{ij-1})))l)}{\Gamma_{\nu}(k,\alpha l)}. \end{align*}

The joint moments of $(N_{1}(t),N_{2}(t))$ also are of practical interest for the applications. They are obtained in the following theorem.

Remark 2.5. Based on (2.6) and (2.5), it is clear that $Cov\left( N_{1}(t_{1}),N_{2} (t_{2})\right) \gt 0$, which shows that $N_{1}(t_{1})$ and $N_{2}(t_{2})$ always are positively correlated, whatever the parameters of the BPGGP are, and whatever the times t ₁ and t ₂ are. Also, $\rho_{\left( N_{1} (t_{1}),N_{2}(t_{2})\right)} $ is non decreasing with respect to $\left( t_{1},t_{2}\right) $, with

\begin{align*} \lim_{\left( t_{1},t_{2}\right) \rightarrow\left( 0,0\right) ^{+}} \rho_{\left( N_{1}(t_{1}),N_{2}(t_{2})\right) } & =0,\\ \lim_{\left( t_{1},t_{2}\right) \rightarrow\left( \infty,\infty\right) }\rho_{\left( N_{1}(t_{1}),N_{2}(t_{2})\right) } & =1. \end{align*}

Now, in the following proposition, some important properties of the BPGGP will be stated. For this, we employ the same notations as those in Cha and Giorgio [Reference Cha and Giorgio8]. Denote by ${\mathcal{H}}_{Pt-} \equiv\{M(u),~0\leq u \lt t\}$ the history of the pooled process in $[0,t)$. Define $M(t-)$ as the total number of events in $[0,t)$ and T_i as the time from 0 until the arrival of the ith event in $[0,t)$ of the pooled process $\{M(t),~t\geq0\}$. Then ${\mathcal{H}} _{Pt-}$ can equivalently be defined in terms of $M(t-)$ and the sequential arrival points of the events $0\leq T_{1}\leq T_{2} \leq\cdots\leq T_{M(t-)} \lt t$ in $[0,t)$. Similarly, define the marginal histories of the marginal processes ${\mathcal{H}}_{it-} \equiv\{N_{i}(u),~0\leq u \lt t\}$, $i=1,2 $. Then, ${\mathcal{H}}_{it-} \equiv\{N_{i}(u),~0\leq u \lt t\}$ can also be defined in terms of $N_{i}(t-)$ and the sequential arrival points of the events $0\leq T_{i1}\leq T_{i2} \leq\cdots\leq T_{iN_{i}(t-)} \lt t$ in $[0,t)$, $i=1,2$, where $N_{i}(t-)$ is the total number of events of type i point process in $[0,t)$, $i=1,2$. In the following, we also use the definition of “p(t)-thinning” in Cha and Giorgio [Reference Cha and Giorgio8]. Also, we define univariate Poisson generalized Gamma process to define the marginal process.

See Cha and Mercier [Reference Cha and Mercier9] for various properties of the PGGP.

For convenience, we now introduce the following notations: $N_{ui}(t)\equiv N_{i}(u+t)-N_{i}(u)$, $\Lambda_{i}(t) \equiv\int_{0}^{t}\lambda_{i}(u)du$, $i=1,2$, $\lambda(t)\equiv\lambda_{1}(t)+\lambda_{2}(t)$, $\Lambda(t)=\int _{0}^{t}\lambda(u)du=\Lambda_{1}(t)+\Lambda_{2}(t)$, and $p_{i}(t)\equiv \lambda_{i}(t)/\lambda(t)$, $i=1,2$.

Properties (iii) and (iv) state about the conditional and unconditional restarting properties of the BPGGP. For more details on the restarting property, see Cha [Reference Cha7] and Cha and Giorgio [Reference Cha and Giorgio8].

3. Further properties of bivariate Poisson generalized gamma process

An efficient characterization for a multivariate point process can be done through the stochastic intensity approach (see Cox and Lewis [Reference Cox and Lewis10], Cha and Giorgio [Reference Cha and Giorgio8]). As mentioned in Cha and Giorgio [Reference Cha and Giorgio8], a marginally regular bivariate process can be specified by the following complete intensity functions:

(3.1)\begin{align} \lambda_{1t} & \equiv\lim_{\Delta t \rightarrow0} \frac{P\left( N_{1}(t,t+\Delta t)\geq1|{\mathcal{H}}_{1t-};{\mathcal{H} }_{2t-}\right) }{\Delta t}=\lim_{\Delta t \rightarrow0} \frac{P\left( N_{1}(t,t+\Delta t)=1|{\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-}\right) }{\Delta t}, \cr \lambda_{2t} & \equiv\lim_{\Delta t \rightarrow0} \frac{P\left( N_{2}(t,t+\Delta t)\geq1|{\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-}\right) }{\Delta t}=\lim_{\Delta t \rightarrow0} \frac{P\left( N_{2}(t,t+\Delta t)=1|{\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-}\right) }{\Delta t}, \cr \lambda_{12t} & \equiv\lim_{\Delta t \rightarrow0} \frac{P\left( N_{1}(t,t+\Delta t)N_{2}(t,t+\Delta t)\geq1|{\mathcal{H}}_{1t-};{\mathcal{H} }_{2t-}\right) }{\Delta t}, \end{align}

where $N_{i}(t_{1},t_{2})$, $t_{1} \lt t_{2}$, represents the number of events in $[t_{1} , t_{2})$, $i=1,2$, respectively (see Cox and Lewis [Reference Cox and Lewis10]). For a regular process, $\lambda_{12t}=0$, and it is sufficient to specify just $\lambda_{1t}$ and $\lambda_{2t}$ in (3.1) in order to define a regular process.

Theorem 3.1. The complete intensity functions of the BPGGP with the set of parameters $(\lambda_{1}(t),\lambda _{2}(t), \nu, k, \alpha, l)$ are given by

(3.2)\begin{align} \lambda_{it}=\frac{1}{(\alpha+\Lambda_{1}(t)+\Lambda_{2} (t))}\frac{\Gamma_{\nu} (k+N_{1}(t-)+N_{2}(t-)+1, (\alpha+\Lambda _{1}(t)+\Lambda_{2}(t))l)}{\Gamma_{\nu} (k+N_{1}(t-)+N_{2}(t-), (\alpha +\Lambda_{1}(t)+\Lambda_{2}(t))l)}\lambda_{i}(t),~i=1,2. \end{align}

Proof. Observe that

\begin{align*} \lambda_{1t} & =\lim_{\Delta t \rightarrow0} \frac{P\left( N_{1}(t,t+\Delta t)=1|{\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-}\right) }{\Delta t} \\ & =E_{(\Phi| {\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-})} \left[ \lim_{\Delta t \rightarrow0} \frac{P\left( N_{1}(t,t+\Delta t)=1|\Phi; {\mathcal{H}}_{1t-};{\mathcal{H} }_{2t-}\right) }{\Delta t}\right], \end{align*}

where $E_{(\Phi| {\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-})}[~\cdot~]$ stands for the expectation with respect to the conditional distribution of $(\Phi| {\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-})$ and

\begin{align*} \lim_{\Delta t \rightarrow0} \frac{P\left( N_{1}(t,t+\Delta t)=1|\Phi; {\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-}\right) }{\Delta t}=\Phi\lambda _{1}(t). \end{align*}

Thus, $\lambda_{1t}=E_{(\Phi| {\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-})} [\Phi\lambda_{1}(t)]$.

Similar to the procedure described in Cha [Reference Cha7], the conditional distribution of $(\Phi| {\mathcal{H}}_{it-}={\mathbf{h}}_{it-},~i=1,2)$, where ${\mathbf{h} }_{it-}\equiv(t_{i1},t_{i2},\cdots,t_{in_{i}},n_{i})$ is the realization of ${\mathcal{H}}_{it-}$, $i=1,2$, respectively, is given by

\begin{align*} & \Bigg(\phi^{n_{1}+n_{2}} \exp\left\{ -\phi\int_{0}^{t} \sum_{i=1} ^{2}\lambda_{i}(x)dx\right\} f(\phi) \Bigg) \\ & \times \Bigg(\int_{0}^{\infty}\phi^{n_{1}+n_{2}} \exp\left\{ -\phi\int_{0}^{t} \sum_{i=1}^{2}\lambda_{i}(x)dx\right\} f(\phi) d\phi\Bigg)^{-1}, \end{align*}

where we recall that f stands for the probability density function of Φ. Then,

\begin{align*} \lambda_{1t} & =E_{(\Phi| {\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-})} [\Phi\lambda_{1}(t)]\\ & = \frac{\int_{0}^{\infty} \phi^{n_{1}+n_{2}+1} \exp\left\{ -\phi\int_{0}^{t} \sum_{i=1}^{2}\lambda_{i}(x)dx\right\} f(\phi) d\phi}{\int_{0}^{\infty}\phi^{n_{1}+n_{2}} \exp\left\{ -\phi\int_{0}^{t} \sum_{i=1}^{2}\lambda_{i}(x)dx\right\} f(\phi) d\phi} \lambda_{1}(t), \end{align*}

which extends Proposition 4.1 in Grandell [Reference Grandell16] to a bivariate and non homogeneous mixed Poisson process. Now, one can check that

\begin{align*} & \int_{0}^{\infty}\phi^{n} \exp\left\{ -\phi\int_{0}^{t} \sum_{i=1} ^{2}\lambda_{i}(x)dx\right\} f(\phi) d\phi\\ & =\frac{\alpha^{k-\nu}} {(\alpha+\Lambda_{1}(t)+\Lambda_{2}(t))^{k+n-\nu}} \frac{\Gamma_{\nu} (k+n, (\alpha+\Lambda_{1}(t)+\Lambda_{2}(t))l)}{\Gamma_{\nu} (k, \alpha l)}. \end{align*}

Therefore,

\begin{align*} \lambda_{1t}=\frac{1}{(\alpha+\Lambda_{1}(t)+\Lambda_{2}(t))}\frac{\Gamma _{\nu} (k+N_{1}(t-)+N_{2}(t-)+1, (\alpha+\Lambda_{1}(t)+\Lambda_{2} (t))l)}{\Gamma_{\nu} (k+N_{1}(t-)+N_{2}(t-), (\alpha+\Lambda_{1} (t)+\Lambda_{2}(t))l)}\lambda_{1}(t). \end{align*}

The intensity function $\lambda_{2t}$ can be obtained symmetrically.

The result of Proposition 3.2 is also clear from the definition of the BPGGP and the characterization of the BVGPP in Theorem 2 of Cha and Giorgio [Reference Cha and Giorgio8].

It can be shown that

\begin{align*} \eta(k)\equiv\frac{\Gamma_{\nu} (k+1, \alpha l)}{\Gamma_{\nu} (k, \alpha l)}, \end{align*}

is increasing in k > 0 for any $\nu\geq0, \alpha, l \gt 0$ (see the proof of Proposition 8 in Cha and Mercier [Reference Cha and Mercier9]). Thus, the complete intensity functions in (3.2) are increasing in $N_{1}(t-)+N_{2}(t-)$. This implies that the proneness to the future event occurrence in each marginal process is increasing with the number of events, previously occurred in the pooled process.

4. Multivariate Poisson generalized gamma process

In this section, we study the multivariate Poisson generalized gamma process (MPGGP), by extending the results obtained in the previous sections. As in Cha and Giorgio [Reference Cha and Giorgio8], let $\{{\mathbf{N}}(t),~t\geq0\}$, where ${\mathbf{N} }(t)=(N_{1}(t),N_{2}(t),\cdots,N_{m}(t))$, be a multivariate process and define the corresponding “pooled” point process $\{M(t),~t\geq0\}$, where $M(t)=N_{1}(t)+N_{2}(t)+\cdots+N_{m}(t)$. Also, define the marginal point processes $\{N_{i}(t),~t\geq0\}$, $i=1,2,\cdots,m$, and the corresponding marginal histories of the marginal processes: ${\mathcal{H}}_{it-}$, $i=1,2,\cdots,m$. The MPGGP can be defined by generalizing Definition 2.1.

To state the properties of the MPGPP, we define $\lambda (t)=\sum_{i=1}^{m}\lambda_{i}(t)$, $\Lambda(t)=\int_{0}^{t}\lambda(v)dv$, $p_{i}(t)=\lambda_{i}(t)/\lambda(t)$, and $N_{ui}(t)\equiv N_{i} (u+t)-N_{i}(u)$, $i=1,2,\cdots,m$.

In addition to the basic properties stated in Proposition 4.2, it can be shown that the complete intensity functions are given by:

\begin{align*} \lambda_{it} &\equiv\lim_{\Delta t \rightarrow0} \frac{P\left( N_{i}(t,t+\Delta t)=1|{\mathcal{H}}_{1t-};{\mathcal{H}}_{2t-};\cdots ;{\mathcal{H}}_{mt-}\right) }{\Delta t} \\ & =\frac{1}{(\alpha+\sum_{i=1} ^{m}\Lambda_{i}(t))}\frac{\Gamma_{\nu} (k+\sum_{i=1}^{m}N_{i}(t-)+1, (\alpha+\sum_{i=1}^{m}\Lambda_{i}(t))l)}{\Gamma_{\nu} (k+\sum_{i=1}^{m} N_{i}(t-), (\alpha+\sum_{i=1}^{m}\Lambda_{i}(t))l)}\lambda_{i}(t), i=1,2,\cdots,m. \end{align*}

Furthermore, it can also be shown that

\begin{align*} & P(N_{i}(u_{i2})-N_{i}(u_{i1})=n_{i},~i=1,2,\cdots,m)\\ &= \Bigg[\prod_{i=1}^{m} \frac{(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1}))^{n_{i}}}{{n_{i}}!} \Bigg]\frac{\alpha^{k-\nu} }{(\alpha+\sum_{i=1}^{m}(\Lambda_{i}(u_{i2} )-\Lambda_{i}(u_{i1})))^{k+\sum_{i=1}^{m}n_{i}-\nu}}\\ & \times\frac{\Gamma_{\nu} (k+\sum_{i=1}^{m}n_{i}, (\alpha+\sum_{i=1}^{m}(\Lambda_{i}(u_{i2})-\Lambda_{i}(u_{i1})))l)} {\Gamma_{\nu} (k, \alpha l)}. \end{align*}

In addition to the above results, other properties obtained in the previous sections could be extended to the multivariate case in similar ways.

5. Comparison and monotony results with respect to the multivariate likelihood ratio ordering

The results of the previous sections were dependent on the specific properties of the generalized gamma distribution. We here provide some other results of MPGGPs based on the notion of multivariate likelihood ratio ordering, which hold in the more general setting of a multivariate mixed Poisson process which we now define.

In order to write down the different properties bellow, we now recall the notion of multivariate likelihood ratio ordering. We refer to Karlin and Rinott [Reference Karlin and Rinott18] for more details.

Definition 5.3. Let X and Y be two random vectors on $\mathbb{R}^{n}$ with respective density $f_{\mathbf{X}}$ and $f_{\mathbf{Y}}$ with respect to a common product measure $\sigma\left( d\mathbf{x}\right) $. (For instance X and Y may be both absolutely continuous or both discrete). Then, X is said to be smaller than Y in the multivariate likelihood ratio ordering (written ${\mathbf{X}}\prec_{\mathbf{lr}}{\mathbf{Y}}$) as soon as

\begin{equation*} f_{\mathbf{X}}\left( \mathbf{x}\right) f_{\mathbf{Y}}\left( \mathbf{y} \right) \leq f_{\mathbf{X}}\left( \mathbf{x}\wedge\mathbf{y}\right) f_{\mathbf{Y}}\left( \mathbf{x}\vee\mathbf{y}\right), \end{equation*}

for all x, $\mathbf{y}\in\mathbb{R}^{n}$, where the minimum $\wedge$ and the maximum $\vee$ are taken componentwise.

In the univariate case, the likelihood ratio order is denoted by $\prec_{lr} $.

Definition 5.4. Let X be a random vector on $\mathbb{R}^{n}$ with density $f_{\mathbf{X}}$ with respect to a product measure $\sigma\left( d\mathbf{x}\right) $. Then, X is said to be Multivariate Totally Positive property of order 2 (MTP2) as soon as

\begin{equation*} f_{\mathbf{X}}\left( \mathbf{x}\right) f_{\mathbf{X}}\left( \mathbf{y} \right) \leq f_{\mathbf{X}}\left( \mathbf{x}\wedge\mathbf{y}\right) f_{\mathbf{X}}\left( \mathbf{x}\vee\mathbf{y}\right), \end{equation*}

for all x, $\mathbf{y}\in\mathbb{R}^{n}$, which is just equivalent to ${\mathbf{X}}\prec_{\mathbf{lr}}{\mathbf{X}}$.

We first provide a result that extends Proposition 4 in Cha and Mercier [Reference Cha and Mercier9] to the multivariate setting.

Proof. Let $0\leq t_{1} \lt t_{2}$. Let us show that ${\ \mathbf{N}}(t_{1})\prec_{\mathbf{lr}}{\mathbf{N}}(t_{2})$.

Let

\begin{align*} g_{j}\left( n_{1},\cdots,n_{m}|\phi\right) & \equiv P\left( N_{i} (t_{j})=n_{i},~i=1,2,\cdots,m|\Phi=\phi\right) \\ & =\prod_{i=1}^{m}P\left( N_{i}(t_{j})=n_{i}|\Phi=\phi\right), \end{align*}

for $j=1,2$ and $n_{i}\in\mathbb{N}$, $i=1,2,\cdots,m$, where $P\left( N_{i}(t_{j})=\cdot|\Phi=\phi\right) $ is the Poisson distribution with parameter $\Lambda_{i}\left( t_{j}\right) ~\phi$. As this distribution increases with respect to $\Lambda_{i}\left( t_{j}\right) $ in the likelihood ratio ordering and as $\Lambda_{i}\left( t_{1}\right) \leq \Lambda_{i}\left( t_{2}\right) $, we derive that $P\left( N_{i} (t_{1})=\cdot|\Phi=\phi\right) \prec_{lr}P\left( N_{i}(t_{2})=\cdot |\Phi=\phi\right) $, and next that $g_{1}\left( \cdots|\phi\right) \prec_{\mathbf{lr}}g_{2}\left( \cdots|\phi\right) $, as the multivariate likelihood ratio ordering is stable through conjunction (see Shaked and Shanthikumar [Reference Shaked and Shanthikumar24], Theorem 6.E.4(a) page 299).

Using that $\Phi\prec_{lr}\Phi$, we derive from Theorem 2.4 in Karlin and Rinott [Reference Karlin and Rinott18] that

\begin{equation*} \int_{\mathbb{R}_{+}}g_{1}\left( \cdots|\phi\right) f_{\Phi}\left( \phi\right) d\phi\prec_{\mathbf{lr}}\int_{\mathbb{R}_{+}}g_{2}\left( \cdots|\phi\right) f_{\Phi}\left( \phi\right) d\phi, \end{equation*}

which is just equivalent to ${\mathbf{N}}(t_{1})\prec_{\mathbf{lr} }{\mathbf{\ N}}(t_{2})$, and allows to conclude.

We next show that the marginal increments in a MMPP (taken at possibly different times for each margin) exhibit the MTP2 property.

Proof. The distribution of $\left( N_{i}(u_{i2} )-N_{i}(u_{i1}),~i=1,2,\cdots,m\right) $ can be written as

\begin{equation*} P\left( N_{i}(u_{i2})-N_{i}(u_{i1})=n_{i},~i=1,2,\cdots,m\right) =\int_{\mathbb{R}_{+}}\left( \prod_{i=1}^{m}g_{i}\left( n_{i}|\phi\right) \right) f_{\Phi}\left( \phi\right) d\phi, \end{equation*}

for all $n_{i}\in\mathbb{N}$, $i=1,2,\cdots,m$, where

\begin{equation*} g_{i}\left( \cdot|\phi\right) =P\left( N_{i}(u_{i2})-N_{i}(u_{i1} )=\cdot|\Phi=\phi\right), \end{equation*}

stands for the Poisson distribution with parameter $\left( \Lambda_{i}\left( u_{i2}\right) -\Lambda_{i}\left( u_{i1}\right) \right) \phi$. As this distribution increases in the likelihood ordering with ϕ, the result follows from Property 7.2.18 in Denuit et al. [Reference Denuit, Dhaene, Goovaerts and Kaas12].

We recall that the MTP2 property is a strong positive dependence property, which entails for instance conditional increasingness in sequence and positive association, see e.g., Belzunce et al. [Reference Belzunce, Martínez-Riquelme and Mulero4]. Hence, such properties are fulfilled by the marginal increments in a MPGGP.

In a similar way, the concept of positive upper orthant dependent multivariate process (PUODMP) was defined in Cha and Giorgio [Reference Cha and Giorgio8]. We refer to this paper for more details. The following positive dependence result now is a direct consequence of Proposition 5.6.

Corollary 5.7. A MMPP is a positive upper orthant dependent multivariate process:

\begin{equation*} P(N_{i}(t_{i2})-N_{i}(t_{i1}) \gt n_{i},i=1,2,\cdots,m)\geq\Pi_{i=1}^{m} P(N_{i}(t_{i2})-N_{i}(t_{i1}) \gt n_{i}), \end{equation*}

for all $t_{i2} \gt t_{i1}$ and n_i, $i=1,2,\cdots,m$

We now come to the extension of Proposition 5 in Cha and Mercier [Reference Cha and Mercier9] to the multivariate setting.

Proposition 5.8. Let $\{{\mathbf{N}}(t),~t\geq0\}$ and $\{\bar{\mathbf{N}}(t),~t\geq0\}$ be two MMPPs with respective sets of parameters $\left( \lambda_{i}(t),i=1,2,\cdots,m,\Phi\right) $ and $\left( \bar{\lambda}_{i}(t),i=1,2,\cdots,m,\bar{\Phi}\right) $. Assume that $\Phi\prec_{\mathbf{lr}}\bar{\Phi}$. Also, let $0\leq u_{i1}\leq u_{i2}$, $i=1,2,\cdots,m$, such that

(5.1)\begin{equation} \Lambda_{i}\left( u_{i2}\right) -\Lambda_{i}\left( u_{i1}\right) \leq \bar{\Lambda}_{i}\left( u_{i2}\right) -\bar{\Lambda}_{i}\left( u_{i1}\right), \end{equation}

for all $i=1,2,\cdots,m$.

Then, we have the following result:

(5.2)\begin{equation} \left( N_{i}(u_{i2})-N_{i}(u_{i1}),~i=1,2,\cdots,m\right) \prec_{ \mathbf{lr}}\left( \bar{N}_{i}(u_{i2})-\bar{N}_{i}(u_{i1}),~i=1,2,\cdots ,m\right), \end{equation}

for all $0\leq u_{i1}\leq u_{i2}$, $n_{i}\in\mathbb{N}$, $i=1,2,\cdots,m$, where lr refers to the multivariate likelihood ratio ordering.

Proof. We use some arguments from the proof of Theorem 3.8 in Belzunce et al. [Reference Belzunce, Mercader, Ruiz and Spizzichino5]. See also Theorem 2.7 in Khaledi and Shaked [Reference Khaledi and Shaked19] (which is written only for absolutely continuous random variables and hence cannot be applied here).

Let us first write

(5.3)\begin{align} g\left( n_{1},\cdots,n_{m}|\phi\right) & \equiv P\left( N_{i} (u_{i2})-N_{i}(u_{i1})=n_{i},~i=1,2,\cdots,m|\Phi=\phi\right) \nonumber\\ & =\prod_{i=1}^{m}P\left( N_{i}(u_{i2})-N_{i}(u_{i1})=n_{i}|\Phi =\phi\right) , \end{align}

where $P\left( N_{i}(u_{i2})-N_{i}(u_{i1})=\cdot|\Phi=\phi\right) $ is the Poisson distribution with parameter $\left( \Lambda_{i}\left( u_{i2}\right) -\Lambda_{i}\left( u_{i1}\right) \right) \phi$.

Using similar arguments as for the proof of Proposition 5.5, it is clear that

\begin{equation*} g\left( \cdots|\phi\right) \prec_{\mathbf{lr}}\bar{g}\left( \cdots |\phi\right) , \end{equation*}

where $\bar{g}$ is defined in a similar way as (5.3) for the process $\{\bar{\mathbf{N}}(t),~t\geq0\}$, because $\Lambda_{i}\left( u_{i2}\right) -\Lambda_{i}\left( u_{i1}\right) \leq\bar{\Lambda}_{i}\left( u_{i2}\right) -\bar{\Lambda}_{i}\left( u_{i1}\right) $.

As the Poisson distribution with parameter $\left( \Lambda_{i}\left( u_{i2}\right) -\Lambda_{i}\left( u_{i1}\right) \right) \phi$ also increases with respect to ϕ in the likelihood ratio ordering, we derive that $g\left( n_{1},\cdots,n_{m}|\phi\right) $ is MTP2 in $\left( n_{1},\cdots,n_{m},\phi\right) $.

Then, based on Theorem 2.4 in Karlin and Rinott [Reference Karlin and Rinott18], we get that

\begin{equation*} \int_{\mathbb{R}_{+}}g\left( \cdots|\phi\right) f_{\Phi}\left( \phi\right) d\phi\prec_{\mathbf{lr}}\int_{\mathbb{R}_{+}}\bar{g}\left( \cdots|\bar{\phi }\right) f_{\bar{\Phi}}\left( \bar{\phi}\right) d\bar{\phi}, \end{equation*}

which is just equivalent to (5.2) and achieves the proof.

As a by-product of the previous proposition, considering $u_{i1}=0$ and $u_{i2}=t$ for all $i=1,2,\cdots,m$, one can see that, if all parameters are fixed except from one, then ${\mathbf{N}}(t)$ increases in the likelihood ordering when k increases or $\Lambda\left( t\right) $ increases, and when α or ν decreases.

We next explore the conditions given in Proposition 5.8 to derive the comparison result on a simple example (BPGGP).

For a given t, we know from Theorem 2.3 that the joint pdf of $\left( N_{1}(t),N_{2}(t)\right) $ is given by

\begin{align*} g_{12}\left( x_{1},x_{2}\right) & =\mathbb{P}\left[ \left( N_{1}(t)=x_{1},N_{2}(t)=x_{2}\right) \right] \\ & =\frac{\lambda_{1}^{x_{1}}\lambda_{2}^{x_{2}}t^{x_{1}+x_{2}}}{{x}_{1} !x_{2}!}\frac{\alpha^{k-\nu}}{\left[ \alpha+t\left( \lambda_{1}+\lambda _{2}\right) \right] ^{k+x_{1}+x_{2}-\nu}}\frac{\Gamma_{\nu}\left[ k+x_{1}+x_{2},\alpha+t\left( \lambda_{1}+\lambda_{2}\right) \right] }{\Gamma_{\nu}(k,\alpha)}, \end{align*}

for all $x_{1},x_{2}\in\mathbb{N}$, with a similar expression for the joint pdf of $\left( \bar{N}_{1}(t),\bar{N}_{2}(t)\right) $ (denoted by$\ \bar {g}_{12}$).

We set

\begin{align*} G\left( x_{1},x_{2},y_{1},y_{2}\right) & =g_{12}\left( x_{1}\wedge y_{1},x_{2}\wedge y_{2}\right) \bar{g}_{12}\left( x_{1}\vee y_{1},x_{2}\vee y_{2}\right) -g_{12}\left( x_{1},x_{2}\right) \bar{g}_{12}\left( y_{1},y_{2}\right) ,\\ \bar{G}\left( x_{1},x_{2},y_{1},y_{2}\right) & =\bar{g}_{12}\left( x_{1}\wedge y_{1},x_{2}\wedge y_{2}\right) g_{12}\left( x_{1}\vee y_{1},x_{2}\vee y_{2}\right) -g_{12}\left( x_{1},x_{2}\right) \bar{g} _{12}\left( y_{1},y_{2}\right), \end{align*}

for all $x_{1},x_{2},y_{1},y_{2}\in\mathbb{N}$. Then $\left( N_{1} (t),N_{2}(t)\right) \prec_{\mathbf{lr}}\left[ \succ_{\mathbf{\ lr}}\right] \left( \bar{N}_{1}(t),\bar{N}_{2}(t)\right) $ if and only if $G\left[ \bar{G}\right] $ remains non negative on $\mathbb{R}_{+}^{4}$.

We take t = 5, $\lambda_{1}=1 \lt $ $\bar{\lambda}_{1}=2$, $\lambda_{2} =1.5 \lt \bar{\lambda}_{2}=3$, ν = 1, k = 1, α = 1, $\bar{\nu}=1.25$, $\bar{k}=0.5$, $\bar{\alpha}=0.75$.

Then Condition (5.1) on the $\Lambda_{i}$’s and $\bar{\Lambda}_{i}$’s given in Proposition 5.8 is true. However Φ and $\bar{\Phi}$ are not comparable with respect to lr ordering. Indeed, the quotient of their respective pdfs (f and $\bar{f}$) is not monotonic, as can be seen in Figure 1.

Figure 1. The quotient $f/\bar{f}$ of the respective pdfs of Φ and $\bar{\Phi}$.

The functions $G(x_{1},x_{2},y_{1},y_{2})$ and $\bar{G}(x_{1},x_{2},y_{1},y_{2})$ are next plotted in Figure 2 with respect to $\left( x_{2},y_{2}\right) $ for $\left( x_{1},y_{1}\right) =\left( 3,10\right) $ and $\left( x_{1} ,y_{1}\right) =\left( 8,2\right) $, respectively. As can be seen, G and $\bar{G}$ both change sign on $\mathbb{R}_{+}^{4}$ and consequently, $\left( N_{1}(t),N_{2}(t)\right) $ and $\left( \bar{N}_{1}(t),\bar{N}_{2}(t)\right) $ are not comparable with respect to the bivariate likelihood ratio ordering (for t = 5).

Figure 2. The functions $G(x_{1},x_{2},y_{1},y_{2})$ and $\bar{G}(x_{1},x_{2},y_{1},y_{2})$ with respect to $\left( x_{2},y_{2}\right) $ for $\left( x_{1},y_{1}\right) =\left( 3,10\right) $ and $\left( x_{1},y_{1}\right) =\left( 8,2\right) $, respectively.

Based on the previous example, Condition (5.1) on the $\Lambda_{i}$’s and $\bar{\Lambda}_{i}$’s is not sufficient to derive the comparison result in Proposition 5.8, and some additional comparison assumption between Φ and $\bar{\Phi}$ is required.

We finally come to the comparison between the points of two MMPPs with different parameters. To begin with, we consider the case where the two processes share the same λ_i’s, $i=1,\dots,m$.

Proposition 5.10. Let $\{{\mathbf{N}}(t),~t\geq0\}$ and $\{\bar{\mathbf{N}}(t),~t\geq0\}$ be two MMPPs which share the same $\left( \lambda_{i}(t),i=1,2,\cdots,m\right) $ with different mixture distributions Φ and $\bar{\Phi}$, respectively. Assume that $\Phi \prec_{\mathbf{lr}}\bar{\Phi}$. For $i=1,2,\cdots,m$ and $n\in\mathbb{N}^{\ast}$, let T_in (resp. $\bar{T}_{in}$) be the n-th point of $\left\{ N_{i}\left( t\right) ,t\geq0\right\} $ (resp. $\left\{ \bar{N}_{i}\left( t\right) ,t\geq0\right\} $).

Then

\begin{equation*} \left( \bar{T}_{in_{i}},~i=1,2,\cdots,m\right) \prec_{\mathbf{lr}}\left( T_{in_{i}},~i=1,2,\cdots,m\right), \end{equation*}

for all $n_{i}\in\mathbb{N}^{\ast}$ and all $i=1,2,\cdots,m$.

Proof. Our aim is to use Theorem 3.8 in Belzunce et al. [Reference Belzunce, Mercader, Ruiz and Spizzichino5]. Let $i\in\left\{ 1,2,\cdots,m\right\} $ be fixed. Let us first show that $\left[ -T_{in_{i}}|\Phi=\phi\right] $ increases with respect to ϕ in the likelihood ratio ordering.

For that, we will use Theorem 3.7 in Belzunce et al. [Reference Belzunce, Lillo, Ruiz and Shaked3]. Let $\phi_{1} \leq\phi_{2}$. Then, the ratio of the respective density functions of $\left[ T_{i1}|\Phi=\phi_{1}\right] $ and $\left[ T_{i1}|\Phi=\phi_{2}\right] $ (first points in the i-th marginal processes) is

\begin{equation*} \frac{\phi_{1}\lambda_{i}\left( t\right) e^{-\phi_{1}\Lambda_{i}\left( t\right) }}{\phi_{2}\lambda_{i}\left( t\right) e^{-\phi_{2}\Lambda _{i}\left( t\right) }}=\frac{\phi_{1}}{\phi_{2}}e^{\left( \phi_{2}-\phi _{1}\right) \Lambda_{i}\left( t\right) }, \end{equation*}

and it is increasing. This implies that $\left[ T_{i1}|\Phi=\phi_{2}\right] \prec_{lr}\left[ T_{i1}|\Phi=\phi_{1}\right] $.

Also, the ratio of the corresponding cumulative hazard rate functions is

\begin{equation*} \frac{\phi_{1}\Lambda_{i}\left( t\right) }{\phi_{2}\Lambda_{i}\left( t\right) }=\frac{\phi_{1}}{\phi_{2}}, \end{equation*}

which is constant and hence non decreasing. Based on Theorem 3.7 in Belzunce et al. [Reference Belzunce, Lillo, Ruiz and Shaked3], we derive that $\left[ T_{in_{i}}|\Phi=\phi_{2}\right] \prec_{lr}\left[ T_{in_{i}}|\Phi=\phi_{1}\right] $, or equivalently that $\left[ -T_{in_{i}}|\Phi=\phi_{1}\right] \prec_{lr}\left[ -T_{in_{i}} |\Phi=\phi_{2}\right] $.

Hence, $\left[ -T_{in_{i}}|\Phi=\phi\right] $ increases with respect to ϕ in the likelihood ratio ordering.

Also, based on our assumptions, we know that $\Phi\prec_{lr}\bar{\Phi}$.

The result now is a direct consequence of Theorem 3.8 in Belzunce et al. [Reference Belzunce, Mercader, Ruiz and Spizzichino5].

A natural question now is: Is it possible to compare the points in two MMPPs which share the same Φ with different $\lambda_{i}(t)$’s for $i=1,2,\cdots,m$? The answer to this question is explored in next example.

Example 5.11. Let $\{{\mathbf{N}}(t),~t\geq0\}$ and $\{\bar{\mathbf{N}}(t),~t\geq0\}$ be two MPGGPs (which are specific MMPPs), which share $\left( \nu,k,\alpha,1\right) =\left( 0,1,1,1\right) $, so that Φ is exponentially distributed with mean 1. Let $\lambda_{i}(t)=\lambda_{i}$ and $\bar{\lambda}_{i}(t)=\bar {\lambda}_{i}$, $i=1,2$, be the corresponding constant baseline intensity functions of the NHPPs, respectively. Then, it is easy to check that the joint density function of $\left( T_{11},T_{21}\right) $ is

\begin{align*} f_{12}\left( x_{1},x_{2}\right) & =\int_{\mathbb{R}_{+}}\lambda_{1} \phi~e^{-\lambda_{1}\phi x_{1}}~\lambda_{2}\phi~e^{-\lambda_{2}\phi x_{2} }~e^{-\phi}d\phi\\ & =\frac{2\lambda_{1}\lambda_{2}}{\left( \lambda_{1}x_{1}+\lambda_{2} x_{2}+1\right) ^{3}}, \end{align*}

with a similar expression for the joint density function $\bar{f}_{12}$ of $\left( \bar{T}_{11},\bar{T}_{21}\right) $. The point is to see whether $\left( T_{11},T_{21}\right) $ and $\left( \bar{T}_{11},\bar{T} _{21}\right) $ are comparable with respect to the bivariate likelihood ordering. Let

\begin{align*} H\left( x_{1},x_{2},y_{1},y_{2}\right) & =f_{12}\left( x_{1}\wedge y_{1},x_{2}\wedge y_{2}\right) \bar{f}_{12}\left( x_{1}\vee y_{1},x_{2}\vee y_{2}\right) -f_{12}\left( x_{1},x_{2}\right) \bar{f}_{12}\left( y_{1},y_{2}\right) ,\\ \bar{H}\left( x_{1},x_{2},y_{1},y_{2}\right) & =\bar{f}_{12}\left( x_{1}\wedge y_{1},x_{2}\wedge y_{2}\right) f_{12}\left( x_{1}\vee y_{1},x_{2}\vee y_{2}\right) -f_{12}\left( x_{1},x_{2}\right) \bar{f} _{12}\left( y_{1},y_{2}\right), \end{align*}

for all $(x_{1},x_{2},y_{1},y_{2})\in\mathbb{R}_{+}^{4}$.

Then $\left( T_{11},T_{21}\right) \prec_{\mathbf{lr}}\left[ \succ _{\mathbf{\ lr}}\right] \left( \bar{T}_{11},\bar{T}_{21}\right) $ if and only if $H\left[ \bar{H}\right] $ remains non negative on $\mathbb{R} _{+}^{4}$. The functions $H(x_{1},x_{2},y_{1},y_{2})$ and $\bar{H}(x_{1} ,x_{2},y_{1},y_{2})$ are plotted in Figure 3 with respect to $\left( x_{2},y_{2}\right) $ for $\left( x_{1},y_{1}\right) =\left( 2,0.01\right) $ and $\left( x_{1},y_{1}\right) =\left( 0.01,1\right) $, respectively, with $\lambda_{1}=\bar{\lambda}_{1}=\lambda_{2}=1 \lt \bar{\lambda}_{2}=6$. As can be seen, H and $\bar{H}$ both change sign on $\mathbb{R}_{+}^{4}$ and consequently, $\left( T_{11},T_{21}\right) $ and $\left( \bar{T}_{11},\bar{T}_{21}\right)$ are not comparable with respect to the bivariate likelihood ratio ordering.

Based on this simple example (with constant λ_i’s and $\bar{\lambda }_{i}$’s such that $\lambda_{1}=\bar{\lambda}_{1}=\lambda_{2} \lt \bar{\lambda }_{2}$), it seems that there is no hope to find conditions under which the points in two MMPPs with different $\lambda_{i} (t)$’s could be comparable with respect to the multivariate likelihood ratio ordering.

Figure 3. The functions $H(x_{1},x_{2},y_{1},y_{2})$ and $\bar{H}(x_{1},x_{2},y_{1},y_{2})$ with respect to $\left( x_{2},y_{2}\right) $ for $\left( x_{1},y_{1}\right) =\left( 2,0.01\right) $ and $\left( x_{1},y_{1}\right) =\left( 0.01,1\right) $, respectively.

However, it is possible to get comparison results with respect to the weaker usual stochastic ordering. We recall that given two random vectors X and Y on $\mathbb{R}^{n}$, then X is said to be smaller than Y in the usual stochastic ordering (written ${\ \mathbf{X}}\prec_{\mathbf{sto}}{\mathbf{Y}}$) as soon as

\begin{equation*} E\left[ \varphi\left( \mathbf{X}\right) \right] \leq E\left[ \varphi\left( \mathbf{Y}\right) \right], \end{equation*}

for all non-decreasing function $\varphi:\mathbb{R}^{n}\longrightarrow \mathbb{R}$ such that the expectations exist. In the univariate setting (written $X\prec_{sto}Y$), it is equivalent to $\bar{F}_{X}\left( t\right) \leq\bar{F}_{Y}\left( t\right), $ for all $t\geq0$.

The multivariate likelihood ratio ordering is known to imply the usual stochastic ordering. See Shaked and Shanthikumar [Reference Shaked and Shanthikumar24] for more details.

We now come to the comparison result.

Proposition 5.12. Let $\{{\mathbf{N}}(t),~t\geq0\}$ and $\{\bar{\mathbf{N}}(t),~t\geq0\}$ be two MMPPs with sets of parameters $\left( \lambda_{i}(t),i=1,2,\cdots,m,\Phi\right) $ and $\left( \bar{\lambda}_{i}(t),i=1,2,\cdots,m,\bar{\Phi}\right) $, respectively. Assume that $\Phi\prec_{\mathbf{lr}}\bar{\Phi}$ and $\Lambda_{i}\left( t\right) \leq\bar{\Lambda}_{i}\left( t\right), $ for all $t\geq0$ and all $i=1,2,\cdots,m$. Using the notations of Proposition 5.10, we have:

\begin{equation*} \left( \bar{T}_{in_{i}},~i=1,2,\cdots,m\right) \prec_{\mathbf{sto}}\left( T_{in_{i}},~i=1,2,\cdots,m\right), \end{equation*}

for all $n_{i}\in\mathbb{N}^{\ast}$ and all $i=1,2,\cdots,m$.

Proof. Our aim is to use Theorem 3.1 in Belzunce et al. [Reference Belzunce, Mercader, Ruiz and Spizzichino5]. We already know from the proof of Proposition 5.10 that $\left[ -T_{in_{i}} |\Phi=\phi\right] $ increases with respect to ϕ in the likelihood ratio ordering, and hence also in the usual stochastic ordering.

Also, based on $\Lambda_{i}\left( t\right) \leq\bar{\Lambda}_{i}\left( t\right), $ for all $t\geq0$ and all $i=1,2,\cdots,m$, it is easy to check that the conditional survival functions of $\bar{T}_{i1}$ given $\bar{\Phi }=\phi$ and $T_{i1}$ given $\Phi=\phi$ fulfill

\begin{equation*} \bar{F}_{\bar{T}_{i1}|\bar{\Phi}=\phi}\left( t\right) \equiv\mathbb{P} \left( \bar{T}_{i1} \gt t|\bar{\Phi}=\phi\right) =e^{-\phi\bar{\Lambda} _{i}\left( t\right) }\leq e^{-\phi\Lambda_{i}\left( t\right) }=\bar {F}_{T_{i1}|\bar{\Phi}=\phi}\left( t\right), \end{equation*}

for all $t\geq0$, which means that $\left[ \bar{T}_{i1}|\bar{\Phi} =\phi\right] \prec_{sto}\left[ T_{i1}|\Phi=\phi\right] $. We derive from Theorem 3.1 in Belzunce et al. [Reference Belzunce, Lillo, Ruiz and Shaked3] that $\left[ \left( \bar{T}_{i1},\bar {T}_{i2},\dots,\bar{T}_{in_{i}}\right) |\bar{\Phi}=\phi\right] \prec _{sto}\left[ \left( T_{i1},T_{i2},\dots,T_{in_{i}}\right) |\Phi =\phi\right] $ and next that $\left[ \bar{T}_{in_{i}}|\bar{\Phi} =\phi\right] \prec_{sto}\left[ T_{in_{i}}|\Phi=\phi\right] $ (as the usual stochastic ordering is stable through marginalization), or equivalently that $\left[ -T_{in_{i}}|\Phi=\phi\right] \prec_{sto}\left[ \bar{T}_{in_{i} }|\bar{\Phi}=\phi\right] $.

Finally, based on $\Phi\prec_{lr}\bar{\Phi}$, we can derive that

\begin{equation*} \left( -T_{in_{i}},~i=1,2,\cdots,m\right) \prec_{\mathbf{sto}} \left( -\bar{T}_{in_{i}},~i=1,2,\cdots,m\right), \end{equation*}

from Theorem 3.1 in Belzunce et al. [Reference Belzunce, Mercader, Ruiz and Spizzichino5], which allows to conclude.