Bell-generated family of mixture cure rate models with application to credit default

1 Introduction

Survival analysis is widely applied in finance, statistics, reliability engineering, and medical research to model the likelihood that a subject (or item) survives beyond a specified period of time [1]. A common assumption in standard survival analysis is that all subjects under study will eventually experience the event of interest if followed for a sufficiently long duration. However, in practice, this assumption is often violated, as some subjects may never experience the event, even under indefinite observation [2]. For instance, in credit risk, particularly in time-to-default analysis, a proportion of borrowers may never default. Similarly, in clinical trials, a fraction of patients may never experience the event of interest (e.g., disease recurrence) after receiving treatment. Such individuals are not merely right-censored in the traditional sense but are instead considered non-susceptible (or “cured”) with respect to the event of interest [3].

Traditional survival models do not directly account for the presence of a cured or non-susceptible subpopulation, thereby necessitating the development of cure rate models. Cure rate model, otherwise known as cure fraction model was first introduced by Boag [4] and subsequently extended by Berkson and Gage [5]. Since then, extensive studies on long-term survivor models have been conducted in the biological sciences (see [6, 7] and [8]). More recently, cure rate models have gained prominence in economics and finance literature, particularly for modeling the time-to-default of loans (see [3, 9], and [10]).

To some extent, cure rate models addresses the limitations of standard survival models by explicitly incorporating the fraction of individuals who are not at risk of experiencing the event in the modeling framework [11]. In credit risk applications, cure rate models can be used to analyze time-to-default, where a substantial proportion of observations are right-censored because some borrowers may never default during the lifetime of the loan. This type of data structure has been extensively studied in the literature [12].

Despite the growing interest in cure rate models for survival and credit risk analysis, most existing studies rely on classical baseline lifetime distributions or their conventional extensions, which often lack the flexibility required to adequately model highly skewed data, heavy-tailed distributions, and non-monotonic hazard functions commonly observed in time-to-default data (see [13, 14]). Thus, the development of distributions that provide a more accurate model time-to-event or failure time of events with and without the presence of competing risk factors is imperative. Such improved models can be found in Abdullahi et al. [15] and Abdullahi et al. [16]. Consequently, there remain the need for more flexible cure rate models that can simultaneously accommodate complex hazard behaviors, high censoring proportions, and heterogeneous borrower populations. The primary objective of this study is therefore to propose a new family of mixture cure rate models based on the Bell-generated (Bell-G) class of distributions for modeling time-to-default data. Specifically, we incorporate four widely used baseline distributions–Weibull, Burr XII, log-logistic, and lognormal–within the Bell-G framework to obtain flexible mixture cure rate models capable of representing a broad range of survival patterns. The Bell distribution [17] is particularly advantageous for handling overdispersion, which is a common feature in many real-world data. Additionally, the Bell distribution demonstrates robustness in the presence of excess zeros, making it well-suited for modeling data with high frequencies of zero counts [18]. In the context of cure rate models, the Bell-G contruction presents a more realistic inclusion of both homogeneous and heterogeneous populations by separating long-term survivors from short-terms ones. The proposed models are estimated using the maximum likelihood method, extended to incorporate covariate effects through a logistic link for the cure fraction, and their finite-sample performance is evaluated via Monte Carlo simulation. To the best of our knowledge, this is the first study to integrate the Bell-G family with mixture cure rate modeling using these classical baseline distributions. The practical utility of the proposed framework is further demonstrated through applications to two real-world credit default datasets, where the models are shown to provide improved goodness-of-fit compared with several existing cure rate models.

The remainder of this paper is organized as follows. Section 2 presents the methodology, including the formulation of the Bell-G mixture cure rate models, maximum likelihood estimation, some statistical properties and simulation studies. Section 3 illustrates the application of the proposed models to credit risk data. Finally, Section 4 concludes the paper and highlights its key contributions and limititions.

2 Methodology

2.1 Bell-G family of distributions

