Parameter Estimation Problem for Doubly Geometric Process with the Gamma Distribution and Some Applications

1. Introduction

In the statistical literature, a data set with occurrence times of successive events generally can be modeled by using a counting process (CP). To determine a suitable stochastic CP model, it has to be tested whether the data set has a trend or not. Pekalp and Aydogdu [33] compared the monotonic trend tests for some counting processes. If the successive interarrival times are independent and identically distributed (iid) (there is no trend), the data set may be modeled by a renewal process (RP). However, in real-life examples, the successive inter-arrival times contain a monotone trend because of the aging effect and the accumulated wear [8]. In this case, this trend can be modeled by a nonhomogeneous Poisson process (NHPP) [2,3,4,5,6,7,8,9,10,11]. The data set having a monotone trend can also be analyzed by a geometric process (GP). GP is one of the widely used and known models for the monotone trend. Lam [18,19] first introduced the GP process. Furthermore, the process is used as a model in many areas in the reliability context. For details, see [37]. Lam [23] presented the GP theory and its applications. The real datasets with a monotone trend are modeled by GP in [40].

Although the GP is known as the most commonly used model, it has some limitations that can cause difficulties. The two limitations can be given as follows: Using the GP for non-monotone interarrival times with distributions with varying shape parameters is not suitable. The other limitation is that GP only allows logarithmic growth or explosive growth [7]. The GP model causes difficulties in the applications. Therefore, it can be said that the GP could be unsuitable for the mentioned cases.

To overcome the above difficulties, some stochastic models were developed by Wu and Scarf [37], Wu [38], and Wu and Wang [39]. One of the important models is the DGP for such models.

Wu [38] compares the DGP with the GP and exhibits the advantages and the preferability of the DGP. The definitions of CP, RP, GP, and DGP are presented as follows.

Definition 1 (Counting Process). 

$N\left(t\right)=Sup\left(n{:S}_n\le t\right)$ is the number of events that occurred in the interval (0,t], then {N(t),t ≥ 0} is called a CP where,  $S_0=0,S_n={\sum }_{k=1}^n{X}_{k}n=\mathrm{1,2},\dots$ be the occurrence time of $thenth$ event,  $X_k=S_k-S_{k-1}$ inter-arrival time be the time of $the(k-1)th$ and $kth$ event;

Definition 2 (Renewal Process). 

If {Xk,k = 1,2, … } are (iid) random variables with cumulative distribution function (cdf)F,a CP {N(t),t ≥ 0} corresponds to an RP.

Definition 3 (Geometric Process). 

Let’s assume that  $\left\{N\right(t),t\ge 0\}$ is a CP and  $\{X_k,k=\mathrm{1,2},\dots \}$ is the interarrival time of a CP $\left\{N\right(t),t\ge 0\}$ . If there exists a positive constant value  $a$ defined as a ratio parameter such that $Y_k{=a}^{k-1}{X}_{k}k=\mathrm{1,2},$ …, the CP corresponds to a GPThe expected value and the variance of a GP are given as $E\left(X_k\right)={\displaystyle \frac{\mu }{a^{k-1}}},Var\left(X_k\right)={\displaystyle \frac{{\sigma }^2}{a^{2\left(k-1\right)}}},k=\mathrm{1,2},\dots$ It can be easily seen that the parameters $a,\mu ,and{\sigma }^2$ uniquely determine the expected value and variance of $X_k\text{'}s$ by the formulas. Therefore it is clear that; the cdf of $X_1$ uniquely determines the cdf of  $X_k\text{'}s$ that is;  $F_k\left(x\right)=F\left(a^{k-1}x\right)$ , $k=\mathrm{1,2},3,\dots$ . As a result, the important role of the parameter estimation problem of $a,\mu ,{\sigma }^2$  in GP is seen.

Monotonicity has an important role in the theory of stochastic processes. The monotonicity properties can be defined as follows for a GP. If $a<(>)1$, then $\{X_k,k=\mathrm{1,2},\dots \}$ is defined as stochastically increasing (decreasing), and if $a=1$ then the GP corresponds to the RP.

Definition 4 (Doubly Geometric Process). 

Let’s assume that  $\left\{N\left(t\right),t\ge 0\right\}$ , is a CP, $\left\{{andX}_k,k=\mathrm{1,2},3,\dots \right\}$ , is the interarrival time of a CP $\left\{N\left(t\right),t\ge 0\right\}$ . If there exists a positive constant $a$ , defined as a ratio parameter, the CP corresponds to a DGP such that $Y_k{=a}^{k-1}{X_k}^{h\left(k\right)}k=\mathrm{1,2},$ … whereF is the cdf of the $X_1$ , $h\left(k\right)$ is a positive function of $k$ with $k=\mathrm{1,2},\dots h\left(1\right)=1.$

The DGP is considered with different $h\left(k\right)\text{'}s$ which are ${h\left(k\right)(1+\mathit{log}\left(k\right))}^b{,h\left(k\right)=b}^{k-1},{h\left(k\right)=b}^{log\left(k\right)}$, and $h\left(k\right)=1+b\mathrm{l}\mathrm{o}\mathrm{g}\left(k\right)$ by Wu (33) for ten real data sets in Lam (2007). In the study, the preferable performance of the DGP for $h\left(k\right)={\left(1+log\left(k\right)\right)}^b$, where $b$ is a positive constant value, and $log\left(k\right)$ denotes the logarithm function value under the base 10. Therefore, Wu [38] determines $h\left(k\right)$ as $h\left(k\right)={\left(1+log\left(k\right)\right)}^b$.

The probability density function (pdf), expected value, and the variance of $X_k$,

$k=\mathrm{1,2},\dots$. are given as follows.

$$f_{X_k}\left(x\right)=f_{Y_k}\left(a^{k-1}{x}^{h\left(k\right)}\right)a^{k-1}{x}^{h\left(k\right)-1}h\left(k\right),$$

$$E\left(X_k\right)=a^{\left(1-k\right)h^{-1}\left(k\right)}{\lambda }_{1k},Var\left(X_k\right)=a^{2\left(1-k\right)h^{-1}\left(k\right)}\left({\lambda }_{2k}-{{\lambda }_{1k}}^2\right)$$

where;

$${\lambda }_{1k}={\int }_0^{\infty }{x}^{h^{-1}\left(k\right)}f\left(x\right)dx$$

$${\lambda }_{2k}={\int }_0^{\infty }{x}^{{2h}^{-1}\left(k\right)}f\left(x\right)dx$$

Wu [33] obtains the monotonicity properties of the DGP given as follows;

i) If $P\left(X_1>1\right)=1$ and $b<0$, $0<a<1$ or if $P\left(0<X_1<1\right)=1$ and $0<b<4.898226,0<a<1$, $\{X_k,k=\mathrm{1,2},\dots \}$ increases stochastically.

