EUROPEAN JOURNAL OF PURE AND APPLIED MATHEMATICS Vol. 15, No. 3, 2022, 1280-1300 ISSN 1307-5543 – ejpam.com Published by New York Business Global A New Sushila Distribution: Properties and Applications Somchit Boonthiem1, Adisak Moumeesri2, Watcharin Klongdee2, Weenakorn Ieosanurak2,∗ 1 Mathematics and Statistics Program, Sakon Nakhon Rajabhat University, Sakon Nakhon, Thailand 2 Department of Mathematics, Faculty of Science, Khon Kaen University, Thailand Abstract. In this paper, we introduce a new continuous distribution mixing exponential and gamma distributions, called new Sushila distribution. We derive some properties of the distribution include: probability density function, cumulative distribution function, expected value, moments about the origin, coefficient of variation (C.V.), coefficient of skewness, coefficient of kurtosis, moment generating function, and reliability measures. The distribution includes, a special cases, the Sushila distribution as a particular case p = 1 2 (θ = 1). The hazard rate function exhibits increasing. The parameter estimations as the moment estimation (ME), the maximum likelihood estimation (MLE), nonlinear least squares methods, and genetic algorithm (GA) are proposed. The application is presented to show that model to fit for waiting time and survival time data. Numerical results compare ME, MLE, weighted least squares (WLS), unweighted least squares (UWLS), and GA. The results conclude that GA method is better performance than the others for iterative methods. Although, ME is not the best estimate, ME is a fast estimate, because it is not an iterative method. Moreover, The proposed distribution has been compared with Lindley and Sushila distributions to a waiting time data set. The result shows that the proposed distribution is performing better than Lindley and Sushila distribution. 2020 Mathematics Subject Classifications: 60-08, 62-08, 49M99 Key Words and Phrases: Exponential distribution, Gamma distribution, Least square method, Maximum likelihood estimation, Method of moment, Genetic algorithm, Estimator 1. Introduction Finite mixture models are a widespread tool to heterogeneous data analytic for across a broad number of fields including agriculture, botany, bioinformatics, cell biology, genetics, genomics, genealogy, paleontology, zoology, fisheries research, economics, machine learn- ing, medicine, psychology, insurance, physics, engineering and chemistry, among many ∗Corresponding author. DOI: https://doi.org/10.29020/nybg.ejpam.v15i3.4420 Email addresses: somchit math@snru.ac.th (S. Boonthiem), adisak.m@kkumail.com (A. Moumeesri), kwatch@kku.ac.th (W. Klongdee), weenie@kku.ac.th (W. Ieosanurak) https://www.ejpam.com 1280 © 2022 EJPAM All rights reserved. W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1281 others. Because of their flexibility, finite mixture models can be used to cluster data, classify data, estimate densities, and increasingly. An enormous number of research have presented some approaches into finite mixture models and many methods to fitting lifetime data for improve performance of parameters are used. The Bayesian approach is used into [17–20, 22, 23, 25]. The Maximum Likelihood Estimation approach is used into [7, 8, 10, 21, 23–25]. The Least Square Estimator and Weighted Least Square Estimator are used into [7, 10, 21]. The Lagrangian multiplier approach is used into [5] and the ARMA approach is used into [9]. There are focus on the parameter estimation which is important to bring into the next step on application. A continuous random variable X which is called to have an exponential distribution with one parameter λ > 0, often called rate parameter, denoted by X ∼ Exp(λ). The probability density function (PDF) of the distribution is given by f(x;λ) = λe−λx, x ≥ 0, and the cumulative distribution function (CDF), has been obtained as F (x) = 1− e−λx, x ≥ 0. A continuous random variable X is called to have a gamma distribution with two parameters α > 0, often called shape parameter, and λ > 0, often called rate parameter, denoted by X ∼ Gamma(α, λ). The PDF of the distribution is given by f(x;α, λ) = λαxα−1e−λx Γ(α) , x > 0, where Γ(α) is gamma function. We can see that Gamma(1, λ) = Exp(λ). Its CDF is given by F (x) = 1 Γ(α) γ(α, λx), x > 0. where γ(α, λx) is the lower incomplete gamma function. There are studies that have used exponential, Lindley and gamma distributions to be used for modeling lifetime data. A two-parameter Lindley distribution is proposed by modifying the mixing weight of exponential and gamma [6, 11–16]. They proposed a two-parameter distribution by modifying the mixing weight of Exp( θα) and Gamma(2, θ α). A continuous random variable X is called to have Sushila distribution (SD) [16] with two parameters θ > 0 and α > 0, denoted by X ∼ SD(α, θ). The PDF of the distribution is given by f(x;α, θ) = θ2 α(θ + 1) (1 + x α )e− θ α x, x > 0. Its CDF is given by F (x) = 1− α(θ + 1) + θx α(θ + 1) e− θ α x, x > 0. W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1282 The PDF of Sushila distribution can be shown as a mixture of Exp( θα) and Gamma(2, θ α) as follows: f(x, α, θ) = pf1(x) + (1− p)f2(x), where p = θ θ + 1 , f1(x) = θ α e− θ α x, and f2(x) = θ2 α2 xe− θ α x. In this paper, we focus on the Sushila distribution. The distribution is mixed by exponential and gamma distributions. Sushila distribution therefore pays attention to the weight value for exponential distribution. Notice that the gamma distribution has two parameters that provide a flexible model to fit the reliability data. For this reason, we pay attention to the weight value for gamma distribution. The main aim of this paper is to study a two-parameter (θ > 0 and α > 0) continuous distribution is introduced as follows: f(x, α, θ) = (1− p)f1(x) + pf2(x). Its CDF, expected value, first four moments, and some of the related measures have been proposed. The paper is organized as follows. We define new Sushila distribution and its properties in Section 2. The moments, generating function, and related measures are derived in Section 3. The reliability measures, like survival function, hazard function, and mean residual life function for new Sushila distribution in Section 4. Five methods of estimation are discussed in Section 5. The numerical simulation and comparison are given to the proposed methods in Section 6. Conclusion is discussed in Section 7. 2. Definition and Properties of the new Sushila Distribution In this section, we introduce a continuous distribution with two parameters α and θ, called new Sushila distribution. The PDF of new Sushila distribution can be shown as a mixture of Exp ( θα) and Gamma(2, θ α) distributions as follows: f(x;α, θ) = (1− p)f1(x) + pf2(x), (1) where p = θ θ + 1 , f1(x) = θ α e− θ α x, and f2(x) = ( θ α )2 xe− θ α x. Definition 1. A continuous random variable X is said to be a new Sushila distribution random variable with two parameters α and θ, its PDF is given by f(x;α, θ) = θ(α+ θ2x)e− θ α x α2(θ + 1) , for all x > 0, α > 0, θ > 0. (2) Probability plots of new Sushila distribution are given in Figure 1 for particular values of α and θ. From Figure 1, we can see that new Sushila distribution is a light-tailed distribution and the tail approaches to zero at a faster rate for lower values of α. It is W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1283 Figure 1: Plots of the PDF of new Sushila distribution for various parameter values. evident from the density plots of new Sushila distribution that the distribution is most likely to appropriately model those data sets where the observations of lower magnitude appear more frequently than that of the higher magnitude values. Such data sets are waiting time, claim and lifetime. The first derivative of Equation (2) with respect to x is obtained as d dx f(x;α, θ) = θ2(αθ − α− θ2x)e− θ α x α3(θ + 1) . (3) Now, based on Equation (3), we obtain, 1. f ′(x) = 0 gives x = α(θ − 1) θ2 as the critical point. For θ > 1, x0 = α(θ − 1) θ2 is the unique critical point at which f(x) is maximum. 2. for θ < 1, f ′(x) ≤ 0 i.e. f(x) is decreasing in x. Therefore, the mode of the distribution (2) is given by Mode =  α(θ − 1) θ2 if θ > 1, 0 otherwise. Next, we obtain the cumulative distribution function of new Sushila distribution as follows. Theorem 1. The cumulative distribution function (CDF) of new Sushila distribution is given by F (x) = 1− (αθ + α+ θ2x)e− θ α x α(θ + 1) , x > 0, θ > 0, α > 0. (4) W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1284 Proof. For all x > 0, α > 0 and θ > 0, consider F (x) = P (X ≤ x) = 1− P (X ≥ x) = 1− ∫ ∞ x θ(α+ θ2t)e− θ α t α2(θ + 1) dt = 1− (αθ + α+ θ2x)e− θ α x α(θ + 1) . The CDF plots of new Sushila distribution are given in Figure 2 for particular values of α and θ. Figure 2: Plots of the CDF of new Sushila distribution for various parameter values. 3. Moments, Generating function and Related Measures In this section, we derive the rth moment about origin, moment generating function, the coefficient of variation, coefficient of skewness, and coefficient of kurtosis. Theorem 2. The rth moment about the origin of new Sushila distribution is defined by, µ′ r = r!αr[(r + 1)θ + 1] θr(θ + 1) , r = 1, 2, 3, . . . . Proof. By definition, for r = 1, 2, . . ., the moment about origin, µ′ r is obtained as, µ′ r = E[Xr] = ∫ ∞ 0 xr θ(α+ θ2x)e− θ α x α2(θ + 1) dx W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1285 = θ α2(θ + 1) [∫ ∞ 0 αxre− θ α x dx+ θ3 ∫ ∞ 0 xr+1e− θ α x dx ] = θ α2(θ + 1) [ α ∫ ∞ 0 ( xθ α )r (α θ )r e− θ α x dx+ θ3 ∫ ∞ 0 ( xθ α )r+1 (α θ ) e− θ α x dx ] = 1 α(θ + 1) ( αr+1 θr Γ(r + 1) + αr+1 θr−2 Γ(r + 2) ) = αr+1 α(θ + 1) ( r! θr + 1 θr−2 (r + 1)! ) = r!αr[(r + 1)θ + 1] θr(θ + 1) . The first four moments about origin of new Sushila distribution can be obtained as: µ′ 1 = α(2θ + 1) θ(θ + 1) , µ′ 2 = 2α2(3θ + 1) θ2(θ + 1) , µ′ 3 = 6α3(4θ + 1) θ3(θ + 1) , µ′ 4 = 24α4(5θ + 1) θ4(θ + 1) . In particular, the 1st moment about origin take the form, µ = E[X] = α(2θ + 1) θ(θ + 1) . Using the relationship between moments about mean and the moments about the origin, we have the moment about mean as the following theorem. Theorem 3. The first four moments about mean of new Sushila distribution is defined by, µ2 = α2(2θ2 + 4θ + 1) θ2(θ + 1)2 , µ3 = 2α3(2θ3 + 6θ2 + 6θ + 1) θ3(θ + 1)3 , µ4 = 3α4(8θ4 + 32θ3 + 44θ2 + 24θ + 3) θ4(θ + 1)4 . Proof. For x > 0, α > 0 and θ > 0. Consider µ2 = E[(X − µ)2] = E[X2]− µ2 = 2α2(3θ + 1) θ2(θ + 1) − ( α(2θ + 1) θ(θ + 1) )2 W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1286 = α2(2θ2 + 4θ + 1) θ2(θ + 1)2 , µ3 = E[X3]− 3µE[X2] + 2µ3 = 2α2(3θ + 1) θ2(θ + 1) − 3 α(2θ + 1) θ(θ + 1) [ 2α2(3θ + 1) θ2(θ + 1) ] + 2 [ α(2θ + 1) θ(θ + 1) ]3 = 2α3(2θ3 + 6θ2 + 6θ + 1) θ3(θ + 1)3 , µ4 = E[(X − µ)4] = E[X4]− 4µE[X3] + 6µ2E[X2]− 3µ4 = 24α4(5θ + 1) θ4(θ + 1) − 4 ( α(2θ + 1) θ(θ + 1) )( 6α3(4θ + 1) θ3(θ + 1) ) + 6 ( α(2θ + 1) θ(θ + 1) )2(2α2(3θ + 1) θ2(θ + 1) ) − 3 ( α(2θ + 1) θ(θ + 1) )4 = 3α4(8θ4 + 32θ3 + 44θ2 + 24θ + 3) θ4(θ + 1)4 . In particular, the 2nd moment about the mean is variance that is, σ2 = µ2 = α2(2θ2 + 4θ + 1) θ2(θ + 1)2 . The coefficient of variation (C.V.), coefficient of skewness ( √ β1), and coefficient of kurtosis (β2) of new Sushila are given by C.V. = σ µ′ 1 = √ 2θ2 + 4θ + 1 2θ + 1 ,√ β1 = µ3 µ 3/2 2 = 4θ[θ(θ + 3) + 3] + 2 θ6(θ + 1)6[2θ(θ + 2) + 1]3/2 , β2 = 3(8θ4 + 32θ3 + 44θ2 + 24θ + 3) (2θ2 + 4θ + 1)2 . Now, we derive the moment generating function (MGF) of new Sushila distribution as the following. Theorem 4. The moment generating function MX(t) of new Sushila distribution is de- fined by MX(t) = E[etX ] = θ(αt+ θ2 − θ) (θ + 1)(θ − αt)2 , θ α > t. W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1287 Proof. MX(t) = E[etX ] = θ α2(θ + 1) ∫ ∞ 0 (α+ θ2x)e(t− θ α )x dx = θ α2(θ + 1) lim b→∞ [ α(α2t+ αθ + αθ(θ(xt− 1)− 1)) + θ3(−x) (θ − αt)2 e(t− θ α )x ]b x=0 = θ(αt+ θ2 − θ) (θ + 1)(θ − αt)2 , θ α > t. 4. Reliability Measures In this section, we derive expressions for the reliability measures, like survival function, hazard function, and mean residual life function for new Sushila distribution. Theorem 5. For x > 0, θ > 0, α > 0. The survival function S(x), the hazard rate function h(x), and mean residual life function m(x) of new Sushila distribution is defined by, S(x) = (αθ + α+ θ2x)e− θ α x α(θ + 1) , h(x) = θ(α+ θ2x) α(αθ + α+ θ2x) , and m(x) = α2(θ + 1)(2αθ + α+ θ2x) αθ + α+ θ2x . Proof. By definition of the survival function, we have S(x) = 1− F (x) = (αθ + α+ θ2x)e− θ α x α(θ + 1) , x > 0, θ > 0, α > 0. By definition of the hazard rate function, we have h(x) = lim ∆x→0 P (X < x+∆x|X > x) ∆x = f(x) 1− F (x) = θ(α+ θ2x) α(αθ + α+ θ2x) . By definition of the mean residual life function, we have m(x) = E[X − x|X > x] W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1288 = 1 1− F (x) ∫ ∞ x (1− F (t)) dt = α2(θ + 1)(2αθ + α+ θ2x) αθ + α+ θ2x . It can be easily verified that h(0) = θα θ + α2 = f(0) and m(0) = µ′ 1. The derivative h′(x) > 0 for all x, so h(x) is an increasing function of x, α and θ, whereas m(x) is a decreasing function of x, α and θ. 5. Estimation of Parameters In this section, we proposed five methods for estimating parameters. The first two methods are the moment estimation (ME), the maximum likelihood estimates (MLE). These methods are widely used for estimation because they are simple. The next two methods are the weighted and unweighted least squares methods via the CDF. These methods are more complicated than the first two methods, but the least squares are better than in some situations. The last method is the genetic algorithm which is a heuristic search that mimics the process of natural evolution. This method is an efficient and effective technique. 5.1. The Moment Estimation (ME) The new Sushila distribution has two parameters to be estimated. Using the first moment about origin, we have E[X] = α(2θ + 1) θ(θ + 1) , E[X] = X̄. We obtain, X̄ = α(2θ + 1) θ(θ + 1) . So α̂ = θ(θ + 1)X̄ (2θ + 1) . Using the second moment about mean, we have µ2 = E[X2]− µ2 = α2(2θ2 + 4θ + 1) θ2(θ + 1)2 , W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1289 and E[X2]− µ2 = 1 n n∑ i=1 x2i − X̄2, so 1 n n∑ i=1 x2i − X̄2 = α2(2θ2 + 4θ + 1) θ2(θ + 1)2 , 1 n n∑ i=1 x2i − X̄2 = ( θ(θ + 1)X̄ (2θ + 1) )2 (2θ2 + 4θ + 1) θ2(θ + 1)2 . Solving this equation for θ, we get θ̂ = ( 2X̄2 − 1 n n∑ i=1 x2i ) ± √√√√X̄2 ( X̄2 − 1 2n n∑ i=1 x2i ) 2 n n∑ i=1 x2i − 3X̄2 , 2 n n∑ i=1 x2i − 3X̄2 ̸= 0, where α̂, θ̂ are the estimators of the parameter θ and α, respectively. This method will be referred to as Method 1. 5.2. The Maximum Likelihood Estimates (MLE) Let x1, x2, . . . , xn be a random sample of size n from new Sushila distribution and let fx be the observed frequency in the sample corresponding to X = x (x = 1, 2, ..., k) such that ∑k x=1 fx = n, where k is the largest observed value having non-zero frequency. The likelihood function, L of new Sushila distribution is given by L = ( θ α2(θ + 1) )n k∏ x=1 ( α+ θ2x )fx e− θ α (nx̄), and so the log likelihood function is obtained as logL = n log θ − 2n logα− n log(θ + 1) + k∑ x=1 fx log(α+ θ2x)− nθx̄ α . (5) By differentiating Equation (5) with respect to θ and α, is obtained by: ∂ logL ∂θ = 0, ∂ logL ∂α = 0. (6) W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1290 The equation (6) do not seem to be solved directly. However, the Fisher’s scoring method can be applied to solve these equations. The following equations for θ̂ and α̂ can be solved ∂2 logL ∂θ2 ∂2 logL ∂θ∂α ∂2 logL ∂α∂θ ∂2 logL ∂α2   θ̂ − θ0 α̂− α0  =  ∂ logL ∂θ ∂ logL ∂α  , where θ0 and α0 are the initial values of θ and α respectively. These equations are solved iteratively till sufficiently close estimates of θ̂ and α̂ are obtained. This method will be referred to as Method 2. 5.3. The Unweighted Least Squares (UWLS) Method via the CDF Let X1, X2, . . . , Xn be n independent random variables having the new Sushila distri- bution with parameters α and θ. Without loss of generality, suppose that X1 < X2 < · · · < Xn are the order statistics and by taking the natural logarithm on the both sides of Equation (4), we obtain that log(F (x)) = log ( 1− (αθ + α+ θ2x)e− θx α α(θ + 1) ) , x, α, θ > 0 = log ( (θ + 1)− (θ + 1 + θ2x α )e− θx α ) − log(θ + 1). Let 0 < x1 < x2 < · · · < xn be n orders observations, then log(F (xi)) = log ( (θ + 1)− (θ + 1 + θ2xi α )e− θxi α ) − log(θ + 1), x, α, θ > 0. Let the empirical distribution function of X be denoted by Fn(x), By reference [4], the estimator of F (xi) can be considered Fn(xi) = i− d n− 2d+ 1 , i = 1, 2, . . . , n, (7) for some real number d, 0 ≤ d ≤ 1. So, we choose four popular expressions (see also [1]) that are often used as estimators of F (xi) uik =  i n+1 , k = 1 i−0.3 n+0.4 , k = 2 i−0.375 n+0.25 , k = 3 i−r n , k = 4, W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1291 i = 1, 2, 3, . . . , n and for some real number r ∈ (0, 1). We estimate α and θ by the unweighted least squares method, by using minimizing function Ek(α, θ) = n∑ i=1 ( log(uik)− log ( (θ + 1)− (θ + 1 + θ2xi α )e− θxi α ) + log(θ + 1) )2 , (8) k = 1, 2, 3, 4. By solving ∂ ∂α Ek(α, θ) = 0 ∂ ∂θ Ek(α, θ) = 0, k = 1, 2, 3, 4. We denoted by A(xi, α, θ) = xi(α+ θ2xi) log ( −e− θxi α ( θ2xi α + θ + 1) + θ + 1 ) , Bk(xi, α, θ) = xi(α+ θ2xi) log(uik), C(xi, α, θ) = α(θ + 1)(e θxi α − 1)− θ2xi, Dk(xi, α, θ) = xi[log(uik) + log(θ + 1)− log(−e− θxi α ( θ2xi α + θ + 1) + θ + 1)] (θ + 1)[α(θ + 1)(e θxi α − 1)− θ2xi , where α̂ is the estimator of the parameter α with this method. Using some algebraic manipulations, θ̂ satisfies the following equation: θ̂ = exp  ∑n i=1 A(xi, α, θ) Ck(xi, α, θ) − ∑n i=1 Bk(xi, α, θ) C(xi, α, θ)∑n i=1 xi(α+ θ2xi) C(xi, α, θ) − 1, and α̂ = n∑ i=1 (θ + 1)θ2xiDk(xi, α, θ) Dk(xi, α, θ) , for k = 1, 2, 3, 4. These four UWLS methods via uik, k = 1, 2, 3, 4, will be referred to as Methods 3 through 6, respectively. 5.4. The Weighted Least Squares (WLS) Method via the CDF Using the same weight for all datum can be erroneous [2], therefore, we weigh the values in Equation (8) via a weighting factor. To find one of these factors we may use a W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1292 variance stabilizing transformation. This approach, proposed by Bickel and Doksum [3] relates the variance of log(uik), denoted by V ar(log(uik)), to the uncertainty of uik and that we denote by V ar(uik). We obtain that V ar(log(uik)) = ( ∂ log(uik) ∂uik )2 V ar(uik), which gives V ar(log(uik)) = ( 1 ∂uik )2 V ar(uik). Therefore, we use wik as a weighting factor, denoted by wik = u2ik and minimizing function Ek(α, θ) = n∑ i=1 wik ( log(uik)− log ( (θ + 1)− (θ + 1 + θ2xi α )e− θxi α ) + log(θ + 1) )2 , k = 1, 2, 3, 4. By solving ∂ ∂α Ek(α, θ) = 0 ∂ ∂θ Ek(α, θ) = 0. k = 1, 2, 3, 4. We denoted by A(xi, α, θ) = xi(α+ θ2xi) log ( −e− θxi α ( θ2xi α + θ + 1) + θ + 1 ) , Bk(xi, α, θ) = xi(α+ θ2xi) log(uik), C(xi, α, θ) = α(θ + 1)(e θxi α − 1)− θ2xi, Dk(xi, α, θ) = xi[log(uik) + log(θ + 1)− log(−e− θxi α ( θ2xi α + θ + 1) + θ + 1)] (θ + 1)[α(θ + 1)(e θxi α − 1)− θ2xi , where α̂ is the estimator of the parameter α with this method. Using some algebraic manipulations, θ̂ satisfies the following equation: θ̂ = exp  ∑n i=1wik A(xi, α, θ) Ck(xi, α, θ) − ∑n i=1wik Bk(xi, α, θ) C(xi, α, θ)∑n i=1wik xi(α+ θ2xi) C(xi, α, θ) − 1, and α̂ = n∑ i=1 (θ + 1)θ2xiwikDk(xi, α, θ) wikDk(xi, α, θ) , for k = 1, 2, 3, 4. These four WLS methods via uik, k = 1, 2, 3, 4 and weighting factors wik, will be referred to as Methods 7 through 10 respectively. W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1293 5.5. Genetic algorithm In the genetic algorithm process is as follows: Step 1. Determine the number of chromosomes, generation, and mutation rate and crossover rate value. The number of chromosomes is 2 (α and θ), the number of generations is 10000, and the mutation rate is 0.3. Step 2. Generate chromosome-chromosome number of the population, and the initialization value of the genes chromosome-chromosome with a random value Step 3. Process steps 4-7 until the number of generations is met Step 4. Evaluation of fitness value of chromosomes by calculating objective function. The fitness values is defined by: fi = 1 i , where i is current chromosomes. Step 5. Chromosomes selection. We used the roulette wheel selection. The probability of choosing individual i is equal to pi = fi∑N i=1 fi where fi is the fitness value of i and N is the size of the current generation. Step 6. Crossover. A point on both parents’ chromosomes is picked by chromosome 2. Step 7. Mutation. The chromosome i for i = 1, 2 are mutated as follows: we random ui ∈ (0, 1), if ui < Mutation rate, then we mutated the chromosome i. Step 8. Solution (Best Chromosomes) This method will be referred to as Method 11. 6. A simulated Study In this section. We use two data sets to find the good parameter estimation. The data are waiting time (in minutes) of 100 bank customers [16] as Table 1, and (in days) of 72 guinea pigs infected with virulent tubercle bacilli [14] as Table 2. The comparison of new Sushila distribution with Lindley and Sushila distribution is proposed in Table 5. It can be observed from Table 3 that the estimates obtained via UWLS (Method 3-6) are very same, also the estimates obtained via WLS (Method 7-10), although the estimates obtained via ME (Method 1) is not the best estimates, this method is good estimate, directly calculate. Moreover, it is not an iterative method like the others. The best estimate can be obtained using GA (Method 11). W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1294 Table 1: Waiting time (in minutes) of 100 bank customers. Waiting time (minutes) Observed frequency 0 - 4.9 30 5 - 9.9 32 10 - 14.9 19 15 - 19.9 10 20 - 24.9 5 25 - 29.9 1 30 - 34.9 2 35 - 39.9 1 Total 100 Table 2: Data of survival time (in days) of 72 guinea pigs infected with virulent tubercle bacilli. Survival time (days) Observed frequency 0 - 79 8 80 - 159 30 160 - 239 18 240 - 319 8 320 - 399 4 400 - 479 3 480 - 559 1 Total 72 It can be observed from Table 4 that the estimates obtained via UWLS (Method 3-6) and WLS (Method 7-10) are very same, although the estimates obtained via ME (Method 1) have the largest χ2, this estimates is directly calculated and it is not an iterative method like the others. MLE can not calculate the parameters. The best performance can be obtained using GA (Method 11). Table 5 shows comparison parameter estimation for the waiting time data set. The estimator of the Lindley distribution has been obtained by the moment estimation (ME). The estimators of the Sushila distribution have been obtained by the moment estimation (ME) and the maximum likelihood estimates (MLE). The estimator of new Sushila dis- tribution has been obtained by the genetic algorithm (GA). The results show that new Sushila distribution gives better performance than the Lindley and Sushila distributions for waiting time data set. The survival time data set in Table 2, the estimator of the Sushila distribution is obtained by the moment estimation (ME). We found that the parameter θ̂ < 0, the θ̂ does not satisfy the definition of Sushila distribution. Therefore we can not find the estimators for Sushila distribution. However, the new Sushila distibution can find the estimators for this data as the following Table 4. W. Ieosanurak et al. / Eur. J. Pure Appl. Math, 15 (3) (2022), 1280-1300 1295 Table 3: ME, MLE, UWLS, WLS, and GA estimators of α and θ for waiting time Method α θ χ2 1 28.099565 5.366360 2.317336 2 0.000398 0.000032 10.140794 3 2.977747 0.408350 7.017687 4 5.230973 0.772538 4.996613 5 3.334206 0.437684 6.528164 6 5.470875 0.796321 4.878338 7 4.055956 0.644292 7.033345 8 10.583034 1.606832 3.627089 9 2.576107 0.303251 7.314763 10 1.670607 0.187810 7.806489 11 22 4 2.083757 Table 4: ME, MLE, UWLS, WLS, and GA estimators of α and θ for survival time Method α θ χ2 1 389.388353 0.267269 27.842259 2 - - - 3 319.472049 3.114969 14.267841 4 322.349246 3.125277 14.271553 5 321.338306 3.121736 14.269368 6 318.985074 3.113425 14.267277 7 318.245084 3.110803 14.267468 8 316.905403 3.106044 14.268873 9 315.503637 3.101055 14.271810 10 318.441567 3.111500 14.267377 11 9438.901726 103 7.936444 7. Conclusion A two-parameter continuous distribution, called new Sushila distribution has been presented. Some properties of the distribution such as CDF, expected value, the rth moment, and estimation of parameters by nonlinear least squares methods, the maximum likelihood estimation (MLE), the moment estimation (ME), and genetic algorithm (GA) have been discussed. We found that the distribution contains the Sushila distribution as a particular case p = 1 2 (θ = 1). Finally, the two popular methods to estimate the parameters of a probability distri- bution are maximum likelihood (MLE) and the method of moments (ME), and the other method of least squares (UWLS and WLS) with various uik and genetic algorithm (GA). We compare the 11 methods to estimate the parameters of the new Sushila distribution. In experiment, we observed that the best parameter estimation is obtained via GA. Besides, if we want a fast method to estimate the parameters, we introduce ME for estimating parameters; moreover, we present an application of new Sushila distribution for fitting the REFERENCES 1296 Table 5: Comparing the estimation for waiting time Waiting time (minutes) Observed frequency Lindley (ME) Sushila (ME) Sushila (MLE) New Sushila (GA) 0 -4.9 30 29.8841 29.8295 32.3716 29.7298 5-9.9 32 30.2563 30.3055 21.7393 29.2572 10-14.9 19 18.8852 18.9084 14.5916 18.7491 15-19.9 10 10.0637 10.0653 9.7891 10.3586 20-24.9 5 4.9523 4.9468 6.5642 5.3034 25-29.9 1 2.3226 2.3168 4.3997 2.5920 30-34.9 2 1.0547 1.0506 2.9477 1.2277 35-39.9 1 0.4680 0.4655 1.9740 0.5685 Total 100 97.8870 97.8884 94.3772 97.7864 Paramters θ̂ = 0.1907 θ̂ = 17.9493 θ̂ = 0.1863 θ̂ = 4 α̂ = 213.1197 α̂ = 0.9758 α̂ = 22 χ2 2.3074 2.3173 10.1382 2.0813 waiting time data and survival time. Furthermore, We compare new Sushila distribution with Lindley and Sushila distributions for the waiting time data. The results show that new Sushila gives better fit than both the Lindley and Suhila distribution. Besides, the Sushila distribution cannot find the estimator for the survival time data. In future work, we will adjust the genetic algorithm for estimating parameter. On the other hand as there are other well established meta-heuristics, such as Bee and Ant Colony. We can use them for a better estimation. Acknowledgements This research is funded by the Young Researcher Development Project of Khon Kaen University Year 2022. The authors would like to thank the Editor and anonymous review- ers for their constructive comments that improved the final version of the paper. References [1] V Barnett. Probability plotting methods and order statistics. Journal of the Royal Statistical Society. Series C (Applied Statistics), 24(1):95–108, 1975. [2] B Bergman. Estimation of weibull parameters using a weight function. Journal of Materials Science Letters, 5(6):611–614, 1986. [3] P Bickel and K Doksum. Mathematical statistics: Basic ideas and selected topics. 2d. ed. vol. 1 prentice hall. Upper Saddle River, NJ, 2001. [4] R D’Agostino and M Stephens. Goodness-of-Fit Techniques. Marcel Dekker, Inc., USA, 1986. REFERENCES 1297 [5] Quentin Gemine, Bertrand Cornélusse, Mevludin Glavic, Raphael Fonteneau, and Damien Ernst. A gaussian mixture approach to model stochastic processes in power systems. In 2016 Power Systems Computation Conference (PSCC), pages 1–7. IEEE, 2016. [6] D V Lindley. Fiducial distributions and bayes’ theorem. Journal of the Royal Statis- tical Society, 20(31):102–107, 1985. [7] Showkat Ahmad Lone, Sadia Anwar, Tabassum Naz Sindhu, and Fahd Jarad. Some estimation methods for mixture of extreme value distributions with simulation and application in medicine. Results in Physics, 37:105496, 2022. [8] Showkat Ahmad Lone, Tabassum Naz Sindhu, Anum Shafiq, and Fahd Jarad. A novel extended gumbel type ii model with statistical inference and covid-19 applications. Results in Physics, 35:105377, 2022. [9] Shuwei Miao, Kaigui Xie, Hejun Yang, Rajesh Karki, Heng-Ming Tai, and Tao Chen. A mixture kernel density model for wind speed probability distribution estimation. Energy Conversion and Management, 126:1066–1083, 2016. [10] Anum Shafiq, Tabassum Naz Sindhu, and Naif Alotaibi. A novel extended model with versatile shaped failure rate: Statistical inference with covid-19 applications. Results in Physics, 36:105398, 2022. [11] R Shanker and AG Amanuel. A new quasi lindley distribution. International Journal of Statistics and systems, 8(2):143–156, 2013. [12] R Shanker, H Fesshaye, and S Sharma. On two parameter lindley distribution and its applications to model lifetime data. Biometrics & Biostatistics International Journal, 3(1):9–15, 2016. [13] R Shanker and A Mishra. A quasi lindley distribution. African Journal of Mathe- matics and Computer Science Research, 6(4):64–71, 2013. [14] R Shanker, S Sharma, and R Shanker. A two-parameter lindley distribution for modeling waiting and survival times data. 4:363–368, 2013. [15] R Shanker, S Sharma, U Shanker, and R Shanker. Janardan distribution and its application to waiting times data. Indian Journal of Applied Research, 3:500–502, 10 2011. [16] R Shanker, S Sharma, U Shanker, and R Shanker. Sushila distribution and its applica- tion to waiting times data. International Journal of Business Management, 3(2):1–11, 2013. [17] Tabassum Naz Sindhu, Muhammad Aslam, and Zawar Hussain. A simulation study of parameters for the censored shifted gompertz mixture distribution: A bayesian approach. Journal of statistics and management systems, 19(3):423–450, 2016. REFERENCES 1298 [18] Tabassum Naz Sindhu, Navid Feroze, and Muhammad Aslam. Bayesian analysis of the kumaraswamy distribution under failure censoring sampling scheme. International Journal of Advanced Science and Technology, 51:39–58, 2013. [19] Tabassum Naz Sindhu, Navid Feroze, and Muhammad Aslam. A class of improved informative priors for bayesian analysis of two-component mixture of failure time dis- tributions from doubly censored data. Journal of Statistics and Management Systems, 20(5):871–900, 2017. [20] Tabassum Naz Sindhu and Zawar Hussain. Predictive inference and parameter esti- mation from the half-normal distribution for the left censored data. Annals of Data Science, pages 1–15, 2020. [21] Tabassum Naz Sindhu, Zawar Hussain, Naif Alotaibi, and Taseer Muhammad. Es- timation method of mixture distribution and modeling of covid-19 pandemic. AIMS Mathematics, 7(6):9926–9956, 2022. [22] Tabassum Naz Sindhu, Zawar Hussain, and Muhammad Aslam. On the bayesian analysis of censored mixture of two topp-leone distribution. Sri Lankan Journal of Applied Statistics, 19(1):13, 2019. [23] Tabassum Naz Sindhu, Anum Shafiq, and Qasem M Al-Mdallal. Exponentiated trans- formation of gumbel type-ii distribution for modeling covid-19 data. Alexandria En- gineering Journal, 60(1):671–689, 2021. [24] Tabassum Naz Sindhu, Anum Shafiq, and Qasem M Al-Mdallal. On the analysis of number of deaths due to covid- 19 outbreak data using a new class of distributions. Results in Physics, 21:103747, 2021. [25] TN Sindhu, M Saleem, and M Aslam. Bayesian estimation for topp leone distribution under trimmed samples. Journal of Basic and Applied Scientific Research, 3(10):347– 360, 2013. APPENDIX 1299 Appendix The structure of genetic algorithm is presented, where the notation is used: Nc is the number of generation, Np is the number of population, Mr is the mutation rate, Cr is crossover rate, Gen is a number of the population, [a (i) k , b (i) k ] is the i current chromosomes and k is the current generation, f(X) is the objective function, S∗ is the best chromosomes, Pn(m, :) is the mth row of matrix Pn Initial the parameter (Nc,Np,Mr,Cr,Gen) and Set k = 1; While k ≤ Gen do: Generate a random number of chromosomes X (i) k = [ a (i) k b (i) k f(a (i) k , b (i) k ) ] and calculate f(X (i) k ) for each i = 1, 2, . . . , Np and defined Pk at k iteration as follows: Pk = [ X (1) k X (2) k . . . X (i) k X (Np) k ]′ Np×3 Calculate i∗ = argimin{f(a(i)k , b (i) k )}. Set g∗ = f(a (i∗) k , b (i∗) k ) and P ∗ = [ a (i∗) k b (i∗) k f(a (i∗) k , b (i∗) k ) ] . Calculate the fitness value f (i) = 1 i for each population i = 1, 2, . . . , Np. Calculate the roulette wheel selection pi = f (i)∑Np j=1 f (j) for i = 1, 2, . . . , Np. if (f(a (i∗) k , b (i∗) k ) < g∗) g∗ = f(a (k) n , b (k) n ); k = k + 1; end while j ≤ Np/2 % Selection for l = 1 : 2 generate a random number, ul ∈ (0, 1). if ( ∑m i=1 pi < ul < ∑m+1 i=1 pi) for m = 1, 2, . . . , Np Parent(l, :) = Pk(m, :) end end %Crossover offspring(1, :) = [Parent(1, 1) Parent(2, 2)] offspring(2, :) = [Parent(2, 1) Parent(1, 2)] end Pk+1(2j − 1, :) = [offspring(1, :) f(offspring(1, :))] Pk+1(2j, :) = [offspring(2, :) f(offspring(2 :))] %Mutation APPENDIX 1300 while i ≤ Np do: Perform mutation operation on Pk+1 with Mr. Replaced the bad population with the best population X (i∗) k in Pk+1. end end