Since the time-to-occurrence of many events of interest is a random variable, it is essential to identify an appropriate probability distribution to model such occurrences. Numerous probability distributions have been proposed in the literature to describe survival times. However, many of these distributions lack the flexibility required to adequately model non-monotonic hazard functions, which are commonly observed in credit risk analysis. To address this limitation, we consider the compounded Bell-G class of distributions, primarily due to its flexibility and ability to handle overdispersion. The Bell-G class of distributions introduced by Alsadat et al. [19] offers notable advantages over well-known extensions such as the Poisson-G class of distributions proposed in Gomes et al. [20] and Maurya and Nadarajah [21]. In particular, the Bell-G distribution provides greater flexibility in accommodating various shapes of hazard functions, making it suitable for a wide range of survival data applications.

Definition. Let G(t) and g(t) denote the cumulative distribution function (CDF) and probablity density function (pdf), respectively of any continuous random variable, t. The CDF of the Bell-G distribution by Alsadat et al. [19] is given as:

F0(t)=1-exp(-eη[1-e-ηG(t)])-exp(1-eη)1-exp(1-eη),(1)

with survival function (SF):

S0(t)=exp(-eη[1-e-ηG(t)])-exp(1-eη)1-exp(1-eη),(2)

where t (time-to-event) > 0 and η (rate parameter) >0

The pdf can be derived as:

f0(t)=ηg(t)exp(η[1-G(t)])exp(-eη[1-e-ηG(t)])1-exp(1-eη),(3)

while the hazard rate function (hrf) is given by

h0(t)=ηg(t)exp(η[1-G(t)])exp(-eη[1-e-ηG(t)])exp(-eη[1-e-ηG(t)])-exp(1-eη).(4)

The Bell-G class of distribution has the following advantages over some existing generalized distributions [19]:

• Its original form contains the feature of Bell numbers;
• It adds additional parameter(s), thereby yielding a better fit as compared to the well-known Poisson-G;
• It provides better goodness-of-fit measures for highly skewed data.

2.2 Mixture cure rate models

Formally introduced by Berkson and Gage [5], the mixture cure rate model was further extended by Maller and Zhou [6] by simply adding the proportion of subjects that are not susceptible to the event of interest to the traditional survival function as:

Spop(t)=τ+(1-τ)S0(t),τ∈(0,1),(5)

where S₀(t) is defined in Equation 2 and τ is the proportion of cured subjects in the population.

The corresponding CDF, pdf and hrf are respectively, given by:

Fpop(t)=(1-τ)F0(t),(6)

fpop(t)=(1-τ)f0(t),(7)

and

hpop(t)=(1-τ)f0(t)τ+(1-τ)S0(t),(8)

where, F₀(t), S₀(t) and f₀(t) are the CDF, pdf and hrf of the susceptible group, defined in Equations 1–3, respectively.

2.3 Special submodels

To develop our proposed family of mixture cure rate models, we assume that the baseline survival function, S₀(t) in Equation 5 follows the proposed Bell-Weibull, Bell-Burr XII, Bell-Loglogistic and Bell-Lognormal distributions.

2.3.1 The Bell Weibull mixture cure rate (BWMcr) model

Given the CDF and pdf of the two-parameter Weibull distribution respectively as:

Gw(t)=1-exp(- (tξ)ω),t>0,(9)

gw(t)=ωξ(tξ)ω-1exp(- (tξ)ω),t>0,(10)

where, ω>0 and ξ>0 are the shape and scale parameters respectively [22]. For specification, G(t) in Equation 1 is same as G_*w*(t) in Equation (9) and g(t) in Equation 3 is the same as g_*w*(t) in Equation 10.

We then, derive the SF, CDF, pdf and hrf of the BWMcr model by substituting Equations 9, 10 into Equations 1–4, respectfully as:

Sbw(t)=τ +(1-τ)(exp(-eη[1-e-η(1-Aξω)])-exp(1-eη))1-exp(1-eη), t>0,(11)

Fbw(t)=(1-τ)(1-exp(-eη[1-e-η(1-Aξω)]))1-exp(1-eη),t>0,(12)

fbw(t)=ηωξ(tξ)ω-1Aξω(1-τ)exp(ηAξω)exp(-eη[1-e-η(1-Aξω)])1-exp(1-eη),t>0,(13)

