On the Maintenance Modeling of a Hybrid Model with Exponential Repair Efficiency

In the reliability literature, maintenance efficiency is usually dealt with as a fixed value. Since repairable systems are subject to different degrees and types of repair, it is more convenient to regard a random variable for maintenance efficiency. This paper is devoted to the statistical study of a general hybrid model for repairable systems working under imperfect maintenance. For both failure improvement and virtual age reduction of the system, maintenance efficiency is assumed to be random, with an exponential distribution as a probability density function. The likelihood function of this model is provided, and the estimation of the model parameters is computed by considering the maximization likelihood procedure. Obtained results were tested and applied to simulated and real data sets. To construct confidence intervals, the bias-corrected accelerated bootstrap method has been used. KeywordsReliability, Failure rate, Virtual age, Likelihood estimation, Parametric statistical test.


Introduction
Modeling the maintenance efficiency on a repairable system has been widely studied in the literature, according to this efficiency, models for the failure process are classified as perfect, minimal and imperfect maintenance. Perfect maintenance establishes the system to an as good as new state and the failure intensity restarts at zero. The minimal maintenance brings back the system on when a failure occurs without changing the failure intensity or improving the state of the system. Imperfect maintenance is between these two extreme cases which is more realistic in practical engineering.
In most imperfect models, the maintenance efficiency is represented by a parameter. Evaluating the efficiency then amounts to estimating this parameter. Many models have been proposed. Brown and Proschan (1983) assumed that the system is perfectly repaired with probability and minimally repaired with probability 1 − , while Block et al. (1985) proposed a model where the probability of maintenance is conditioned by the failure times. The most known models for imperfect maintenance are virtual age models, which are characterized by a set of random variables { } ≥1 called reduction factor in the virtual age and represent the maintenance efficiency. Kijima proposed and discussed Kijima's type 1 and type 2 models (Kijima et al., 1988;Kijima, 1989), which are both practical in different kinds of system and maintenance activities. Another known class of imperfect maintenance model are the failure rate models. Nakagawa (1986Nakagawa ( , 1988 assumes that the preventive maintenance effect is expressed by a change in the slope of the failure intensity. Chan and Shaw (1993) quantified the maintenance efficiency by a reduction in the failure rate proportional to the value of the failure intensity. Doyen and Gaudoin (2004) considered an arithmetic reduction of intensity or reduction of age to express the effect of the maintenance. Later Doyen et al. (2017) studied the geometric reduction of intensity or reduction of age. Ghnimi et al. (2017) analyzed the memory effect order of the arithmetic reduction for age classes using likelihood ratio test. Zantek et al. (2015) proposed a stochastic process, where virtual age of the system is reduced by preventive maintenance (PM). Whereas, corrective maintenance (CM) shifts the failure intensity parameters. Lin et al. (2000) presented an imperfect PM model by combining age reduction PM model and failure rate adjustment PM model. Remy et al. (2013) presented a case study where they considered several probabilistic models by taking the effect of both PM and CM on the system age. Fuqing and Kumar (2012) generalized Kijima type 1 and type 2 by expressing the maintenance efficiency with a function depending on time . Random variables were incorporated in failure rate and virtual age models by Wu and Clements-Croome (2005). Peng et al. (2015) optimized the PM policy in a sequential hybrid PM model addressing random factors for adjustment in intensity and reduction in age. Later, Duan et al. (2018) investigated a selective maintenance problem considering a multi-component system with stochastic maintenance efficiency. Gasmi (2019) introduced statistical inference for a system with a uniform distributed degree of repair and developed simultaneous confidence regions for the parameters.
In reliability literature, there exists a lot of work that focuses on maintenance policies. This work considers known parameters for the failure intensity of systems, which is not necessarily the case in practice, statistical studies are then needed to estimate the model parameters. However, estimating the maintenance efficiency is a primordial and crucial step in modeling repairable systems, especially when these latter are subject to imperfect maintenance with different degrees of repair. Therefore, the development of theoretical models, as well as the conception of empirical methods, are necessary in conducting reliability analysis. This paper is intended for the study of a general hybrid reliability model, which combines the failure rate and the virtual age classes. The discussed model considers the random feature for maintenance efficiency.