ii) If $P\left(0<X_1<1\right)=1$ and $b<0,a>1$ or if$P\left(X_1>1\right)=1$ and $0<b<4.898226,a>1$, $\{X_k,k=\mathrm{1,2},3\dots \}$ decreases stochastically.

iii) If $(\mathit{log}\left(y\right)-\mathrm{k}\mathrm{l}\mathrm{o}\mathrm{g}\left(a\right))+{\left(1+\mathit{log}\left(k\right)\right)}^{-b}(\left(k-1\right)\log \left(a\right)-\mathit{log}\left(y\right)){\left(1+\mathit{log}\left(k+1\right)\right)}^{-b}$ alters between positive and negative values, then the sequence $\{X_k,k=\mathrm{1,2},\dots \}$ corresponds to a non-monotonous set over $k\text{'}s$ where $y$ shows $X_k$'s possible values for $k=\mathrm{1,2},\dots$.

The parameter estimation problem naturally arises in DGP. The parameter estimation problem for a DGP contains the parameters $a,b,\mu ,$ and${\sigma }^2.$ These parameters determine the mean and the variance of the first inter-arrival time $X_1$. Therefore, the parameter estimation problem for DGP is a very important issue.

Lam and Chan [21], Chan et al. [8], Aydoğdu et al. [3], and Kara et al. [14] used the lognormal, gamma, Weibull, and inverse Gaussian distributions, respectively, for the $X_1$ interarrival time to estimate the parameters for a GP. Kara et al. [17] consider the parameter estimation problem for the gamma geometric process.

Pekalp and Aydogdu [31] considered the power series expansions for the probability distribution, mean value, and variance function of a geometric process with gamma interarrival times. Pekalp and Aydogdu [34] considered the parameter estimation problem for the mean value and variance functions in GP.

Aydogdu and Altındag [5] computed the mean value and variance functions in a geometric process. Altındag[1] evaluated the multiple process data in a geometric process with exponential failures.

Yılmaz et al [41], Yılmaz [42] used the Bayesian inference with Lindley distribution and generalized exponential distribution by a GP, respectively.

Pekalp and Aydogdu [27] studied the integral equation for the second moment function in a GP. Pekalp and Aydogdu [28], Pekalp et al. [29] obtained the asymptotic solution of the integral equation for the second moment function of a GP by discriminating the lifetime distributions of the ten real data sets used in [22].

Pekalp et.al [30,32,35] estimate the parameters of a DGP by using the ML method under the exponential, Weibull, and lognormal distribution assumptions, respectively, for the first inter-arrival time $X_1$. Eroglu Inan [12] considered the parameter estimation problem for a DGP under the assumption that the first interarrival time has an inverse Gaussian distribution.

Although there are some studies in GP for the gamma distribution, which is used and is an important model in reliability theory, to the best of our knowledge, no study has been done yet about DGP, which removes the lack of GP model definitions. Therefore, the parametric statistical inference problem for DGP under the gamma distribution assumption must be studied.

Additionally, the estimators of the parameters can be obtained by the Modified Moment Method (MM) proposed by Saada et al. [36] for GP. Lam [20] introduces a least squares (LSE) method for nonparametric inference in GP. Aydoğdu and Kara [4] and Kara et al [16] considered a nonparametric estimation approach for α-series processes and GP. Wu [38] presented both MLE and LSE methods for DGP. Jasim and Qazaz [13] proposed the NP inference methods for DGP.

In this study, the parameter estimation problem for a DGP is considered under the assumption that the first inter-arrival time $X_1$ has a gamma distribution. Furthermore, to show the significance of the obtained maximum likelihood estimators in the real data sets, $\widehat{X_k}$‘s are obtained by ML and MM methods. The MM method is used to estimate the values of the parameters.

For these predictions, the MSE and MPE criteria given as follows are used to compare the methods.

$${MSE}^{*}={\displaystyle \frac{1}{n}}{\sum }_{k=1}^n{\left(X_k-\widehat{X_k}\right)}^2$$

$$MPE=\underset{1\le k\le n}{\underbrace{max}}{\displaystyle \frac{\left\{\left|S_k-\widehat{S_k}\right|\right\}}{S_k}}$$

The rest of this paper is organized as follows: In Section 2, the gamma distribution and its probabilistic properties are given. In Section 3, the log-likelihood function, first derivatives of the function, and the asymptotic joint distributions of the estimators are obtained. In Section 4, an extensive simulation study is performed, and the simulated mean, bias, and MSE values are calculated, and the small sample performances of the estimators are evaluated for various values of parameters and a certain repetition number. In Section 5, for a DGP, two illustrative examples are presented by using the ML and MM methods with two real data sets called the coal mining disaster data and the main propulsion diesel engine failure data. Finally, in Section 6, the results are discussed.

2. Gamma Distribution

The gamma distribution is one of the most commonly used models in the reliability and life testing areas for asymmetric data, for details, [9,24] referred to in [15].

The probability density function (pdf) of the gamma distribution is given by:

$$f_X\left(x,\alpha ,\beta \right)={\displaystyle \frac{1}{\Gamma \left(\alpha \right){\beta }^{\alpha }}}{x}^{\alpha -1}{e}^{-{\displaystyle \frac{x}{\beta }}}x>0,\alpha >0,\beta >0$$

(1)

where

$\alpha :$ shape parameter

$\beta$: scale parameter

If $\alpha =1$, the gamma distribution corresponds to the exponential distribution.

Let’s $X$ have a gamma distribution with $\alpha and\beta$parameters. In this case, it is written that $X~gamma\left(\alpha ,\beta \right)$. Some characteristics of the gamma distribution, such as expected value, variance, moment generating function, skewness, and kurtosis of $X~gamma\left(\alpha ,\beta \right)$ are given as $E\left(X\right)=\alpha \beta ,Var\left(X\right)=\alpha {\beta }^2,{\mu }_X\left(t\right)={(1-\beta t)}^{-\alpha },k=1,\dots n,{\gamma }_1={\displaystyle \frac{2}{\sqrt{\alpha }}},{\gamma }_2={\displaystyle \frac{6}{\alpha }}$ respectively

3. Statistical Inference for Gamma Distribution

Let’s assume that {N(t),t ≥ 0} is a DGP with the ratio parameter $a$, $h\left(k\right)={\left(1+logk\right)}^b,$ and the first inter-arrival time distribution $X_1$ has a $gamma\left(\alpha ,\beta \right)$ distribution. $\{X_1,X_2,\dots ,X_n\}$ is a realization of the DGP.

The random variables $Y_k=a^{k-1}{X}_k^{h\left(k\right)}$ are iid with the $gamma\left(\alpha ,\beta \right)$ distribution.

Then the PDF of $X_k$ is obtained as

$$f_{X_k}\left(x\right)=f\left(a^{k-1}{X}_k^{h\left(k\right)}\right)x_k^{h\left(k\right)-1}h\left(k\right)a^{k-1}$$

(2)

Hence, the likelihood function of the realization $\left\{X_1,\dots ,X_n\right\}$ is that

$$L\left(\left(a,b,\alpha ,\beta \right)\right)={\left({\displaystyle \frac{1}{\Gamma \left(\alpha \right){\beta }^{\alpha }}}\right)}^n{a}^{{\displaystyle \frac{\alpha n\left(n-1\right)}{2}}}{\prod }_{k=1}^{n}h\left(k\right)\prod _{k=1}^n{x_k}^{\alpha h\left(k\right)-1}{e}^{-{\sum }_{k=1}^n\left({\displaystyle \frac{a^{k-1}{x_k}^{h\left(k\right)}}{\beta }}\right)}$$

(3)

The log-likelihood function is clear that

$$lnL\left(a,b,\alpha ,\beta \right)=-n\left(\mathit{log}gamma\left(\alpha \right)+\alpha \mathit{log}\left(\beta \right)\right)+{\displaystyle \frac{\alpha n\left(n-1\right)}{2}}lna+{\sum }_{k=1}^{n}ln\left(h\left(k\right)\right)$$

$$+{\sum }_{k=1}^n\left(\alpha h\left(k\right)-1\right)ln\left(x_k\right)-{\sum }_{k=1}^n\left({\displaystyle \frac{a^{k-1}{x_k}^{h\left(k\right)}}{\beta }}\right)$$

(4)

Firstly, the first derivatives of the log-likelihood function defined as Equation (4) are taken concerning parameters $a,b,\alpha ,\beta$, respectively, and then the likelihood equations are obtained as follows by equating the derivatives to zero.

$${\displaystyle \frac{\partial lnL}{\partial a}}=\alpha {\displaystyle \frac{n\left(n-1\right)}{2a}}-{\displaystyle \frac{1}{\beta }}\sum _{k=1}^n\left(k-1\right)a^{k-2}{x_k}^{h\left(k\right)}=0$$

$${\displaystyle \frac{\partial lnL}{\partial b}}={\sum }_{k=1}^{n}ln\left(1+log\left(k\right)\right)+{\sum }_{k=1}^n\alpha ln\left(x_k\right){\left(1+log\left(k\right)\right)}^{b}ln\left(1+log\left(k\right)\right)$$

$$-{\displaystyle \frac{1}{\beta }}\sum _{k=1}^n{a}^{k-1}{x_k}^{h\left(k\right)}\mathit{ln}\left(x_k\right){\left(1+\mathit{log}\left(k\right)\right)}^{b}ln\left(1+\mathit{log}\left(k\right)\right)=0$$

$${\displaystyle \frac{\partial lnL}{\partial \alpha }}={\displaystyle \frac{n\left(n-1\right)}{2}}lna-n{\displaystyle \frac{\partial \Gamma \left(\alpha \right)}{\partial \alpha }}-nln\left(\beta \right)+\sum _{k=1}^{n}h\left(k\right)\mathit{ln}\left(x_k\right)=0$$

$${\displaystyle \frac{\partial lnL}{\partial \beta }}=-{\displaystyle \frac{n\alpha }{\beta }}+\sum _{k=1}^n{\displaystyle \frac{a^{k-1}{x_k}^{h\left(k\right)}}{{\beta }^2}}=0$$

In this study, firstly, it is seen that these likelihood equations can’t be solved analytically. Otherwise, it is seen that the values of the ML estimators of the $a,b,\alpha ,and\beta$ can be estimated numerically. For this purpose, the values are estimated numerically by using maximization of the log-likelihood function (4) with the Nmaximize subroutine, which is a constrained nonlinear optimization technique of the Mathematica package program in the simulation study. Then the MM method is used to estimate the values of the parameters.

In the second step, the Fisher information matrix and its inverse are obtained to investigate the asymptotic joint distribution, the asymptotic properties, and the diagonal elements, which give the asymptotic variance for each parameter, respectively.

The negative secondary derivatives of the $lnL$ function create the elements of the information matrix. Thus, the elements of the information matrix are shown by the notation given as follows.

$$I_{ij}=-E\left({\displaystyle \frac{{\partial }^{2}lnL\left(\theta \right)}{\partial {\theta }_i\partial {\theta }_j}}\right),i,j\in \left\{1,2,3,4\right\},\theta =({\theta }_1,{\theta }_2,{\theta }_3,{\theta }_4)$$

The diagonal elements of the Fisher information matrix are obtained as follows:

$I_{11}=-E\left({\displaystyle \frac{{\partial }^{2}lnL}{\partial a^2}}\right)={\displaystyle \frac{n\left(n-1\right)\alpha }{2a^2}}+{\displaystyle \frac{\alpha }{a^2}}{\sum }_{k=1}^n\left(k^2-3k+2\right)$

$I_{12}=-E\left({\displaystyle \frac{\partial lnL}{\partial a}}{\displaystyle \frac{\partial lnL}{\partial b}}\right)$

$={\sum }_{k=1}^n{\displaystyle \frac{\left(\mathit{ln}\left(1+log(k\right)\right)\left(k-1\right)}{a\beta }}\left[E\left(Y_k\mathit{ln}\left(Y_k\right)\right)-\left(k-1\right)lna\alpha \beta \right]$

$={\sum }_{k=1}^n{\displaystyle \frac{\left(\mathit{ln}\left(1+log(k\right)\right)\left(k-1\right)}{a\beta }}\left[6.4637-\left(k-1\right)lna\alpha \beta \right]$

$I_{13}=-E\left({\displaystyle \frac{\partial lnL}{\partial a}}{\displaystyle \frac{\partial lnL}{\partial \alpha }}\right)=-\left({\displaystyle \frac{n\left(n-1\right)}{2a}}\right)$

$I_{14}=-E\left({\displaystyle \frac{\partial lnL}{\partial a}}{\displaystyle \frac{\partial lnL}{\partial \beta }}\right)=-\left({\displaystyle \frac{\left(n^2-n\right)\alpha }{2a\beta }}\right)$

$I_{22}=-E\left({\displaystyle \frac{{\partial }^{2}lnL}{\partial b^2}}\right)=-\left[\alpha {\left(\mathit{ln}\left(1+log(k\right)\right)}^2\left[E\left(\mathit{ln}\left(Y_k\right)\right)-\left(k-1\right)lna\right]\right]$

$+\left[\left({\displaystyle \frac{1}{\beta }}\right){\sum }_{k=1}^n\left(\left(E\left(\left(Y_k\right)ln\left(Y_k\right)\right)-\left(k-1\right)lna\alpha \beta \right){\left(\mathit{ln}\left(1+log(k\right)\right)}^2\right)\right]$

$+\left[{\displaystyle \frac{{\left(\mathit{ln}\left(1+log(k\right)\right)}^2}{\beta }}\left({\sum }_{k=1}^{n}E\left(\left(Y_k\right){\left(ln\left(Y_k\right)\right)}^2\right)+\alpha \beta {\left(\left(k-1\right)lna\right)}^2\right.\right.$

$\left.-2E\left(\left(Y_k\right)ln\left(Y_k\right)\right)\left(k-1\right)lna\right]$ =

$-\left[\alpha {\left(\mathit{ln}\left(1+log(k\right)\right)}^2\left[1.1159-\left(k-1\right)lna\right]\right]$

$+\left({\displaystyle \frac{1}{\beta }}\right){\sum }_{k=1}^n\left(\left(6.4637-\left(k-1\right)lna\alpha \beta \right){\left(\mathit{ln}\left(1+log(k\right)\right)}^2\right)$

$+{\displaystyle \frac{{\left(\mathit{ln}\left(1+log(k\right)\right)}^2}{\beta }}{\sum }_{k=1}^n\left(12.0247\right)+\alpha \beta {\left(\left(k-1\right)lna\right)}^2$

$-2\left(6.4637\right)\left(k-1\right)lna$

$I_{23}=-E\left({\displaystyle \frac{\partial lnL}{\partial b}}{\displaystyle \frac{\partial lnL}{\partial \alpha }}\right)={\sum }_{k=1}^n\left(\mathit{ln}\left(1+log(k\right)\right)\left(E\left(\mathit{ln}\left(Y_k\right)\right)-\left(k-1\right)lna\right)$

$=-\left({\sum }_{k=1}^n\left(\mathit{ln}\left(1+log(k\right)\right)\left(1.1159-\left(k-1\right)lna\right)\right)$

$I_{24}=-E\left({\displaystyle \frac{\partial lnL}{\partial b}}{\displaystyle \frac{\partial lnL}{\partial \beta }}\right)$

$=-{\displaystyle \frac{1}{{\beta }^2}}{\sum }_{k=1}^n\left(E\left(\left(Y_k\right)ln\left(Y_k\right)\right)-\left(k-1\right)lna\alpha \beta \right)\left(\mathit{ln}\left(1+log(k\right)\right)$

$=-\left({\displaystyle \frac{1}{{\beta }^2}}{\sum }_{k=1}^n\left(\left(6.4637\right)-\left(k-1\right)lna\alpha \beta \right)\left(\mathit{ln}\left(1+log(k\right)\right)\right)$

$I_{33}=-E\left({\displaystyle \frac{{\partial }^{2}lnL}{\partial {\alpha }^2}}\right)=nPolyGamma\left[1,-1+\alpha \right]$

$I_{34}=E\left({\displaystyle \frac{{\partial }^2\mathit{ln}L}{\partial \alpha \partial \beta }}\right)={\displaystyle \frac{n}{\beta }}$

$I_{44}=-E\left({\displaystyle \frac{{\partial }^2\mathit{ln}L}{\partial {\beta }^2}}\right)={\displaystyle \frac{n\alpha }{{\beta }^2}}$

It is shown that the $ML$ estimators have an asymptotically normal $\left(AN\right)$ joint distribution for the $n$ sample size values that go to infinity by [6].

$$\widehat{\theta }~AN\left(\theta ,{V=I}^{-1}\right)$$

(5)

where $\theta$ mean vector and ${V=I}^{-1}$ variance-covariance matrix are given as

$$\theta =\left(\begin{array}{c}a\\ b\\ \alpha \\ \beta \end{array}\right),\widehat{\theta }=\left(\begin{array}{c}\widehat{a}\\ \widehat{b}\\ \widehat{\alpha }\\ \widehat{\beta }\end{array}\right),{V=I}^{-1}=\left[\begin{array}{cc}\begin{array}{cc}{v}_{11}& v_{12}\\ *& v_{22}\end{array}& \begin{array}{cc}{v}_{13}& v_{14}\\ v_{23}& v_{24}\end{array}\\ \begin{array}{cc}*& *\\ *& *\end{array}& \begin{array}{cc}{v}_{33}& v_{34}\\ *& v_{44}\end{array}\end{array}\right]$$

Thus, it is seen that the asymptotic distributions of each estimator can be given as follows,

$$\widehat{a}~AN\left(a,v_{11}\right)$$

$$\widehat{b}~AN\left(b,v_{22}\right)$$

$$\widehat{\mu }~AN\left(\mu ,v_{33}\right)$$

$$\widehat{\lambda }~AN\left(\lambda ,v_{44}\right)$$

where

${\mathit{v}}_{\mathit{i}\mathit{i}}:$ the diagonal element of matrix $V.i=\mathrm{1,2},\mathrm{3,4}$As a result, it can be said that the asymptotic variance values decrease as the sample size n goes to infinity. Therefore, it is seen that the estimators $\widehat{a},\widehat{b,}\widehat{\alpha },\widehat{and\beta }$ are consistent and asymptotically unbiased.

4. Simulation Method

In this section, a Monte Carlo simulation study is performed to compare the performance of the ML estimators for each parameter with various $a,b,\alpha ,and\beta$parameter values and sample sizes by using the Bias and MSE comparison criteria. Throughout the study, Mathematica software is used for calculations.

The simulation is evaluated for a specified number of repetitions, step by step:

1. Generate a sample$\{Y_1,Y_2,\dots ,Y_n\}$ from a gamma distribution with the parameters $\mu$and$\lambda .$

2. Calculate a realization$\left\{X_1,X_2,\dots ,X_n\right\}$ of the DGP with the parameters $a,b$ by using $X_k={\left({\displaystyle \frac{Y_k}{a^{k-1}}}\right)}^{1/{(1+\mathrm{l}\mathrm{o}\mathrm{g}k)}^b},k=\mathrm{1,2},3\dots ,n$.

3. Calculate the mean, bias, and mean square error (MSE) values for each parameter based on the $\left\{X_1,X_2,\dots ,X_n\right\}$ data set with a specified number of repetitions.

Lam [23] considers the application of GP with the 0.95, 0.99, 1.01, and 1.05 values of the parameter a. Therefore, these values are used throughout the study.

The monotonic property of the DGP is altered by the parameter $b$’s positive and negative values. The possible effect of positive and negative values of the $b$ is considered by taking the same constant -2, 2 in this study.

The parameter$\alpha$ and$\beta$ values are chosen as 1 and 2, respectively. The sample size values are chosen as 30, 50, $and10,\mathrm{w}\mathrm{h}\mathrm{i}\mathrm{c}\mathrm{h}\mathrm{a}\mathrm{r}\mathrm{e}$commonly used values in applications. The simulation study is repeated for 1000 trials. Table 1, Table 2, Table 3 and Table 4 include the corresponding simulated mean, biases, and MSE values for each estimator.

It is known that the ML estimators are asymptotically consistent, unbiased, and efficient estimators. The asymptotic unbiasedness, consistency, and efficiency statistical properties of the estimators are evaluated by Table 1.

The MSE gets small values, and the absolute bias values get closer to zero for increasing $n$ values. This expected case shows the fact that the estimators $a,b,\alpha ,and\beta$are asymptotically consistent.

Table 5-7-9-11 shows the efficiency statistical properties of the estimators. The diagonal elements of the $I^{-1}$give the lower bounds of the variances called minimum variance bounds (MVB) of the parameters$a,b,\alpha ,and\beta$. If the MVB values compare the corresponding simulated variance values, and it is seen that the values get close to each other for increasing $n$ values, it can be said that the estimators are highly efficient estimators.

Table 5-7-9-11 includes the simulated variance values and lower variance bounds for $a=\mathrm{0.95,0.99,1.01,1.05};b=2,-2,\alpha =\mathrm{1,2},$and $\beta =1,$2. Table 6-8-10-12 shows the differences between the simulated variance values and the MVB values. When the tables are investigated, it is seen that the values get close to each other, all of the difference values decrease for increasing $n$ values, and this is an expected case. This case supports the fact that the ML estimators are highly efficient. Furthermore, it is seen in the table that the sufficiency level of $n$ sample size is 30 for the desired level of convergence of the simulated variance values to MVB values.

Table 5. Comparing simulation results about the variances and MVB results. $\left(\mathit{\alpha }=1,\mathit{\beta }=2,\mathit{b}=2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.0004	0.0383	0.0746	0.3333	0.0005	0.0640	0.0000	0.1207
	50	0.0001	0.0157	0.0734	0.3323	0.0001	0.0304	0.0000	0.1151
	100	0.00003	0.0065	0.0180	0.1021	0.00005	0.0118	0.0000	0.0374
0.99	30	0.0002	0.0378	0.0786	0.3548	0.0002	0.0621	0.0000	0.1164
	50	0.00005	0.0166	0.0726	0.3410	0.00006	0.0291	0.0000	0.1123
	100	0.00000	0.0060	0.0180	0.1064	0.00000	0.0111	0.0000	0.0382
1.01	30	0.0001	0.0361	0.0756	0.3217	0.00017	0.0610	0.0000	0.1132
	50	0.00003	0.0170	0.0745	0.3309	0.00003	0.0284	0.0000	0.1125
	100	0.00000	0.0060	0.0152	0.0978	0.00000	0.0107	0.0000	0.0381
1.05	30	0.00016	0.0348	0.0832	0.3355	0.00013	0.0587	0.0000	0.1125
	50	0.00003	0.0170	0.0752	0.3357	0.00003	0.0284	0.0000	0.1166
	100	0.00001	0.0060	0.0172	0.1036	0.00001	0.0097	0.0000	0.0387

Table 5. Comparing simulation results about the variances and MVB results. $\left(\mathit{\alpha }=1,\mathit{\beta }=2,\mathit{b}=2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.0004	0.0383	0.0746	0.3333	0.0005	0.0640	0.0000	0.1207
	50	0.0001	0.0157	0.0734	0.3323	0.0001	0.0304	0.0000	0.1151
	100	0.00003	0.0065	0.0180	0.1021	0.00005	0.0118	0.0000	0.0374
0.99	30	0.0002	0.0378	0.0786	0.3548	0.0002	0.0621	0.0000	0.1164
	50	0.00005	0.0166	0.0726	0.3410	0.00006	0.0291	0.0000	0.1123
	100	0.00000	0.0060	0.0180	0.1064	0.00000	0.0111	0.0000	0.0382
1.01	30	0.0001	0.0361	0.0756	0.3217	0.00017	0.0610	0.0000	0.1132
	50	0.00003	0.0170	0.0745	0.3309	0.00003	0.0284	0.0000	0.1125
	100	0.00000	0.0060	0.0152	0.0978	0.00000	0.0107	0.0000	0.0381
1.05	30	0.00016	0.0348	0.0832	0.3355	0.00013	0.0587	0.0000	0.1125
	50	0.00003	0.0170	0.0752	0.3357	0.00003	0.0284	0.0000	0.1166
	100	0.00001	0.0060	0.0172	0.1036	0.00001	0.0097	0.0000	0.0387

Table 6. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0001	0.0257	0.0746	0.2172
	50	0	0.0147	0.0734	0.2126
	100	0.00002	0.0053	0.0180	0.0647
0.99	30	0	0.0243	0.0786	0.2384
	50	0.00001	0.0125	0.0726	0.2287
	100	0	0.0051	0.0180	0.0682
1.01	30	0.00007	0.0249	0.0756	0.2184
	50	0	0.0114	0.0745	0.2085
	100	0	0.0047	0.0152	0.0597
1.05	30	0.00003	0.0239	0.0832	0.2230
	50	0	0.0114	0.0752	0.2191
	100	0	0.0037	0.0172	0.0649

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 6. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0001	0.0257	0.0746	0.2172
	50	0	0.0147	0.0734	0.2126
	100	0.00002	0.0053	0.0180	0.0647
0.99	30	0	0.0243	0.0786	0.2384
	50	0.00001	0.0125	0.0726	0.2287
	100	0	0.0051	0.0180	0.0682
1.01	30	0.00007	0.0249	0.0756	0.2184
	50	0	0.0114	0.0745	0.2085
	100	0	0.0047	0.0152	0.0597
1.05	30	0.00003	0.0239	0.0832	0.2230
	50	0	0.0114	0.0752	0.2191
	100	0	0.0037	0.0172	0.0649

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 7. Comparing simulation results about the variances and MVB results $\left(\mathit{\alpha }=1,\mathit{\beta }=2,\mathit{b}=-2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.0020	0.0000	0.0671	0.3346	0.0043	0.0446	0.0000	0.1195
	50	0.0008	0.00000	0.0432	0.2184	0.0017	0.0163	0.0000	0.0725
	100	0.0002	0.00000	0.0172	0.0972	0.0003	0.0030	0.0000	0.0390
0.99	30	0.0034	0.0000	0.0770	0.3257	0.0040	0.0571	0.0000	0.1144
	50	0.0020	0.0000	0.0403	0.2071	0.0018	0.0259	0.0000	0.0743
	100	0.0011	0.0000	0.0185	0.0997	0.0005	0.0078	0.0000	0.0386
1.01	30	0.0047	0.0000	0.0639	0.3244	0.0033	0.0616	0.0000	0.1151
	50	0.0029	0.0000	0.0398	0.1988	0.0014	0.0300	0.0000	0.0722
	100	0.0022	0.0000	0.0174	0.1030	0.0003	0.0109	0.0000	0.0377
1.05	30	0.0083	0.0000	0.0699	0.3282	0.0014	0.0636	0.0000	0.1173
	50	0.0065	0.0000	0.0348	0.1866	0.0002	0.0309	0.0000	0.0741
	100	0.0055	0.0000	0.0178	0.1046	0.0000	0.0106	0.0000	0.0384

Table 7. Comparing simulation results about the variances and MVB results $\left(\mathit{\alpha }=1,\mathit{\beta }=2,\mathit{b}=-2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.0020	0.0000	0.0671	0.3346	0.0043	0.0446	0.0000	0.1195
	50	0.0008	0.00000	0.0432	0.2184	0.0017	0.0163	0.0000	0.0725
	100	0.0002	0.00000	0.0172	0.0972	0.0003	0.0030	0.0000	0.0390
0.99	30	0.0034	0.0000	0.0770	0.3257	0.0040	0.0571	0.0000	0.1144
	50	0.0020	0.0000	0.0403	0.2071	0.0018	0.0259	0.0000	0.0743
	100	0.0011	0.0000	0.0185	0.0997	0.0005	0.0078	0.0000	0.0386
1.01	30	0.0047	0.0000	0.0639	0.3244	0.0033	0.0616	0.0000	0.1151
	50	0.0029	0.0000	0.0398	0.1988	0.0014	0.0300	0.0000	0.0722
	100	0.0022	0.0000	0.0174	0.1030	0.0003	0.0109	0.0000	0.0377
1.05	30	0.0083	0.0000	0.0699	0.3282	0.0014	0.0636	0.0000	0.1173
	50	0.0065	0.0000	0.0348	0.1866	0.0002	0.0309	0.0000	0.0741
	100	0.0055	0.0000	0.0178	0.1046	0.0000	0.0106	0.0000	0.0384

Table 8. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0023	0.0446	0.0671	0.2151
	50	0.0009	0.0163	0.0432	0.1459
	100	0.0001	0.003	0.0172	0.0582
0.99	30	0.0006	0.0571	0.0770	0.2113
	50	-0.0002	0.0259	0.0403	0.1328
	100	-0.0006	0.0078	0.0185	0.0611
1.01	30	-0.0019	0.0616	0.0639	0.2093
	50	-0.0015	0.03	0.0398	0.1266
	100	-0.0014	0.0109	0.0174	0.0653
1.05	30	-0.0069	0.0636	0.0699	0.2109
	50	-0.0063	0.0309	0.0348	0.1125
	100	-0.0055	0.0106	0.0178	0.0662

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 8. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0023	0.0446	0.0671	0.2151
	50	0.0009	0.0163	0.0432	0.1459
	100	0.0001	0.003	0.0172	0.0582
0.99	30	0.0006	0.0571	0.0770	0.2113
	50	-0.0002	0.0259	0.0403	0.1328
	100	-0.0006	0.0078	0.0185	0.0611
1.01	30	-0.0019	0.0616	0.0639	0.2093
	50	-0.0015	0.03	0.0398	0.1266
	100	-0.0014	0.0109	0.0174	0.0653
1.05	30	-0.0069	0.0636	0.0699	0.2109
	50	-0.0063	0.0309	0.0348	0.1125
	100	-0.0055	0.0106	0.0178	0.0662

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 9. Comparing simulation results about the variances and MVB results $\left(\mathit{\alpha }=2,\mathit{\beta }=1,\mathit{b}=2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.0002	0.0359	0.3572	0.0803	0.0002	0.0320	0.0279	0.0197
	50	0.00009	0.0161	0.1798	0.0455	0.00008	0.0152	0.0170	0.0127
	100	0.00002	0.0057	0.0758	0.0221	0.00002	0.0059	0.0086	0.0067
0.99	30	0.0001	0.0359	0.3151	0.0678	0.0001	0.0310	0.0280	0.0196
	50	0.00003	0.0161	0.1805	0.0451	0.00003	0.0145	0.0170	0.0125
	100	0.0000	0.0057	0.0786	0.0221	0.0000	0.0055	0.0086	0.0066
1.01	30	0.00009	0.0334	0.3484	0.0745	0.00008	0.0305	0.0281	0.0203
	50	0.00001	0.0146	0.1476	0.0435	0.00001	0.0142	0.0172	0.0133
	100	0.00000	0.0053	0.0835	0.0236	0.00000	0.0053	0.0086	0.0068
1.05	30	0.00007	0.0309	0.3299	0.0723	0.00006	0.0294	0.0280	0.0199
	50	0.00002	0.0141	0.1613	0.0431	0.00002	0.0134	0.0171	0.0131
	100	0.00001	0.0051	0.0791	0.0219	0.0000	0.0048	0.0086	0.0067

Table 9. Comparing simulation results about the variances and MVB results $\left(\mathit{\alpha }=2,\mathit{\beta }=1,\mathit{b}=2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.0002	0.0359	0.3572	0.0803	0.0002	0.0320	0.0279	0.0197
	50	0.00009	0.0161	0.1798	0.0455	0.00008	0.0152	0.0170	0.0127
	100	0.00002	0.0057	0.0758	0.0221	0.00002	0.0059	0.0086	0.0067
0.99	30	0.0001	0.0359	0.3151	0.0678	0.0001	0.0310	0.0280	0.0196
	50	0.00003	0.0161	0.1805	0.0451	0.00003	0.0145	0.0170	0.0125
	100	0.0000	0.0057	0.0786	0.0221	0.0000	0.0055	0.0086	0.0066
1.01	30	0.00009	0.0334	0.3484	0.0745	0.00008	0.0305	0.0281	0.0203
	50	0.00001	0.0146	0.1476	0.0435	0.00001	0.0142	0.0172	0.0133
	100	0.00000	0.0053	0.0835	0.0236	0.00000	0.0053	0.0086	0.0068
1.05	30	0.00007	0.0309	0.3299	0.0723	0.00006	0.0294	0.0280	0.0199
	50	0.00002	0.0141	0.1613	0.0431	0.00002	0.0134	0.0171	0.0131
	100	0.00001	0.0051	0.0791	0.0219	0.0000	0.0048	0.0086	0.0067

Table 10. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0000	0.0039	0.3293	0.0606
	50	0.00001	0.0009	0.1628	0.0328
	100	0.0000	0.0002	0.0672	0.0154
0.99	30	0.0000	0.0049	0.2871	0.0482
	50	0.0000	0.0016	0.1635	0.0326
	100	0.0000	0.0002	0.07	0.0155
1.01	30	0.00001	0.0029	0.3203	0.0542
	50	0.0000	0.0004	0.1304	0.0302
	100	0.0000	0.0000	0.0749	0.0168
1.05	30	0.00001	0.0015	0.3019	0.0524
	50	0.0000	0.0007	0.1442	0.0300
	100	0.00001	0.0003	0.0705	0.0152

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 10. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0000	0.0039	0.3293	0.0606
	50	0.00001	0.0009	0.1628	0.0328
	100	0.0000	0.0002	0.0672	0.0154
0.99	30	0.0000	0.0049	0.2871	0.0482
	50	0.0000	0.0016	0.1635	0.0326
	100	0.0000	0.0002	0.07	0.0155
1.01	30	0.00001	0.0029	0.3203	0.0542
	50	0.0000	0.0004	0.1304	0.0302
	100	0.0000	0.0000	0.0749	0.0168
1.05	30	0.00001	0.0015	0.3019	0.0524
	50	0.0000	0.0007	0.1442	0.0300
	100	0.00001	0.0003	0.0705	0.0152

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 11. Comparing simulation results about the variances and MVB results. $\left(\mathit{\alpha }=2,\mathit{\beta }=1,\mathit{b}=-2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.00001	0.0051	0.3157	0.0783	0.0000	0.0048	0.0281	0.0201
	50	0.0004	0.0000	0.1481	0.0422	0.0009	0.0092	0.0172	0.0134
	100	0.0001	0.0000	0.0786	0.0224	0.0002	0.0017	0.0086	0.0066
0.99	30	0.0018	0.0000	0.3117	0.0686	0.0017	0.0303	0.0280	0.0195
	50	0.0010	0.0000	0.1762	0.0456	0.0008	0.0143	0.0171	0.0130
	100	0.00051	0.0000	0.0820	0.0233	0.0002	0.0046	0.0086	0.0068
1.01	30	0.0024	0.0000	0.3546	0.0715	0.0012	0.0318	0.0279	0.0191
	50	0.0016	0.0000	0.1675	0.0421	0.0005	0.0157	0.0171	0.0129
	100	0.0011	0.0000	0.0716	0.0217	0.0001	0.0060	0.0086	0.0068
1.05	30	0.0047	0.0000	0.2926	0.0695	0.0002	0.0307	0.0281	0.0201
	50	0.0037	0.0000	0.1730	0.0456	0.00002	0.0144	0.0170	0.0128
	100	0.0030	0.0000	0.0724	0.0209	0.00002	0.0044	0.0086	0.0068

Table 11. Comparing simulation results about the variances and MVB results. $\left(\mathit{\alpha }=2,\mathit{\beta }=1,\mathit{b}=-2\right).$

$a$	$n$	Simulation results about the variances				MVB results
		$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$	$\widehat{a}$	$\widehat{b}$	$\widehat{\alpha }$	$\widehat{\beta }$
0.95	30	0.00001	0.0051	0.3157	0.0783	0.0000	0.0048	0.0281	0.0201
	50	0.0004	0.0000	0.1481	0.0422	0.0009	0.0092	0.0172	0.0134
	100	0.0001	0.0000	0.0786	0.0224	0.0002	0.0017	0.0086	0.0066
0.99	30	0.0018	0.0000	0.3117	0.0686	0.0017	0.0303	0.0280	0.0195
	50	0.0010	0.0000	0.1762	0.0456	0.0008	0.0143	0.0171	0.0130
	100	0.00051	0.0000	0.0820	0.0233	0.0002	0.0046	0.0086	0.0068
1.01	30	0.0024	0.0000	0.3546	0.0715	0.0012	0.0318	0.0279	0.0191
	50	0.0016	0.0000	0.1675	0.0421	0.0005	0.0157	0.0171	0.0129
	100	0.0011	0.0000	0.0716	0.0217	0.0001	0.0060	0.0086	0.0068
1.05	30	0.0047	0.0000	0.2926	0.0695	0.0002	0.0307	0.0281	0.0201
	50	0.0037	0.0000	0.1730	0.0456	0.00002	0.0144	0.0170	0.0128
	100	0.0030	0.0000	0.0724	0.0209	0.00002	0.0044	0.0086	0.0068

Table 12. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0000	0.0003	0.2876	0.0582
	50	0.0005	0.0092	0.1309	0.0288
	100	0.0001	0.0017	0.07	0.0158
0.99	30	0.0001	0.0303	0.2837	0.0491
	50	0.0002	0.0143	0.1591	0.0326
	100	0.0003	0.0046	0.0734	0.0165
1.01	30	0.0012	0.0318	0.3267	0.0524
	50	0.0011	0.0157	0.1504	0.0292
	100	0.0010	0.006	0.063	0.0149
1.05	30	0.0045	0.0307	0.2645	0.0494
	50	0.0036	0.0144	0.156	0.0328
	100	0.0029	0.0044	0.0638	0.0141

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

Table 12. The differences between the simulation results and the MVB results.

$a$	$n$	$\widehat{\Delta a}$	$\widehat{\Delta b}$	$\widehat{\Delta \alpha }$	$\widehat{\Delta \beta }$
0.95	30	0.0000	0.0003	0.2876	0.0582
	50	0.0005	0.0092	0.1309	0.0288
	100	0.0001	0.0017	0.07	0.0158
0.99	30	0.0001	0.0303	0.2837	0.0491
	50	0.0002	0.0143	0.1591	0.0326
	100	0.0003	0.0046	0.0734	0.0165
1.01	30	0.0012	0.0318	0.3267	0.0524
	50	0.0011	0.0157	0.1504	0.0292
	100	0.0010	0.006	0.063	0.0149
1.05	30	0.0045	0.0307	0.2645	0.0494
	50	0.0036	0.0144	0.156	0.0328
	100	0.0029	0.0044	0.0638	0.0141

$\Delta :\left|Simulationvariances-minimumvariancebounds\right|.$

5. Illustrative Examples

In this section, two real data sets of different sizes in various types of areas are used to calculate the estimators of the model parameters and investigate the advantages and disadvantages of ML and MM methods. Let’s assume $Y_k\text{'}s$ are iid with $E\left(Y_k\right)=\mu ,V\left(Y_k\right)={\sigma }^2$, $Y_k=a^{k-1}{X_k}^{h\left(k\right)}$

$${\widehat{a}}_{MM}=\mathit{exp}\left(lna\right)$$

(6)

where

$$lna={\displaystyle \frac{\mathsf{\lambda }\mathrm{*}3\mathrm{n}\mathrm{*}\left(\mathrm{n}-1\right)-6\sum _{k=1}^n\left(k-1\right){\left(1+log\left(k\right)\right)}^{b}ln\left(X_k\right)}{n\left(n-1\right)\left(2n-1\right)}}$$

(7)

$$\mathsf{\lambda }={\displaystyle \frac{\left(4n-2\right)\sum _{k=1}^n{\left(1+log\left(k\right)\right)}^{b}ln\left(X_k\right)-6\sum _{k=1}^n\left(k-1\right){\left(1+log\left(k\right)\right)}^{b}ln\left(X_k\right)}{n\left(n+1\right)}}$$

(8)

$${\widehat{{\sigma }^2}}_{MM}=\left\{\begin{array}{c}{\sum }_{k=1}^n{\displaystyle \frac{{\left(\widehat{y_k}-\overline{\widehat{y_k}}\right)}^2}{n-1}},a\ne 1,h\left(k\right)\ne 1,b\ne 0\\ {\sum }_{k=1}^n{\displaystyle \frac{{\left(x_k-\overline{x}\right)}^2}{n-1}},a=1,h\left(k\right)=1,b=0\end{array}\right.$$

(9)

where

$$\widehat{y_k}={\widehat{a}}_{MM}{x}_k^{{\widehat{h\left(k\right)}}_{MM}},\overline{\widehat{y_k}}={\sum }_{k=1}^n{\displaystyle \frac{\widehat{y_k}}{n}},\overline{x}={\sum }_{k=1}^n{\displaystyle \frac{x_k}{n}}$$

For a data set $\left\{X_1,\dots ,X_n\right\}$, the $\widehat{X_k}$ is calculated as follows,

$$\widehat{X_k}=\left\{\begin{array}{c}\\ {{\widehat{a}}_{ML}}^{\left(1-k\right)h^{-1}\left(k\right)}{\lambda }_{1k},byaDGPwithMLestimators\\ {\left({\displaystyle \frac{\mu }{{{\widehat{a}}_{MM}}^{k-1}}}\right)}^{{\left(1+log\left(k\right)\right)}^{-b_{MM}}},byaDGPwithMMestimators\end{array}\right.$$

(10)

$${\lambda }_{1k}={\int }_0^{\infty }{x}^{h^{-1}\left(k\right)}f\left(x\right)dx$$

(11)

The mean squared error (MSE) and maximum percentage error (MPE) criteria are used to evaluate the performance of the DGPs with the ML and MM estimators, see Lam et al. [22]. The $S_k,\widehat{S_k},\widehat{Y_k},\widehat{\mu },\widehat{{\sigma }^2}$, and MSE, MPE criteria’s equalities, which are used in real-world data applications, are given as follows.

$$S_k=X_1+\dots +X_k,k=\mathrm{1,2},\dots ,n$$

(12)

$$\widehat{S_k}={\sum }_{k=1}^n\widehat{X_k}$$

(13)

$${MSE}^{*}={\displaystyle \frac{1}{n}}{\sum }_{k=1}^n{\left(X_k-\widehat{X_k}\right)}^2$$

(14)

$$MPE=\underset{1\le k\le n}{\underbrace{max}}{\displaystyle \frac{\left\{\left|S_k-\widehat{S_k}\right|\right\}}{S_k}}$$

(15)

The estimations of the parameters $\mu$ and ${\sigma }^2$ are obtained as follows

$${\widehat{\mu }}_{MM}=\overline{\widehat{Y_k}}={\sum }_{k=1}^n{\displaystyle \frac{{\widehat{a}}^{k-1}{X_k}^{\widehat{h\left(k\right)}}}{n}}$$

(16)

$${\widehat{{\sigma }^2}}_{MM}=\left\{\begin{array}{c}{\sum }_{k=1}^n{\displaystyle \frac{{\left(\widehat{y_k}-\overline{\widehat{y_k}}\right)}^2}{n-1}},a\ne 1,h\left(k\right)\ne 1,b\ne 0\\ {\sum }_{k=1}^n{\displaystyle \frac{{\left(x_k-\overline{x}\right)}^2}{n-1}},a=1,h\left(k\right)=1,b=0\end{array}\right.$$

(17)

[13].

The following approach is used to investigate whether the data set can be fitted by a DGP with the gamma distribution or not. Firstly, it is assumed that for the DGP model, the datasets follow a DGP with a gamma distribution. Then the K-S test is used to see whether the DGP model can represent the dataset or not. The ML and LSE estimators’ values are obtained for each parameter.

The prediction values are obtained by the formula ${\widehat{y}}_k={\widehat{a}}_{DGP}^{k-1}{x}_k^{{\left(1+\log \left(k\right)\right)}^{{\widehat{b}}_{DGP}}}k=1,\dots ,n$where

$\widehat{a}:$ the ML estimator of the model parameter $a$

$\widehat{b}:$ the ML estimator of the model parameter $b$

The goodness of fitness measurement values of the ${\widehat{y}}_k\text{'}s$ for gamma distribution with parameters $\alpha$, $\beta$are calculated for each model by the $KS$ test. As a result, if it is seen that the predictions have gamma random variables by the KS test statistic, then it can be said that the DGP can model the data set with a particular gamma distribution.

Example1. Coal mining disasters data.

The first data set, originally studied by Maguire et al. [26] (Data Set No.1), shows the interarrival time between coal mining disasters. The data set is given in Table 13. The results for each method are given in Table 15.

For this data, the KS test statistic value and the corresponding p value are given in Table 14.

It is seen in Table 14 that the gamma distribution is appropriate for the data set 1.

When the data set 1 is modeled by a DGP with the gamma distribution, the ML estimator values of the $a,b,\alpha ,\beta$ parameters, the MM estimates of the parameters $a,b,$MSE, and MPE values, and the estimates of $E\left(Y_k\right)=\mu ,Var\left(Y_k\right)={\sigma }^2$ are given in Table 15.

It can be seen from Table 15 that ML estimators have the smallest MSE and MPE values. In Figure 1, $S_k$ and $\widehat{S_k}$ disaster times are plotted.

Figure 1 is compatible with Table 15.

Example 2. Main propulsion diesel engine failure data.

The second data set is called the U.S.S. Grampus No. 4 main propulsion diesel engine failure data. (Data set No.2) contains the cumulative operating hours until significant maintenance events occur for one of the thirty engines on nine different submarines. The data set was studied by Lee[25]. It was originally studied by the Corporation [10]. The sample size of the data set is given in Table 16. The results for each method are given in Table 18.

For this data, the KS test statistic value and the corresponding p value are given in Table 17.

It is seen in Table 17 that the gamma distribution is an appropriate model for dataset set2.

When the data set 2 is modeled by a DGP with the gamma distribution, the ML estimator values of the $a,b,\alpha ,\beta$ parameters, the MM estimates of the parameters $a,b,$ MSE, and MPE values, and the estimates of $E\left(Y_k\right)=\mu ,Var\left(Y_k\right)={\sigma }^2$ are given in Table 18.

It can be seen from Table 18 that ML estimators have the smallest MSE and MPE values. In Figure 1, $S_k$ and $\widehat{S_k}$ disaster times are plotted.

Figure 2 is compatible with Table 18.

6. Conclusions

In this study, the parameter estimation problem for a DGP is considered under the assumption that the first inter-arrival time has a gamma distribution. The estimators of the parameters $a,b,\alpha ,and\beta$are calculated by the ML method. The asymptotic joint distributions of the ML estimators are obtained. The unbiasedness and consistency statistical properties of the estimators are investigated. The Monte Carlo simulation study was performed to evaluate the performance of the estimators. The simulation study supports the fact that ML estimators are highly efficient and unbiased. Finally, two data sets are considered in applications. The ML and MM estimators’ values, MSE, MPE measures, and the values of the K-S are calculated. As a result, it is seen that the gamma distribution is an appropriate model for these data sets. It can be easily said that ML estimators with smaller MSE and MPE measure values outperform the MM estimators.