and

hbw(t)=ηωξ(tξ)ω−1Aξω(1−τ)exp(ηAξω)exp(−eη[1−e−η(1−Aξω)])τ(1−exp(1−eη))+(1−τ)(exp(−eη[1−e−η(1−Aξω)])(−exp(1−eη)),t>0,(14)

where Aξω=exp(-(tξ)ω) and 0 < τ < 1.

At different parameter values as indicated in the plot, the survival curves as shown in Figure 1a plateaus above zero, while the CDF curve exhibits a reverse case in Figure 1b indicating the presense of a cured proportion. The pdf of BWMcr exhibit right-skewed, left-skewed and approximately symmetric shapes (see Figure 1c). The hazard rate function also displays a variety of shapes including increasing, decreasing and modified upside-down bathtub as shown in Figure 1d. This indicates the model's ability to accomodate heterogeneous and time-varying risk dynamics.

2.3.2 The Bell Burr XII mixture cure rate (BBXIIMcr) model

Given the CDF and pdf of the three parameter Burr XII distribution as

Gbxii(t)=1-(1+(tξ)ω)-φ,t>0,(15)

and

gbxii(t)=ωφξ-ωtω-1(1+(tξ)ω)-φ-1,t>0,(16)

where, ω>0, φ>0 are the shape parameters and ξ>0 is the scale parameter [23]. Here, G(t) in Equation 1 is same as G_*bxii*(t) in Equation (15) while g(t) in Equation 3 is the same as g_*bxii*(t) in Equation 16.

Now, we derive the SF, CDF, pdf and hrf of the BBXIIMcr model by substituting Equations 15, 16 into Equations 1–4, respectfully as:

Sbbxii(t)=τ +(1-τ)(exp(-eη[1-e-η(1-Bξω-φ)])-exp(1-eη))1-exp(1-eη),t>0,(17)

where Bξω=1+(tξ)ω,

Fbbxii(t)=(1-τ)(1-exp(-eη[1-e-η(1-Bξω-φ)]))1-exp(1-eη),t>0,(18)

fbbxii(t)=ηωφξ-ωtω-1Bξω-φ-1(1-τ)exp(ηBξω-φ)exp(-eη[1-e-η(1-Bξω-φ)])1-exp(1-eη),t>0,(19)

hbbxii(t)=ηωφξ−ωtω−1Bξω−φ−1(1−τ)exp(ηBξω−φ)exp(−eη[1−e−η(1−Bξω−φ)])τ(1−exp(1−eη))+(1−τ)(exp(−eη[1−e−η(1−Bξω−φ)]−exp(1−eη)),t>0.(20)

For some parameter values indicated, Figures 2a, b depict different potential survival function and CDF shapes for the BBXIIMcr which confirms its validity as a cure rate model. Figure 2c illustrates skewed and nearly symmetric shapes of the pdf. The hazard rate function depicts monotone and non-monotone shapes of the BBBXIIMcr model showing the model's ability to mimic a variety of datasets (see Figure 2d).

2.3.3 The Bell loglogistic mixture cure rate (BLLMcr) model

Given the CDF and pdf of the standard Log-logistic distribution as

Gll(t)=tωtω+ξω,t>0,(21)

and

gll(t)=(ωξ)(tξ)ω-1[1+(tξ)ω]2,t>0,(22)

where, G(t) in Equation 1 is equal to G_*ll*(t) in Equation 21, g(t) is substituted as g_*ll*(t) in Equation 22, ω>0, ξ>0 are the shape and scale parameters respectively [24]. The derivation of the SF, CDF, pdf, and hrf of the BLLMcr model is done by simple substitution of Equations 21, 22 into Equations 1–4, respectfully as

Sbll(t)=τ+(1-τ)(exp(-eη[1-e-ηtωtω+ξω])-exp(1-eη))1-exp(1-eη),t>0,(23)

Fbll(t)=(1-τ)(1-exp(-eη[1-e-ηtωtω+ξω]))1-exp(1-eη),t>0,(24)

fbll(t)=(1-τ)η(ωξ)(tξ)ω-1exp(η[1-(tωtω+ξω)])exp(-eη[1-e-ηtωtω+ξω])[1+(tξ)ω]2(1-exp(1-eη)),t>0,(25)

hbll(t)=(1−τ)η(ωξ)(tξ)ω−1exp(η[1−tωtω+ξω])exp(−eη[1−e−ηtωtω+ξω])[1+(tξ)ω)2[τ[1−exp(1−eη])+(1−τ))(exp(e−η[1−e−ηtωtω+ξω])−exp(1−eη))],t>0.(26)

The plots of BLLMcr model is given in Figure 3. Here, survival and CDF plots in Figures 3a and b, respectively indicates strong evidence of the presence of non-susceptible proportion of subjects in the model at different parameter values. The pdf shows skewed and almost symmetric shapes whereas the hrf depicts increasing, decreasing and modified upside-down bathub shapes as shown in Figures 3c and d, respectively.

2.3.4 The Bell lognormal mixture cure rate (BLNMcr) model

Given the CDF and pdf of the standard Lognormal distribution as

Gln(t)=12[1+erf(ln t-μσ2)],t>0,(27)

and

gln(t)=1tσ2πexp(-(ln t-μ)22σ2),t>0,(28)

where, G(t) in Equation 1 is equal to G_*ln*(t) in Equation 27 and g(t) in Equation 3 is also equal to g_*ln*(t) in Equation 28, μ∈(−∞, +∞) and σ>0 are the mean and standard deviation, respectively.

Hence, the SF, CDF, pdf and hrf of the proposed BLNMcr model is derived by substituting Equations 27, 28 into Equations 1–4, respectfully as

Sbln(t)=τ+(1-τ)(exp(-eη[1-e-η2Lμσ])-exp(1-eη))1-exp(1-eη), t>0,(29)

Fbln(t)=(1-τ)(1-exp(-eη[1-e-η2Lμσ]))1-exp(1-eη),t>0,(30)

fbln(t)=(1-τ)ηexp(-(ln t-μ)22σ2)exp(η[1-Lμσ2])exp(-eη[1-e-ηLμσ2])tσ2π(1-exp(1-eη)), t>0(31)

hbln(t)=(1−τ)ηexp(−(lnt−μ)22σ2)exp(η[1−Lμ,σ2])exp(−eη[1−e−ηLμ,σ2])tσ2π[τ(1−exp(1−eη))+(1−τ))(exp(−eη[1−e−ηLμ,σ2])−exp(1−eη))],t>0,(32)