Imperfect Repair Model
The failure rate function is considered from the literature (Peng et al., 2015) and the corresponding failure intensity ( ) after the ℎ maintenance is expressed by: where is the improvement factor and represents the reduction in the virtual age of the ℎ maintenance activity. These two quantities are random with distribution functions and respectively. is the time between ( − 1) ℎ and ℎ maintenance activities. The recursive relation of the failure intensity function ( ) can be expressed for ∈ {2, … , } as: When = 1, we have 1 ( ) = ( ), by considering the ( , ) distribution for the initial failure intensity, we obtain: where > 0 is the shape parameter and > 1 is the scale parameter, thus the initial reliability function of the system equals: As maintenance is supposed to reduce the failure intensity in order to increase system reliability, it would be reasonable and more convenient to consider a density function for which 0 < ∫ ( ) +∞ 0 < 1. Also, reduction in the virtual age of the system would as well satisfy 0 < ∫ ( ) +∞ 0 < 1.
As is not an integer and by considering ! = Γ( ), where Γ( ) is the Gamma function defined by Γ( ) = ∫ −1 +∞ 0 exp(− ) , we obtain: the reliability function of the system for = 2 to corresponds to:

Parameter Estimation
Let { 1 , ⋯ , } be a random sample of size of the Exponential Maintenance Efficiency ( ) process, where = ( , , ). The likelihood function of the model is given by: The log-likelihood function can be written as: The maximum likelihood estimators (MLEs) of can be obtained by maximizing the log-likelihood function. Usually, the MLEs of are denoted by ̂ which corresponds to the solution of the system ∇ ln(ℒ( , )) = 0 3 .
Since there is no closed analytical method to solve the resulting nonlinear system, numerical methods using an iterative process are applied to find unknown parameter estimators. Such as the trust region approach (Byrd et al., 2000). We can as well transform the problem into an optimization problem and apply some constraint optimization method as the interior point (Coleman and Li, 1996;Vanderbei and Shanno, 1999). Thus, estimates of the model parameters are obtained.
To approximate confidence intervals, standard errors of the model parameters can be computed by corrected-bias accelerated (BCa) method (Efron, 1987). First, estimators are considered as true values for the parameters, then, a large number of samples is simulated. Next, the MLEs of each bootstrap sample are determined, which we denote by 1 intervals, calculation of the bias-correction parameter 0 and the acceleration parameter is required. A bias-correction 0 is computed as: Where −1 is the inverse of the standard normal cumulative distribution and is the Bootstrap cumulative distribution function of * . Let − * denote the estimate of when we omit the ℎ observation from the original data and let + represent the mean of these partial estimates. Using the Jackknife method, the acceleration parameter is estimated by: By denoting α 2 the upper α 2 -percentage of standard normal distribution, the BCa confidence intervals of the model parameters are then constructed by:

Simulation Study
To test the performance of the estimation method, a simulation study was conducted. were generated from the inverse of the cumulative density function (cdf) evaluated at a uniform distribution on [0,1], as follows: For = 1, the corresponding cdf is given by: The time 1 is simulated by solving the equation -= 0, where is a uniform random number, thus: For = 2 to , we have: ))} , The associated mean squared error (MSE) of the estimators were estimated. The MSE is defined by: By computing the estimators {θ } ∈{1,…, } of simulations, the mathematical expectation (θ) and the variance (θ) were estimated by their unbiased estimators (θ) and 2 (θ) as follows: The criterion MSE is then estimated by: The performance of estimation the model has been assessed by the mean of a simulation. = 100 samples of size = {50, 100, 200, 300, 400, 500} were generated by taking = 4, = 2.5, = 1.5. Obtained results are resumed in Table 1.  Figure 1. Parameter estimation and of the EME model. Figure 1 shows the simulation results of the EME model. Small crosses represent the obtained estimates based on simulated sample sizes, while true values of the parameters are illustrated by red segments. Simulation results showed that by increasing sample size, the decreases and parameter estimates of the EME model are close to their true values. This indicates the MLEs consistency properties.
Coverage probability of the 95% Wald confidence interval was also estimated. The Wald confidence interval is calculated by [̂± 1.96√ (̂)].
The Wald confidence interval proportions (associated with each simulated sample) that include true values of the parameters approximate a nominal value of 95% as desired. These proportions vary between 94% and 96%, which denotes the robustness of the proposed estimation procedure.
In the second part of the simulation study, method performance was assessed. Thus, datasets were simulated under different cases as described in Table 2. Several values have been considered for fixed maintenance efficiency . The parameter α was randomly selected, while β is equal to the scale parameter which characterizes aging. The parameter is chosen such that the mean of the exponential distribution 1 approximates c. These datasets were used to test the null hypothesis of non-random efficiency. When maintenance efficiency is constant (CME) and equals c the corresponding ℎ failure intensity and reliability function are respectively given by: Let = 1000 samples of = {200, 500} maintenance times were generated. The proposed estimation procedure is applied to each of these datasets. Then, the null hypothesis of non-random efficiency is tested, i.e. 0 : maintenance efficiency is constant, against 1 : maintenance efficiency is random.
The proportion of times that and of the model was less than the model is calculated. and are defined by: where, is the number of model parameters. Here equals 3 and is the generated sample size.
The selected model has the smallest value of these criterions. Identical proportions for both and were obtained. These proportions are reported in the last column of Table 2. Figure 2 illustrates the obtained proportions of AIC according to different cases of Table 2. From an overview, and according to the proportions of the AIC, it is noticeable that cases 2, 6 and 10, have the same behavior as cases 4, 8 and 12. Whereas, cases 3, 7 and 11 tend to behave as cases 1, 5 and 9.  In cases where datasets were generated under the null hypothesis, i.e., cases 1, 5 and 9, the estimated proportion of that favored the model was lower than the model as intended. This proportion tends to zero by increasing . When equals 1.4 in case number 11, we have a proportion of 0.19 for = 500. This indicates that the alternative hypothesis still holds, but the decrease in the proportion has just slowed down. For = 1000, we got = 0.15, it indicates that the values taken by affect the speed with which the proportion of increases or decreases. This can be explained by calculating the variance of the density of maintenance efficiency, for = 1.4 the corresponding variance equals 0.51 which is much greater when = 2 and = 3. This gives less harmonious values for the reductions and improvements, which creates a more dynamic process to the system. For cases 4, 8 and 12, datasets were generated under the alternative hypothesis ( ). The estimated proportion of for the model was higher than the model as desired. In cases where only factor improvement is random, prefers the null hypothesis than the alternative one. For these cases, the model distribution function is given by: .
By taking μ = 1 , we obtain the distribution function of the . This behavior can be justified by the strong law of large numbers. Thus, the average reduction in the intensity converges to their theoretical expectation which equals

Application
For illustration, a well-known automobile data studied by Ahn et al. (1998) and Guida and Pulcini (2006) has been analyzed. Considering the set of times between failures of auto 1, the proposed model was compared to three others. Data set is reported in Table 3. The maximum likelihood estimates of parameters , , , and ln (ℒ) were numerically computed for different cases as described in Table 4. Since the sample size is not much larger than the number of parameters ( < 40), it is recommended to apply the bias-corrected ( ). This latter provides better model choices than other model selection criterions. Note that the is computed by: According to results reported in Table 4, one can conclude that a best fit for this data is the model with lowest . This selected model has both factor improvement and virtual age reduction random and follows the exponential density function. It corresponds to the EME model with parameter approximates 2.54. Figure 3 represents a cumulative number of observed failures and mean cumulative number of the estimated model using AMC data. To approximate the bounds of the BCa confidence intervals, B=1000 Bootstrap samples were simulated. The numerical results are reported in Table 5. We notice that the standard deviations are considerably high. This is due to the small size = 16 of the bootstrap samples. Note that with larger sample sizes, the intervals get narrower and more accurate.  It is worth mentioning that the estimation procedure plays a key role in system reliability analysis. This helps to predict the future failure times and optimize the maintenance policies as needed. Therefore, misevaluated maintenance efficiency would lead to bad maintenance policies and thus cause high losses for companies.

Conclusion
This paper focuses on the study of a hybrid model that combines two principal classes of repairable systems, the virtual age model and the failure rate model. Both factor improvement in failure rate and reduction in virtual age are assumed to be random. The two quantities follow a distribution function for which an exponential density function is considered. The likelihood function of the model has been computed and maximized to obtain parameter estimators. Since there is no analytical form of the estimators, numerical optimization tools have been applied. To test of the method performance, a simulation study was conducted. In the first part, we considered samples generated from the same values so that we test consistency of parameter estimates, and robustness of the estimation method. In the second part, we observed the processes behavior resulting from a combination of random and constant maintenance efficiency. As it is tedious to obtain exact confidence intervals, the bias-corrected accelerated bootstrap method was applied to approximate upper and lower bounds of the confidence interval of each parameter. Finally, an application with real data set has been tested. Obtained results showed that the model is the best fit for the AMC ambassador data.

Conflict of Interest
On behalf of all authors, the corresponding author states that there is no conflict of interest.