Bivariate Generalized Half-Logistic Distribution: Properties and Its Application in Household Financial Affordability in KSA

: The generalized half-logistic distribution is ideal to ﬁt the lifetime of some products, such as ball bearings and electrical insulation. In this paper, we aim to extend this scope by creating a motivated bivariate version. We thus introduce the bivariate generalized half-logistic distribution using the Farlie Gumbel Morgenstern (FGM) copula, which is called the FGM bivariate generalized half-logistic distribution (FGMBGHLD for short). In particular, the FGMBGHLD ﬁnds application in describing bivariate lifetime datasets that have weak correlations between variables. Some statistical properties and functions of our new distribution, such as the product moments, moment generating function, reliability function, and hazard rate function, are derived. We discuss the maximum likelihood estimation method of the FGMBGHLD parameters. As an application of the FGMBGHLD in reliability, we consider the stress–strength model when the stress and strength random variables are dependent. We also derive the point and interval estimates of the stress–strength coefﬁcient. Finally, we use the data from the household income and expenditure survey of KSA 2018 for Saudi households by administrative region to demonstrate the practicability of the proposed model. A comparison with a modern bivariate Weibull distribution is performed.


Introduction
In various fields, such as life testing, reliability, and biological and engineering sciences, there is a need for flexible lifetime distributions with various probability density and hazard rate properties. To this end, Mudholkar et al. (1995) [1] introduced the exponentiated Weibull family of distributions, which includes unimodal distributions with bathtub hazard rates as well as a broader class of monotone hazard rates. Alternative distributions have been examined since, presenting slightly different features. Gupta and Kundu (1999) [2] proposed a generalized exponential distribution. Olopade (2008) [3] considered two distributions, named type-I and type-III generalized half-logistic distributions. Kantam et al. (2014) [4] proposed a type-II generalized half-logistic distribution (GHLD-II for short). For the purpose of this paper, a brief presentation of the GHLD-II is necessary. On the mathematical plan, the probability density function (PDF), cumulative distribution function (CDF), and reliability function of the GHLD-II with scale parameter σ and power parameter µ are given by (2) and Thus, the GHLD-II is developed through the exponentiation of the reliability function of the half-logistic distribution (see Balakrishnan (1985) [5]). The flexibility of the GHLD-II is mainly in the mode and tail of the distribution, making it an interesting distribution for the modeling of lifetime phenomena. It is proven to define a better model than the exponential, Weibull, and half-logistic models (see Kantam et al. (2014) [4]).
The first objective of this paper is to derive a comprehensive bivariate generalized half-logistic distribution (BGHLD for short) using the copula approach and study its statistical properties, such as PDF, CDF, product moments, moment generating function, and hazard rate function. Many authors discuss the same idea but other distributions; see Almetwally et al. (2020) [6], Almetwally and Muhammed (2020) [7], and Muhammed et al. (2021) [8]. In view of the impact of the GHLD-II in the recent literature, we derive that bivariate versions have a promising future in terms of modeling and data analysis. Now, in order to detail and motivate the construction of our BGHLD, let us present some basics of the notion of the copula. As a first approach, we can say that a copula is a multivariate CDF for which the marginal distribution of each variable is uniform on the interval (0, 1). It describes the dependence between random variables. The definitions below provide more technical details. Definition 1. Let us consider a random vector (X 1 , . . . , X d ) and the marginal CDFs denoted by F i (x) = P(X i < x), for i = 1, . . . , d. Then, using probability integral transform (PIT) for each component, the distribution of the random vector (U 1 , . . . , U d ) = (F 1 (X 1 ), . . . , F d (X d )) belongs to the (uni f (0, 1)) d family of distributions, and the copula related to (X 1 , . . . , X d ) is defined as the joint CDF of (U 1 , . . . , U d ), i.e., with (u 1 , . . . , u d ) ∈ [0, 1] d .
The Sklar theorem, established by Sklar (1959) [9], is pivotal in copula theory. It states that, for two random variables X 1 and X 2 with marginal CDFs F 1 (x 1 ) and F 2 (x 2 ) and marginal PDFS f 1 (x 1 ) and f 2 (x 2 ), respectively, the CDF and PDF of (X 1 , X 2 ) are given by and respectively, where c(u 1 , u 2 ) denotes the copula density related to C(u 1 , u 2 ), i.e., c(u 1 , u 2 ) = ∂ 2 C(u 1 , u 2 )/(∂u 1 ∂u 2 ). Gumbel (1960) [10] discussed one of the most popular parametric families of copulas, called the Farlie Gumbel Morgenstern (FGM) copula. The FGM copula and its density are specified by and respectively. The parameter θ can be thought of as a dependence parameter that is dependent on the underlying random variables, with the independent case being θ = 0. The FGM copula is thus simple, flexible, and can be adapted when dealing with the construction of bivariate distributions with complicated marginal distributions in terms of functions. It is used in our study to create our BGHLD, which we naturally call the FGMBGHLD. The second objective is to develop the maximum likelihood (ML) estimation method of the FGMBGHLD parameters. Finally, the third goal is to derive the corresponding stress-strength model, but when and how this makes sense: in the dependent case, which can occur in engineering, operations research, quality control, education, economics, and insurance. Domma and Giordano (2013) [11] provided an example. In this paper, we are interested in economics, where X and Y are household income and consumption, and R = P(Y < X) is a measure of household financial affordability.
This paper is organized as follows. In Section 2, the FGMBGHLD is described. In Section 3, we derive some statistical properties of the FGMBGHLD. In Section 4, we exploit the copula approach to take into account the dependence of stress and strength variables in evaluating R. In Section 5, the ML estimation method for the FGMBGHLD parameters is discussed. In Section 6, point and interval estimations for R are elaborated. In Section 7, a Monte Carlo simulation study is performed to study the behavior of different estimates. In Section 8, the estimation of R is applied to KSA data (year 2018) to measure the household financial affordability for Saudi households by administrative region, with comparison to a modern bivariate Weibull distribution. The conclusion of this paper appears in Section 9.