where Lμσ=[1+erf(ln t-μσ2)].

At some parameter values as indicated, plots in Figures 4a and b, respectively show different potential valid cure rate survival function and CDF shapes for the BLNMcr model. Figures 4c, d also demonstrate different pdf and hazard rate functional shapes for the BLNMcr model. These shapes is an indication of the ability of the model to deal with different characteristics of different datasets.

2.4 Inference

In this section, we present some statistical properties of the our proposed models, a systematic procedure to obtain estimates by the maximum Likelihood Estimation (MLE) method and perform simulation studies.

2.4.1 Quantile function

We obtain the quantile function (QF) for the BWMcr which is useful for the generation of random numbers from the model. This is derived as the inverse of CDF of the BWMcr model given in Equation 12 as

Qu=ξ[−log(1 +1ηlog(1+1eηlog(1−u(1−exp(1−eη))1−τ)))]1ω, u∈(0,1).(33)

Proof. By definition, the QF is given as Qu=Fbw-1(t).

Let, F_*bw*(t) = u, thus,

(1-τ)(1-exp(-eη[1-e-η(1-exp(-(tξ))ω)]))1-exp(1-eη)=u(1-τ)(1-exp(-eη[1-e-η(1-exp(-(tξ))ω)]))=u(1-exp(1-eη))

For simplification, let L = (1 − exp(1−e^η)) and A=exp(-(tξ)ω) so that

uL=(1-τ)(1-exp(-eη[1-e-η(1-A)]))

dividing through by (1−τ) gives

uL1-τ=1-exp(-eη[1-e-η(1-A)])⇒exp(-eη[1-e-η(1-A)]) =1-uL1-τ

Taking the logarithm of both sides of the equation and dividing by −e^η gives

1-e-η(1-A)=-e-ηlog(1-uL1-τ)⇒e-η(1-A)=1+e-ηlog(1-uL1-τ)

Taking the logarithm of both sides again and dividing by −η gives

A=1+1ηlog[1+e−ηlog(1−u1−τ)]

substituting L and A back into the equation to obtain:

exp(−(tξ)ω)=1 +1ηlog[1+e−ηlog(1−u(1−exp(1−eη))1−τ)]

Taking the logarithm of both sides and making “t” the subject gives the QF of the BWMcr model as

Qu=ξ[−log(1+1ηlog(1 +1eηlog(1−u(1−exp(1−eη))1−τ)))]1ω, u∈(0,1).

2.4.2 Tail behavior

The tail behavior of the mixture cure rate survival model describes how the probability of surviving decays as time, t increases beyond bounds.

2.4.2.1 Proposition 1

For the proposed Bell-G mixture cure rate model,

limtarrow∞S(t)=τ.

Proof. Substituting S₀(t) in Equation 2 into Equation 5, we can re-write the survival function of the mixture cure rate model as

S(t)=τ+(1-τ)[exp(-eη[1-e-ηG(t)])-exp(1-eη)1-exp(1-eη)](34)

Since G(t) is a proper CDF,

limtarrow∞G(t)=1(35)

Substituting G(t) = 1 into Equation 34 gives

S(t)=τ+(1-τ)[exp(1-eη)-exp(1-eη)1-exp(1-eη)]

which simplifies to

S(t)=τ+(1-τ)(0)

Hence, limtarrow∞S(t)=τ.

2.4.3 Identifiability

Identifiability is an important property of statistical models that ensures that each set of parameter values corresponds to a unique distribution. This implies that if two parameters are different, then their corresponding distribution must also be different. A model is said to be generally identifiable if the set of non-identifiable parameters is of measure zero within the parameter space [25].

2.4.3.1 Proposition 2

The pdf of the Bell-G distribution {f₀(t; η):η>0} is identifiable.

Proof. Given the pdf of the Bell-G family in Equation 3.

Suppose we have two different parameter values of η₁ and η₂, thus, η₁≠η₂.

We consider the ratio,

R(t)=f0(t;η1)f0(t;η2).

Then

R(t)=P⋅exp{(η1−η2)(1−G(t))−eη1(1−e−η1G(t)) +eη2(1−e−η2G(t))}

where

P=η1[1-e1-η2]η2[1-e1-η1]

Now, if η₁≠η₂, the exponents depends on G(t), so R(t) cannot equal 1 for all t.

Hence,

f0(t;η1)≠f0(t;η2)

This proves that the model is identifiable.

2.4.4 Maximum likelihood estimation

Given a random sample size, n with observed pairs (t_*i*, δ_*i*), where δ_*i* represents the censoring indicator on the interval (0, 1) for censored and uncensored observations, respectively for i = 1, 2, …, n and t_*i* denoting the observed default times. The log-likelihood function for the mixture cure rate model is given by

L(θ)=∑i=1nδilog[(1-τ)f0(ti)]+∑i=1n(1-δi)log[τ+(1-τ)S0(ti)](36)

To obtain the log-likelihood functions of the BWMcr, BBXIIMcr, BLLMcr and BLNMcr models, Equations 11, 13, 17, 19, 23, 25, 29, 31 are substituted into Equation 36. Subsequently, the estimates of the parameters of the models are obtained by direct maximization of the log-likelihood function given in Equation 36.

2.4.5 Simulation study

Simulation studies were conducted to evaluate the finite-sample performance and statistical properties of the proposed models. The maximum likelihood estimates via Equation 36 were obtained for the family of mixture cure rate models following the steps below:

• a. Derive the quantile function Q_*u* for one member of family; BWMcr.
• b. Generate a sample of size n of u₁, u₂, ..., u_*n* from U~(0, 1).
• c. Given that τ_*i* is the cure rate parameter,