FGM Bivariate Generalized Half-Logistic Distribution (FGMBGHLD)
Applying the Sklar theorem as stated in Equations (6) and (7) with Equations (1) and (2), and the FGM copula in Equations (8) and (9), we obtain the CDF and PDF of a random vector (Y 1 , Y 2 ) following the FGMBGHLD. They are given by and respectively, with the restrictions of the variables and parameters already mentioned. In order to illustrate the effect of the dependence parameter θ on the shape of these functions, Figure 1 shows the three-dimensional plots of the PDF and CDF with different values of θ (positive and negative). From Figure 1, we see that the variable variations of θ play a significant role; the PDF can take different forms in the space, with various skewness and kurtosis.

Statistical Properties of the FGMBGHLD
Here, we discuss some statistical properties of the FGMBHLD as defined by Equations (10) and (11). The marginal distributions, product moments, moment generating function, conditional distribution, generating random variables, and reliability function are derived.

Marginal PDFs
From a random vector (Y 1 , Y 2 ) following the FGMBHLD, for i = 1, 2, the distribution of Y i has the following PDF: Thus, more concretely, Y 1 has the following PDF: and Y 2 has the following PDF: On the other hand, for i = j with i, j = 1, 2, the general formula for the conditional PDF of Y i given Y j = y j is where F i (y i ) and F j y j denote the CDFs of Y i and Y j , respectively. Similarly, the conditional CDF of Y i given Y j = y j is We omit their analytical expressions for the FGMBHLD for the sake of brevity.

Moment Generating Function
The moment generating function of (Y 1 , Y 2 ) following the FGMBHLD is obtained as where and where 2 F 1 (a, b; c; z) refers to the (generalized) hypergeometric function. The parameters t 1 and t 2 must be selected such that the above quantities exist in the mathematical sense, which is the case for t 1 ≤ 0 and t 2 ≤ 0 among other more technical cases.

Product Moments
To obtain the product moments about the origin of (Y 1 , Y 2 ) following the FGMBHLD, for any positive real numbers r 1 and r 2 , we calculatè where and It is understood that Γ(x) refers to the standard gamma function, with Γ(m + 1) = m! for any integer m. From the product moments, various measures of moment skewness and kurtosis can be presented. On this topic, see, for instance, Almetwally et al. (2020) [6], Almetwally and Muhammed (2020) [7], and Muhammed et al. (2021) [8].

Reliability and Hazard Rate Functions
The reliability function of a bivariate distribution with an associated copula is defined by the copula composed with its marginal reliability functions. See Osmetti and Chiodini (2011) [12]. Hence, based on (Y 1 , Y 2 ) following the FGMBHLD, it is expressed as where R 1 (y 1 ) and R 2 (y 2 ) denote the reliability functions of Y 1 and Y 2 , respectively. According to the FGM copula, we obtain For the FGMBHLD, the reliability function is Moreover, Basu (1971) [13] defined the bivariate hazard rate function as For the FGMBHLD, the hazard rate function is indicated as (31)

Reliability for Dependence Stress-Strength Model
Domma and Giordano (2013) [11] introduced the concept of dependence via the stressstrength model. They calculated the reliability measure under the hypothesis that the bivariate distribution of the stress and strength variables, modeled by the random variables X and Y, is defined by joining their respective marginal CDFs F(x) and G(y) for any copula. In this setting, the measure R for dependent X and Y can be defined as where f (x) and g(y) denote the PDFs of X and Y, respectively, and c(u 1 , u 2 ) the copula density.
Using the FGM copula, we have the following relationship: and Now, we calculate R when X and Y have possibly non-identical GHLD with the CDFs , respectively. Hence, σ is common to the two marginal distributions. In this case, after some integral developments, we obtain