ti′={Qu,ui≤1-τi∞,otherwise

• d. Generate censored samples c₁, c₂, ..., c_*n* from the exponential distribution with a rate equal to 0.2.
• e. Obtain the right-censored samples data as ti=min(ti′,ci).
• f. Obtain the pairs (t_*i*, δ_*i*), where δ_*i* = 1 if t_*i* ≤ c_*i* and 0, otherwise.
• g. Obtain the maximum likelihood estimates of the parameters using Equation 36 and obtain the average estimates (AE), absolute bias (AB) and the root mean square error (RMSE) of the estimates.
• h. Repeat steps (b)–(f) 1,000 times for sample sizes n = 50, 100, 300, 500, 800.
• i. Repeat steps (b)–(g) for parameter set τ, η, ξ, ω = (0.1, 1.2, 0.3, 0.8) at τ_*i*= 10% and 8%, respecively.

The results of the simulation study for the BWMcr model are given in Table 1. From the results presented, it can be observed that all the parameters are consistent as the AB and RMSE decreases as sample size increases for both 10 and 8 cure percentages considered, respectively. Also, the AE values of the parameters get closer to the true parameter values as sample size increases.

Table 1. Simulation results of BWMcr model for (τ, η, ξ, ω)=(0.1, 1.2, 0.3, 0.8).

2.5 Regression

We now present regression models with cure rates. In this case, the survival time-to-event of loan contract is considered as the response variable which is assumed to follow the our proposed distributions. The essence of regression models in this study is to measure the possible effects of covariates on the time-to-default of loans issued to customers. Using the logit function, we relate the cure parameter τ to the covariates given as

τi=eXiTβ1+eXiTβ,i=1,2,…n,(37)

where XiT=(1,xi1,xi2,…,xik) is the ith vector of covariates and β=(β0,β1,β2,…,βk)T is the vetcor of regression coefficients. The structure of the regression model is subsequently given as

log(τi1-τi)=XiTβ,i=1,2,…,n.(38)

The log-likelihood functions for the proposed mixture cure rate regression models are obtained by substituting τ_*i* into Equation 36. The maximum likelihood method is used to obtain the estimates of the regression parameters.

2.6 Randomized quantile residuals (RQRs)

In this study, the randomized quantile residual (RQR) is used to examine the best fit of the cure rate models. For a model to provide best fit to a given data, the RQRs of the model are expected to assume a standard normal distribution. Thus, it is mathematically given as:

ri=Qu(Ft(ti;β^)),(39)

where Q_*u* is the quantile function of the standard normal distribution and β^ is the vector of the estimated parameters of the cure rate models [26].

3 Application

We now apply the proposed models to a time-to-default loan dataset, which includes loan survival times, interest rates, loan amounts, borrower income, and other related variables. All computations were done with R software using the Broyden-Fletcher-Goldfarb-Shanno optimization algorithm with bounds (“L-BFGS-B”) as the maximization algorithm. We adopted this optmizer due to its ability to produce quicker convergence with more stable estimates in medium size models [27]. To ensure that our optimation reached global maximum, we used multiple random initial values that produced same log-likehood values and parameter estimates.

3.1 Time-to-default loan datasets

Loans issued for the purposes of medical care and credit card were extracted from the XYZ lending dataset [28]. The dataset spans from loan contracts issued from January 2011 to December 2015 for funds equal or above $5,000 with repayment period of 60 months.

A total of 1,476 and 53,820 loan contracts were issued for medical care and credit card purposes, respectively. Table 2 presents the summary statistics for the two loan datasets, showing that about 91% and 95% of the observations are censored for medical care and credit card loans, respectively.

Table 2. Summary statistics of time-to-default of loan datasets (in months).

As depicted in Table 2, the two datasets exhibit right-skewness and heavy tails. However, the credit card loan dataset shows slightly more positive assymetry and higher kurtosis as compared to the medical loan data. The evidence of the existence of “non-subjects” of interest is illustrated in Figure 5, where the Kaplan–Meier (K–M) curve reaches a plateau well for both datasets.

3.2 Model selection criteria

The measure for model selection used in this study includes the Akaike information criterion (AIC), Bayesian Information Criterion (BIC) and visual comparison of the fitted survival curves vs. Kaplan–Meier curves. We compared the performance of the proposed family of mixture cure rate models to some well-establsihed models namely; Odds Burr XII-Logistic Mixture Cure Rate (OBXIILMcr) by Hassan et al. [29], Weibull Exponentiated Exponential Mixture Cure Rate (WEEMcr) by Usman et al. [30], Log-exponentiated-Weibull Mixture Cure Rate (LEWMcr) by Cancho et al. [31] and the Exponentiated Weibull Exponential Mixture Cure Rate (EWEMcr) by Omer et al. [32]. The mixture cure rate models proposed in this paper appeared to be relatively superior since it had the least of values of the measurement criteria, presented in Tables 3, 4 for medical and credit card loans, respectively. It is evident that our proposed models produced the lowest criteria values, more specifically the BLNMcr presented a reduction of 4.84 units in AIC as compared to the best competitor model; LEWMcr model. Addionally, when we compare the other three propsed models to that of the other three competitos models, there is reduction in AIC (>10 units) for the medical loan dataset for the medical loan data. This inidcates a strong evidence of improved fit presented by our proposed models. For credit card loan data presented in Table 4, our poposed models presents a reduction in AIC (> 10 units) comparatively on the average. Specifically, the BWMcr models gave the lowest criteria values.

Table 3. Model selection criteria for medical loan data.

Table 4. Model selection criteria for credit card loan data.

3.3 Mixture cure rate models without covariates

The usefulness of the proposed cure rate models are demostrated in this section. All the four mixture cure rate models are fitted to the data as well as the other four competitor mixture cure rate models. Tables 5, 6 present parameter estimates of the fitted models to the medical and credit card loan datasets, respectively with their corresponding standard errors and p values using the maximum likelihood estimation (MLE) technique.

Table 5. Parameter estimates for medical care loan data.

Table 6. Parameter estimates for credit card loan data.

To illustrate the models fit performances, the fitted survival curves of the models are obtained and compared with the K–M survival curve. The proposed models and the four competitor models are overlaid with the K–M survival cuves. Figure 6 shows that the four sub models developed from the family are more closer to the K–M curve than the other four competitor models in the case of the medical loan dataset. Hence, the proposed models demonstrates improved fit to the medical loan dataset. This results can be confirmed in Table 3 with the BLNMcr model emerging relatively better.

As shown in Figure 7, the graphical comparison shows that all four of the proposed models lies closer to the K–M curve.

3.4 Mixture cure rate models with covariates

The applicability of the Bell-G familly of mixture cure rate regression model developed in this study is also illustrated in this section. The data used for the demonstration is the XYZ cooperation loan dataset on loans issued between 2011 to 2015 for the purpose of medical care. The dependent variable is defined as the survival times-to-default while the covariates considered are interest rate, loan amount and annual income of borrower. Thus, the regression model used takes the form:

log(τi1−τi)=β0+β1xi1+β2xi2+β3xi3, i=1,2,…1,476,(40)

where x_*i*1, x_*i*2 and x_*i*3 represents interest rate, loan amount and annual income of the ith borrower respectively. Mixture cure rate regression models for the BWMcr, BBXIIMcr, BLLMcr and BLNMcr are fitted to the data. The performance of these regression models are compared with competitor regression models such OBXIILMcr by Hassan et al. [29], WEEMcr by Usman et al. [30], Sine Weibull Geometric Mixture Cure Rate (SWGMMcr) by Angbing et al. [33] and the Cosine-Frêchet Mixture Cure Rate (CosFMcr) by Tee et al. [10]. Table 7 gives the summary statistics of medical loan data with covariates, while Table 8 show the parameter estimates with their corresponding standard errors and p-values for all the fitted regression models.

Table 7. Summary statistics of the medical loan data with covariates.

Table 8. Parameter estimates for regression models.

It can be observed that the covariates; loan amount and annual income of borrower are not statistically significant for all the regression models fitted at the 5% significance level, with loan amount being highly insignificant with a p-value greater than 0.14. However, interest rate was found to be statistically significant for all the regression models suggesting that the likelihood of non-default in the medical care loan is greatly influenced by the interest rate. This is because interest rate acts as a composite borrower risk signal that is usually set by lenders based on loan amount requested, annual income of borrower, and overall creditworthiness of the borrower. Hence, the marginal explanatory power of loan amount and annual income is absorbed by the interest rate in the regression model. This is consistent with works by Abdymomunov et al. [34] and Jayaram et al. [35]. Again, the results also indicated that all parameter coefficients for interest rate and loan amount were negative, indicating that a one unit increase in interest rate and loan amount will cause a reduction in the log odds of the probability of non-default. Meanwhile, the parameter coeficient for annual income of borrower were postive for all regression models which means that a unit rise in annual income of borrowers will result in a unit increase in the log odds of the probability of non-default.

Table 9 shows the goodness-of-fit measures of the fitted regression models. It is evident that our proposed models demostrated better goodness-of-fit than all the other competing regression models as they have the least of the information criteria values. Overall, the BLNMcr regression model provides the most accurate description of the medical loan data since it recorded the lowest AIC and BIC values.

Table 9. Goodness-of-fit measures of the regression models.

Probability–probability (P–P) plots of the RQRs for the fitted regression models are shown in Figure 8. It can be observed that the plotted points for our proposed regression models are closer to the diagonal line as compared to the competitor regression models. This demonstrates that the emperical probabilities are in close agreement with the theoretical probabilities, which suggests that the proposed cure rate regression models provides a good fit to the medical loan data.

4 Conclusion

In this paper, we propose an alternative to existing survival models with a cured fraction based on the Bell-G class of distributions. Key properties of the proposed family, including the quantile function and the analytical forms of the survival, cumulative distribution, probability density and hazard rate functions were derived. Four special cases were presented: the Bell-Weibull, Bell-Burr XII, Bell-Log-Logistic, and Bell-Lognormal mixture cure rate models. A simulation study was conducted to evaluate the performance of the estimators which produced satisfactory results. The applicability of the proposed models was demonstrated using two loan datasets and compared with other well-known mixture cure rate models. The findings revealed that our proposed family of mixture cure rate models demonstrates improved fit to the datasets. Specifically, the BLNMcr and BWMcr model were found to be the most accuate in modeling medical and credit card loan datasets repectively. Furthermore, regression models were developed using interest rate, loan amount and income of borrower as covariates by the logistic link function. The BBXIIMcr model was found to best performing cure rate regression model where only interest rate was the only statistically significant covariate that could affect the time-to-default. The results indicate that the proposed family provides a competitive and flexible alternative in modeling survival data with a cured fraction as evidenced by the model selection criteria considered. This improvement is attributed to the additional flexibility introduced by the Bell-G construction, which allows the models to capture complex distributional characteristics such as skewness, heavy tails and non-monotonic hazard rate patterns often observed in credit risk data, contributing meaningfully to the existing literature. These results imply that lenders can better identify long-term non-defaulters, leading to improved credit scoring and risk-adjusted pricing. Despite these promising results, some limitations should be noted. First, Bell-G cure rate model introduces additional parameters to the standadard cure rate model, thereby increasing computational complexities. Second, parameter estimation was carried out using the maximum likelihood method, which may encounter numerical difficulties in highly censored datasets. Third, the empirical validation relied on only two datasets; therefore, additional applications across different domains would further establish the robustness and generalizability of the proposed models. Future research may extend the proposed Bell-generated mixture cure rate models by incorporating other covariates such as age and gender of borrowers, exploring Bayesian estimation procedures, and developing diagnostic tools for assessing model adequacy. Such extensions would further enhance the applicability of the proposed framework in credit risk modeling, reliability analysis, and biomedical survival studies where long-term survivors or cured fractions are commonly observed.

Generative AI statement

The author(s) declared that Generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Footnotes

Edited by: Aliyu Ismail Ishaq, Ahmadu Bello University, Nigeria

Reviewed by: Usman Aliyu Abdullahi, Kano State University of Technology, Nigeria

Lazhar Benkhelifa, University of Biskra, Algeria