Estimation Method for the Distribution Parameters
In this section, we present the ML method for estimating the FGMBHLD parameters.
Let (x 1 , y 1 ) . . . (x n , y n ) be a random sample from a random vector (X, Y) following the FGMBHLD with the parameters µ 1 , µ 2 , σ 1 , σ 2 , and θ. Hence, in particular, X follows the GBHLD(µ 1 , σ 1 ) and Y follows the GBHLD(µ 2 , σ 2 ). Elaal and Jarwan (2017) [14] introduced the ML estimation method for bivariate distributions based on copula. The basis consists of constructing the log-likelihood function as where F(x) and G(y) are the CDFs of X and Y, and f (x) and g(y) are their respective PDFs, and c(u 1 , u 2 ) refers to the copula density. The ML estimates (MLEs) of the involved parameters are obtained by maximizing this function with respect to these parameters. Under the setting of the FGMBHLD, we have where φ(x i , µ 1 , σ 1 ) = 1 − 2F(x i ) and η(y i , µ 2 , σ 2 ) = 1 − 2G(y i ). The MLEs of the parameters µ 1 , µ 2 , σ 1 , σ 2 and θ, say µ 1 , µ 2 , σ 1 , σ 2 , and θ, are those maximizing this function. They can be obtained by differentiation. To be more precise, by differentiating the log-likelihood with respect to the distribution parameters, we obtain and ∂Ln L ∂θ where and By setting the above first partial derivatives of Ln L to zero, we obtain µ 1 , µ 2 , σ 1 , σ 2 and θ. Since we cannot obtain a closed form for these estimates, a numerical method must be used.

Estimation of the Stress-Strength Distribution Parameter
In this section, we introduce the MLE for R = P(Y < X). Moreover, we derive a motivated asymptotic confidence interval and a bootstrap confidence interval for it.

Maximum Likelihood Estimate of R
From observed data (x 1 , y 1 ) . . . (x n , y n ), which are taken from a random vector (X, Y) following the FGMBHLD with the parameters µ 1 , µ 2 , σ 1 , σ 2 , and θ, with σ = σ 1 = σ 2 , we consider the MLEs µ 1 , µ 2 , σ and θ of these parameters, respectively. Then, based on Equation (35) and the invariance property, the MLE of R is obtained by substitution as

Asymptotic Confidence Interval (ACI)
We now aim to compute the ACI for R with a large sample. Let Θ = (µ 1 , µ 2 , σ, θ), and Θ i be the i-th component of this vector. First, we construct the Fisher information matrix as follows: where ], i, j = 1, . . . , 4, Θ i refereing to the ith component of Θ. Second, we construct the variance-covariance matrix by replacing the distribution parameters by their MLEs, and we obtain where To obtain the ACI of R, the following two theorems are useful.
where A = V n .
Proof. The theorem can be demonstrated using the asymptotic properties of MLEs of the distribution parameters under regularity conditions and the multivariate central limit theorem.
Proof. The proof is based on Theorem 1 and the application of the delta method.
According to Xu and Long (2007) [15], a 100( where Z α/2 denotes the value providing an area of α 2 in the upper tail of the standard normal distribution, and B = b T A −1 b, where b is defined as Equation (51) with substitution of the unknown parameters by the corresponding MLEs.

Simulation
In this section, a Monte Carlo simulation study is introduced to describe the point and interval estimation of R.
The obtained n pair of values are thus generated values from (Y 1 , Y 2 ) following the FGMBHLD.
(ii) Re-sample the simple random sample (x i , y i ) with replacement. (iii) Obtain the new simple random sample (x * i , y * i ). Use the algorithm in Section 7.2 to generate different sample sizes with n = 30, 50, 70 and 100, with 10,000 replications. All computations are obtained using Mathematica 11.1.

3.
Calculate R MLE according to the methodology in Section 6.1 and the "average R MLE ", say R * MLE , based on all the samples at a fixed size.

4.
Evaluate the ACI and BCI according to the methodology in Sections 6.2 and 7.2.

5.
Study the behavior of R MLE by evaluating the bias defined by the "average of (R MLE − R)" and the mean square error (MSE) indicated as the "average of (R MLE − R) 2 ". 6.
In the context of interval estimation, we compare the ACI and BCI using the asymptotic confidence length (ACL) and converge probability (CP).
The results of the simulation study are presented in Table 1. In general, the length of the ACI becomes smaller than the length of the BCI.
The CP in almost all cases of the ACI is more than the CP in the BCI.
Hence, from the above results, the behavior of the MLEs is good for large samples. Moreover, the ACI is more suitable than the BCI for the stress-strength model.

Application: Household Financial Affordability in KSA 2018
In this section, we introduce a real application of the stress-strength model in an economic data setting, where X and Y represent household income and consumption, respectively. Here, R = P(Y < X) is a household's financial affordability. We use the data from the household income and expenditure survey of KSA 2018. The survey period was from 28 February 2017 to 31 March 2018 in each month. In this study, we are interested in studying the behavior of R when X represents the average household monthly income by administrative region for Saudi households and Y represents the average household monthly consumption expenditure by administrative region for Saudi households, in order to measure the financial affordability for Saudi households by administrative region in 2018. The data are shown in Table 2. Table 3 presents the descriptive statistics for the data.  To achieve our aim, we demonstrate the practicability of our proposed model. The Anderson-Darling (AD) goodness of fit statistic value is used to confirm that the GHLD is suitable for the income and consumption data; the corresponding p-values are almost equal to 1. Moreover, the quantile-quantile (Q-Q) plot is used to confirm this statement, as shown in Figure 2. Now, we evaluate R = P(Y < X) in the following two cases: Case 1: If X and Y are independent with X following the GHLD(µ 1 , σ) and Y following the GHLD(µ 2 , σ), and the dependent parameter θ is set as 0; Case 2: If X and Y are dependent with (X,Y) following the FGMBGHLD.
We calculate, in both cases, the MLEs of the distribution parameters and R, as well as the ACI and ACL. The results are shown in Table 4. 1.
Since θ is estimated as 0.4713, and is therefore positive, then the relation between X and Y is positive, as we see in Figure 3.

2.
The measure of affordability when X and Y are dependent is less than when X and Y are independent, so the case of dependent variables is more realistic. Finally, Figure 4 shows the (estimated) PDF and CDF of the FGMGBHLD with the estimated parameters from the considered data. It can be noted that the PDF seems unimodal (bump effect) with a long two-dimensional tail. With the FGMBGHLD, the equation behind the calculated PDF and CDF can be employed for further modeling.
To conclude this section, in order to show the performance of our new distribution on KSA data, we compare it with the bivariate Weibull distribution (BWD) as presented in Almetwally et al. (2020) [6]. First, we use the goodness of fit test and Q-Q plot to show that the BWD is a good fit to the KSA data. From the AD goodness of fit test, we find that the p-value equals 0.082 and 0.125 for the two considered data sets, respectively. As a result, the BWD fits the KSA data well. Figure 5 illustrates this conclusion. Now, we repeat our application but replace our proposed distribution by the BWD. Table 5 shows the result of the MLEs, R, ACIs, and ACLs of the distribution parameters and stress-strength model in the following two cases: Case 1: If X and Y are independent with X following the Weibull(α 1 , β) and Y following the Weibull(α 2 , β); Case 2: If X and Y are dependent with (X,Y) following the BWD. From the ACL viewpoint, we can compare the performance of our distribution and the BWD on the KSA data. Thus, from Tables 4 and 5, we observe that the ACLs for our proposed distribution are lower than those of the BWD for both cases.
We complete this result by using the AD test for copula-based distributions as described in Genest et al. (2013) [18]. Table 6 shows the p-values of this AD test for our distribution and the BWD (dependent case for both). Table 6. AD test for the proposed distribution and the BWD. The lower p-value is obtained for the FGMBGHLD distribution. Based on the results above, we can confirm that the proposed distribution is more suitable than the BWD for the considered KSA data.

Conclusions
In this paper, we introduced the bivariate distribution using the FGM copula approach, abbreviated as FGMBGHLD. We studied some of its statistical properties, such as the PDF, CDF, product moments, moment generating function, reliability function, and hazard rate function. In a multivariate statistical setting (and bivariate in particular), it is well known that the maximum likelihood estimation method gives unique estimates (under some regularity conditions) and guarantees their asymptotic performance from the unbiased and normality viewpoints. For these reasons, we developed it for the FGMBGHLD. We also applied the FGMBGHLD in a real-life data analysis scenario. We investigated the stressstrength model represented by R when the stress and strength variables are dependent and have the FGMBGHLD as a joint distribution. A simulated study was performed to study the behavior of the maximum likelihood estimate of R. Confidence intervals were constructed using two different techniques. Finally, we provided a real application of the considered (dependent) stress-strength model when X and Y measure the household financial affordability in KSA 2018 for Saudi households by administrative region. The obtained results are quite good and competitive with those of a valuable competitor (the bivariate Weibull distribution as introduced by [6]). Research perspectives include the application of the FGMBGHLD to more different bivariate data types, its multivariate version, and the development of regression model types.