Parameter estimation of lifetime distribution for the meta-action unit with uncertainty failure modes under type-I censored data

Abstract

This paper presents a parameter estimation method for the lifetime distribution of the Meta-Action Unit (MAU) with uncertainty failure modes under type-I censored data. The MAU is regarded as the basic functional unit to accomplish the function of mechanical equipment, and its failure modes are classified according to the abnormal kinematic parameters in Meta-Action (MA), which are more succinct than the traditional mechanical failure modes on parts. However, there is some uncertain information about the failure data and censored data of MAU because of the technology limitations and the space accessibility constraints for monitoring the kinematic parameters of MA, which uncertainty information can impact the parameter estimates of MAU lifetime distribution. In order to avoid the impacts on the estimating accuracy of distribution parameters, the evidential likelihood function based on the belief function theory is constructed in view of the credibility level of the failure data and censored data. In addition, the Evidential Expectation Maximization (E2M) algorithm is proposed to estimate the parameters of the mixed exponential distribution of MAU lifetime under type-I censored data. Finally, an application of an Automatic Pallet Changer (APC) is used to illustrate the validity of the MAU failure modes classification. The simulations of the E2M algorithm are conducted to show that the proposed parameters estimation method can integrate uncertain information in the failure data and the censored data, and obtain more stable results than those based on the conventional Expectation-Maximization (EM).

Introduction

The parameter estimation of lifetime distribution is fundamental work for the reliability assessment or lifetime prediction of the mechanical equipment. With the rapid development of modern mechanical technic, the mechanical equipment become more and more complicated. This complication can contribute to a large number of correlated failure modes appearing on the mechanical parts or components, which aggravates the difficulty of modeling lifetime distribution. In view of the multiple failure modes in mechanical equipment, some equipment-oriented decomposition rules^1,2 are proposed to simplify the analysis of multiple failure modes, which however is usually under the assumption that these failure modes are independent and consistent. The failure physical mechanisms is also important for modeling lifetime distributions for mechanical equipment. ³ Due to the complicated structures and the correlated failure mechanisms, it’s difficult to perform the physical failure analysis for all parts or units for mechanical equipment in engineering practice. To alleviate this difficulty, Zhang et al.^4,5 proposed an equipment-oriented functional decomposition method defined as the Function-Motion-Action (FMA), which is used to obtain the MAU by decomposing the function of mechanical equipment. Up to now, some researches about the MAU have involved in the field of the failure mechanism analysis, ⁶ the failure risk, ⁷ and the failure diagnose, ⁸ etc. Compared with the conventional static failure modes on parts as the analysis object, in these literatures, the action failure modes on MAU is regarded as the analysis object. At this new point of view, the failure behaviors of mechanical equipment are shown through the faulty actions of MAU. Based on the application of MAU failure modes, the reliability allocation technic, ⁹ early failure active elimination method, ¹⁰ and remaining useful life prediction ¹¹ that these researches has been concerned with many scholars. Therefore, the lifetime distribution modeling based on the MAU is an interesting and necessary job. Considering that the MAU is usually used for repairable complex mechanical equipment, these equipment usually have less samples for reliability test because of limited funding and testing sites. The type-I censored data is taken into account for the lifetime distribution modeling of MAU.

On the other hand, as for the reliability test, the failure data are usually imprecise for the unavoidable censored data in the test, which further increases the fitting difficulty of lifetime distribution for mechanical equipment. Lone et al. ¹² proposed the geometric process to analysis the statistical data of accelerated life testing under type-I censored data and, they also considered the occurrence of competing failure causes under type-I and type-II censored data, as well as progressive censoring. ¹³ Although many optimized methods of the single lifetime distribution modeling have been used to tackle the difficulty, the complexity of parameter estimation algorithm for single distributions makes it difficult to apply these methods in practice. Then, many mixed distribution populations^14,15 are proposed to model the lifetime distribution of mechanical equipment. These literatures are devoted to increase the parameter estimation accuracy for the mixed distributions. Generally the mixed distribution represents the probability distribution of each sample in the overall population, ¹⁶ which however discards the correlation between the failure modes and samples. In addition, nonparametric distribution models also can be used to predict reliability, Warr et al. ¹⁷ developed reliability estimates for complex system using component, subsystem, and system failure data with uncertainty estimates based on a Bayesian nonparametric approach. The uncertainty is an unavoidable and fundamental element in the reliability analysis of any complex system. ¹⁸

As an extension of probability theory, the belief function has many advantages in dealing with uncertain information. ¹⁹ In recent years, it has been applied to the classification, ²⁰ decision making, ²¹ failure mode and effects analysis, ²² etc. Zhou and Shi ²³ introduced the E2M algorithm to optimize the parameters estimation when the sample is obtained from Type-II progressively censored data. Different from the commonly adopted methods such as Bayes method or EM for the parameters estimation of mixed distributions, the estimation method based on the evidential likelihood function cooperates the uncertain information of the failure modes which can be used to express the possible degrees of the probability that each sample belongs to a sub-distribution.

The motivation of this paper is the uncertain failure data from the mechanical equipment which would lead to the deviations of parameters estimating for the lifetime distribution. In view of the fact that the MAU failure sample is usually observed at a fixed test termination time, that is, type-I censored test, in this paper, the parameter estimation method of the mixed distribution of MAU with uncertainty failure modes is proposed under the type-I censored data. The novelty of proposed mode in this paper contains two points; the decomposition of mechanical equipment by MA to classify the failure modes of mechanical equipment, the quantifying of uncertainty information in the failure modes of MAU by belief function, and the E2M algorithm alleviates the problem of inaccurate calculation results caused by unreasonable initial value assignment in the traditional EM algorithm. The parameters estimating method in this paper can be used to assess the equipment reliability more accurate and closer to the practical situation where there are some uncertain failure data in filed. Among the censored failure data, for example, the multiple failure modes would lead to the ambiguity of failure time to a lifetime distribution.

Based on the mentioned above, the rest of this paper is organized as follows. Specifically, the concept of MA and the failure modes of MAU are introduced in Section 2. Then, Section 3 proposes the notion of belief function and the evidential likelihood function based on the mixed Exponential distribution, and the parameters of mixed Exponential distribution of MAU are estimated by E2M algorithm. Next in Section 4, an engineering example is used to show the advantages of MAU failure modes and the validity of the proposed estimation method. In addition, considering the assignment problem of the initial values, a simulation example is used to demonstrate the superiority of the E2M algorithm over the traditional EM algorithm under type-I censored data. The conclusion is summarized in Section 5.

Failure modes of MAU

Before figuring out the MAU failure modes, the problems about what the MAU is and why it can be used as an analytical object of mechanical failure modes instead of the parts or other physical bodies need to be known. The MA is regarded as a basic action used to realize the function of mechanical equipment. The mapping relationship between the mechanical equipment function and the MA is illustrated in the Figure 1, which can be regarded as a decomposition rule for the mechanical equipment function. The decomposition rule can be implemented through two steps. Firstly, dividing the function of equipment into the motions of assemblies next, dividing the motion into the actions of parts and these actions are defined as MA. For example, the milling function of the machine tool center is realized through the spindle rotating action, the piston translating action, the rack translating action, etc. The more details about the FMA can be seen in the literature. ²⁴ The basic MA is shown in the Figure 2.

Figure 1.

FMA decomposition rule.

Figure 2.

The translation MA and the rotation MA.

The MAU, defined as a structural unified unit of MA, contains such elements as the input, the output, the supporter, the fastener and the middleware. ²⁴ Additionally, these components are combined together as a functional unit to execute the MA and are assembled with each other such as the gear rotating unit, the workbench translating unit as shown in the Figure 3.

Figure 3.

The rotation MAU (a) and the translation MAU (b).

According to definition of MAU, once the MA fails the function of mechanical equipment can’t be realized, which means the equipment has failed. So the equipment failure data can be collected in terms of the failure modes of MAU. As for a mechanical equipment, compared with the conventional failure modes on the individual part, ²⁵ the MAU failure modes are upon the abnormal actions of the output element in the MAU, and the failure modes of other elements in the MAU can be regarded as the failure reasons of these abnormal actions.

For example the gear rotating unit as a MAU in the Figure 3, the gear is the output element, and the gear rotating is the MA, the failure that the fluctuation excessive of gear rotating can be seen as one of failure modes of the gear rotating unit, and the failure reasons could be the bearing fatigue, the shaft eccentric or other things in the gear rotating unit. For screening the MA failure data conveniently from the field warranty records, the failure modes of MAU have been defined in the literature. ⁶ In this paper, an extended classification of the MAU failure modes is collected in Table 1.

Table 1.

The MA failure modes.

MAU failure modes		Failure criteria
Index	Mode	Translation	Rotation
Speed	No action	$v=0$	$w=0$
	Speed insufficient	$v<[v{]}_{min}$	$w<[w{]}_{min}$
	Speed excessive	$v>[v{]}_{max}$	$w>[w{]}_{max}$
Acceleration	Acceleration insufficient	$\varphi <[\varphi {]}_{min}$	$\beta <[\beta {]}_{min}$
	Acceleration excessive	$\varphi >[\varphi {]}_{max}$	$\beta >[\beta {]}_{max}$
Displacement	Displacement insufficient	$l<[l{]}_{min}$	$\theta <[\theta {]}_{min}$
	Displacement excessive	$l>[l{]}_{max}$	$\theta >[\theta {]}_{max}$
Fluctuation	Fluctuation excessive	$\sigma >[\sigma {]}_{max}$	$\delta >[\delta {]}_{max}$

The values $[.{]}_{min}$ and $[.{]}_{max}$ mean the allowed interval values of the MA performances, which can be obtained from the design specification. From the Table 1, the MA can exhibit eight different failure modes which can be classified into speed, acceleration, displacement, and fluctuation, the four physical failure indexes.

Since some parts in the mechanical equipment are invisible or difficult to detect its performance, it becomes difficult for warranty engineers to judge the primary failure causes (failure modes). For the failure modes of MAU, there may be uncertain information in the four indexes of MAU failure modes. For example, when an action’ speed becomes smaller than that at normal, it may be judged as the “Speed insufficient” failure mode or the “Acceleration insufficient.” Consequently, the heterogeneous failure samples, from failure modes of the MAU, cannot be explicitly classified into the homogeneous sub-distributions, which would block the estimating accuracy of mixed distributions.

Evidential likelihood function of the mixed exponential distribution of MAU

Mixed exponential distribution

The different failure modes usually express different failure mechanisms, which could correspond to different lifetime distributions. ²⁶ There can be more than one failure mode for any MAU, so the single distribution is unsuitable for a MAU lifetime modeling. In the field of lifetime distribution modeling for the mechanical equipment, it’s well known that the Weibull models have been a general distribution to fit the lifetime of mechanical equipment. ²⁷ In view of the multiple failure modes, the mixed Weibull distribution has been widely used to model the lifetime of mechanical equipment or parts,^28,29 and the mixed Exponential distribution is a special form of the Weibull distribution when the shape parameter is constant to 1. In this paper, a finite mixed Exponential distribution is used to fit the lifetime of MAU. The mixed Exponential distribution with $m$ sub-distributions has its PDF and survival function respectively as:

$$f(t;\Theta )=\sum _{i=1}^m{\omega }_i{\lambda }_i\exp (-{\lambda }_{i}t)$$

(1)

$$R(t;\Theta )=\sum _{i=1}^m{\omega }_i\exp (-{\lambda }_{i}t)$$

(2)

Where ${\omega }_i(i=1,\dots ,m)$ is the weighting factor, and $\sum _{i=1}^m{\omega }_i=1$ $(0<{\omega }_i<1)$ , $\Theta ={\left\{({\omega }_i,{\lambda }_i)\right\}}_{i=1}^m$ is called the parameter vector of a mixed Exponential distribution.

Belief function

The theory of belief function is a generalization of the theory of probability, ²³ which can be used to express the support degrees assigned to some subsets (failure modes) in the domain $\Omega$ through the mass function $h:2^{\Omega }\to [0,1]$ . Compared to the probability models that are given to the atomic element of $\Omega$ , the mass function can be given to any subset of $\Omega$ . The belief function $\mathrm{Bel}(A)$ and plausibility function $\mathrm{Pl}(A)$ defined as:

$$\mathrm{Bel}(A)=\sum _{B\subseteq A}h(B)$$

(3)

$$\mathrm{Pl}(A)=\sum _{B\cap A\ne \varnothing }h(B)$$

(4)

The former $\mathrm{Bel}(A)$ denotes as the degree to which the evidence supports the $A$ , and the latter $\mathrm{Pl}(A)$ measures an upper bound on the degree of potential support to the $A$ . The function $\mathrm{pl}:\Omega \to [0,1]$ such that $\mathrm{pl}(t)=\mathrm{Pl}(\{t\})$ is called the plausibility function, which can be used to express the support degrees of failure data belonging to a certain failure mode of the MAU or to a sub-distribution of the MAU mixed lifetime distribution. The theory of belief functions can be considered as having greater expressive power than possibility theory, and more details about the belief function theory can be found in the literature. ³⁰

Evidential likelihood function under type-I censored data

Let $T=\{t_1,t_2,\dots t_n\}$ be a stochastic variable denoted as the failure time of the MAU in the domain ${\Omega }_T$ . In view of the uncertainty information between failure time and failure modes in MAU, as well as the imprecision of censored data, the likelihood function of the mixed distribution can be described as:

$$\begin{array}{cc}\multicolumn{1}{c}{L(\theta ;h)}& =\sum _{i=1}^{r}h(A_i)L(\theta ;A_i)\hfill \\ \multicolumn{1}{c}{}& =\sum _{t\in {\Omega }_T}f(t;\theta )\sum _{t\in A_i,A_i\subseteq {\Omega }_T}h(A_i)\hfill \\ \multicolumn{1}{c}{}& =\sum _{t\in {\Omega }_T}f(t;\theta )\mathrm{pl}(t)\hfill \\ \multicolumn{1}{c}{}& ={\mathrm{E}}_{\theta }[\mathrm{pl}(T)]\hfill \\ \multicolumn{1}{c}{}& =L(\theta ;\mathrm{pl})\hfill \end{array}$$

(5)

Where $L(\theta ;A_i)$ is the likelihood function of imprecise failure data ( $t_i\in A_i$ ). From the formula (5), the index $L(\theta ;\mathrm{pl})$ can be seen as the expectation of $\mathrm{pl}(t)$ with respect to $f(t;\theta )$ , which is often noted as the evidential likelihood function. ³¹

For a MAU with $K$ failure modes, its failure data under type-I censored time (completely observed failure time: $(t_1,t_2,\dots ,t_N)$ and censoring time: $(s_1,s_2,\dots ,s_S)$ ) contain some uncertain information. Its likelihood function can be shown as:

$$L(\Theta |t)=\underset{i=1}{\overset{N}{\Pi }}f(t_i;\Theta )\underset{j=1}{\overset{S}{\Pi }}R(s_j;\Theta )$$

(6)

Let $Z_i=(z_{i1},z_{i2},\dots ,z_{\mathrm{iK}})$ is the vector that the completely observed data $t_i$ belongs to the $k_{\mathrm{th}}$ failure mode, and $R_j=(r_{j1},r_{j2},\dots ,r_{\mathrm{jK}})$ is the vector that the censored data $s_j$ belongs to the $k_{\mathrm{th}}$ failure mode. The symbols $z_{\mathrm{ik}}$ and $r_{\mathrm{jk}}$ are Boolean variables, they satisfy:

$$z_{\mathrm{ik}}=\{\begin{array}{c}1,{\mathrm{t}}_i\in A_k\hfill \\ \multicolumn{1}{c}{0,\mathrm{otherwise}}\end{array};{\mathrm{r}}_{\mathrm{jk}}=\{\begin{array}{c}1,{\mathrm{s}}_j\in A_k\hfill \\ \multicolumn{1}{c}{0,\mathrm{otherwise}}\end{array}$$

(7)

According to the formulas (1) and (2), the likelihood function (6) can be rewritten as:

$$\begin{array}{cc}\multicolumn{1}{c}{}& L(\Theta |t)=\underset{i=1}{\overset{N}{\Pi }}\sum _{k=1}^K{({\omega }_k{\lambda }_k\exp (-{\lambda }_k{t}_i))}^{z_{\mathrm{ik}}}\hfill \\ \multicolumn{1}{c}{}& \underset{j=1}{\overset{S}{\Pi }}\sum _{k=1}^K{({\omega }_k\exp (-{\lambda }_k{s}_j))}^{r_{\mathrm{jk}}}\hfill \end{array}$$

(8)

Then, the evidential likelihood function under type-I censored data for the MAU based on the mixed Exponential distribution can be denoted as

$$\begin{array}{cc}\multicolumn{1}{c}{\log L(\Theta |t)}& =\sum _{i=1}^N\log \sum _{k=1}^K(z_{\mathrm{ik}}){\omega }_k{\lambda }_k\exp (-{\lambda }_k{t}_i)\hfill \\ \multicolumn{1}{c}{}& +\sum _{j=1}^S\log \sum _{k=1}^K(r_{\mathrm{jk}}){\omega }_k\exp (-{\lambda }_k{s}_j)\hfill \end{array}$$

(9)

Parameters estimation of the mixed Exponential distribution of MAU

The EM method is one of the main approaches to estimate the parameters for a mixed distribution.^16,32 In this paper, an E2M algorithm is proposed for modeling the lifetime distribution of the MAU under type-I censored data. The theoretical difference between EM and E2M methods is that the latter takes the uncertainty into the expectation function in the EM. Similarly to the EM, this algorithm is also executed through two steps as:

E-step: the Expectation to the evidential likelihood function of the MAU failure samples. According to the Jensen inequality, the expectation of $\log L(\Theta |t)$ with respect to the $\mathrm{pl}$ can be shown as:

$$\begin{array}{ccc}\multicolumn{1}{c}{}& E[\log L(\Theta |t)|\mathrm{pl};{\Theta }^{(q)}]\ge \hfill & \hfill \\ \multicolumn{1}{c}{}& \sum _{i=1}^N\sum _{k=1}^K[\log {\omega }_k+\log {\lambda }_k-{\lambda }_k{t}_i]E(z_{\mathrm{ik}}|\mathrm{pl};{\Theta }^{(q)})\hfill \\ \multicolumn{1}{c}{}& +\sum _{j=1}^S\sum _{k=1}^K[\log {\omega }_k-{\lambda }_k{s}_j]E(r_{\mathrm{jk}}|\mathrm{pl};{\Theta }^{(q)})=Q(\Theta ,{\Theta }^{(q)})\hfill \end{array}$$

(10)

Where the expectation $E(z_{\mathrm{ik}}|\mathrm{pl};{\Theta }^{(q)})$ can be calculated as follows:

$$\begin{array}{cc}\multicolumn{1}{c}{E(z_{\mathrm{ik}}|\mathrm{pl};{\Theta }^{(q)})}& =p(z_{\mathrm{ik}}=1|\mathrm{pl};{\Theta }^{(q)})\hfill \\ \multicolumn{1}{c}{}& \hfill \\ \multicolumn{1}{c}{}& =\frac{{{\omega }_k}^{(q)}{{\lambda }_k}^{(q)}\exp (-{{\lambda }_k}^{(q)}{t}_i)p{l}_{t_{i,k}}}{\sum _{l=1}^K{{\omega }_l}^{(q)}{{\lambda }_l}^{(q)}\exp (-{{\lambda }_l}^{(q)}{t}_i)p{l}_{t_{i,l}}}={a_{i,k}}^{(q)}\hfill \end{array}$$

(11)

Similarly, we can get

$$E(r_{\mathrm{jk}}|\mathrm{pl};{\Theta }^{(q)})=\frac{{{\omega }_k}^{(q)}\exp (-{{\lambda }_k}^{(q)}{t}_i)p{l}_{s_{j,k}}}{\sum _{l=1}^K{{\omega }_l}^{(q)}\exp (-{{\lambda }_l}^{(q)}{t}_i)p{l}_{s_{j,l}}}={b_{j,k}}^{(q)}$$

(12)

M-step: the Lagrange method can be used to maximize the $Q(\Theta ,{\Theta }^{(q)})$ with $\sum _{k=1}^K{\omega }_k=1$

$$L(\Theta ,\beta )=Q(\Theta ,{\Theta }^{(q)})+\beta (\sum _{k=1}^K{\omega }_k-1)$$

(13)

Then take the partial derivative of equation (13), we can get:

$$\frac{\partial L(\Theta ,\beta )}{\partial {\omega }_k}=\frac{\sum _{i=1}^N{a_{i,k}}^{(q)}+\sum _{j=1}^S{b_{j,k}}^{(q)}}{{\omega }_k}+\beta =0$$

(14)

$$\begin{array}{cc}\multicolumn{1}{c}{}& \frac{\partial L(\Theta ,\beta )}{\partial {\lambda }_k}=\frac{1}{{\lambda }_k}\sum _{i=1}^N{a_{i,k}}^{(q)}-\sum _{i=1}^N{t}_i{a_{i,k}}^{(q)}\hfill \\ \multicolumn{1}{c}{}& -\sum _{j=1}^S{s}_j{b_{j,k}}^{(q)}=0\hfill \end{array}$$

(15)

Further, we can get:

$${{\omega }_k}^{(q+1)}=\frac{{C_k}^{(q)}}{\sum _{k=1}^K{C_k}^{(q)}}$$

(16)

$${{\lambda }_k}^{(q+1)}=\frac{\sum _{i=1}^N{a_{i,k}}^{(q)}}{\sum _{i=1}^N{t}_i{a_{i,k}}^{(q)}+\sum _{j=1}^S{s}_j{b_{j,k}}^{(q)}}$$

(17)

$${C_k}^{(q)}=\sum _{i=1}^N{a_{i,k}}^{(q)}+\sum _{j=1}^S{b_{j,k}}^{(q)}$$

(18)

The E2M needs to be repeated iteratively until:

$$\{\begin{array}{c}\sqrt{\sum _{k=1}^K{({{\omega }_k}^{(q+1)}-{{\omega }_k}^{(q)})}^2}\le \xi \hfill \\ \multicolumn{1}{c}{\sqrt{\sum _{k=1}^K{({{\lambda }_k}^{(q+1)}-{{\lambda }_k}^{(q)})}^2}\le \zeta }\end{array}$$

(19)

Where $\xi$ and $\zeta$ are given small positive constants.

In order to demonstrate the fitted accuracy with different estimation methods for the same test sample, the Root Mean Square Error (RMSE) and Normalized Root Mean Square Error (NRMSE) as two indicators are used to express the fitted accuracy of proposed method as shown in formula (20) and (21)

$$\mathrm{RMSE}=\sqrt{\sum _{i=1}^n{[F_n(t_i)-\stackrel{\wedge }{F_n}(t_i)]}^2/\sum _{i=1}^n{[F_n(t_i)-\stackrel{\wedge }{F_n}(t_i)]}^{2}nn}$$

(20)

$$\mathrm{NRMSE}=\sqrt{\sum _{i=1}^n{[F_n(t_i)-\stackrel{\wedge }{F_n}(t_i)]}^2/\sum _{i=1}^n{[F_n(t_i)-\stackrel{\wedge }{F_n}(t_i)]}^2\sum _{i=1}^n{[F_n(t_i)]}^2\sum _{i=1}^n{[F_n(t_i)]}^2}$$

(21)

Where, $n$ is the sample size, $F_n(t_i)$ is the indicator value of experience distribution, $\stackrel{\wedge }{F_n}(t_i)$ is the indicator value of estimated distribution based on E2M. The smaller the indicators values, the better the fitted accuracy. The parameter estimation flow of the lifetime distribution for the MAU is shown in Figure 4.

Figure 4.

Flow chart of parameters estimation of lifetime distribution for the MAU.

Application

As one of the critical assemblies on the machine tool center, the automatic pallet changer (APC) is used to exchange the position of the parts to be processed and parts processed already, which improves the processing efficiency for machine tool centers. The APC consists of gear shaft, pallet, spindle, lifting piston, end bearing, seal ring, lifting hydraulic cylinder, sleeve, rack, and sliding block, etc. as the Figure 5 shows.

Figure 5.

Structure diagram of APC.

According to the structures diagram of the APC, there are mainly hydraulic circuit components and mechanical transmission parts being used to finish the function of changing position. Specifically, the lifting piston is driven to move up and down by the hydraulic power, and is fixed with the pallet by the bolts. Then, the gear shaft, meshed with the rack, is driven to turn back and forth by the hydraulic power. The gear shaft is fixed with the pallet. The APC function can be decomposed by FMA as shown in the Figure 6.

Figure 6.

FMA diagram of APC.

In view of the coupling of electromechanical and hydraulic in the APC, as well as the instability of hydraulic power, there are usually exhibiting various failure modes on the APC. We collected 50 failure samples about the APC from six machine tool centers at the reliability test field, and found that the most failures occurred when the pallet rotated. Taking the pallet rotating action unit as the object of failure analysis, in this paper, the different descriptions on the failure modes of the parts and the MA are collected in Table 2.

Table 2.

Failure modes of the pallet rotating action unit.

On parts			On MA
No.	Failure modes	Count	No.	Failure modes	Count
1	Pallet can’t rotating	7	1	No action	11
2	Pallet stuck	4
3	Pallet inflexible rotation	13	2	Fluctuation excessive	31
4	Noise when pallet rotating	8
5	Pallet buff failure	7
6	Pallet tilt	3
7	Pallet exchanging position error	6	3	Displacement insufficient	8
8	Valve oil leakage when pallet rotating	2

From the records in the Table 2, the number of failure modes is reduced from 8 to 3. The traditional failure modes on parts are often expressed through subjective methods such as human vision and hearing, which could lead to a situation that the same failure phenomenon is expressed as different failure modes because of the subjectivity of the maintainers. But, based on the MA failure modes taxonomy, in this paper, the failure modes of the pallet rotating action unit are expressed through the objective kinematic parameters values detected from various sensors technics, which greatly reduce the number of failure modes without discarding the time to failure data.

Although the numbers of failure modes of the pallet rotating action unit are decreased to 3, the numbers of time to failures with the completely observed failure data and the partly censored data are all constant as the Table 3.

Table 3.

The failure time and censored time of the pallet rotating action unit.

Failure time (h)	5	6	9	11	12	15	17	18	20	23	24	29	32	34	38
	43	47	52	66	75	86	94	125	151	175	188	199	248	261	300
	352	378	407	451	518	592	648	711	817	905	1012	1236	1412	1994
Censored time (h)	9+			23+		50+	190+			302+	645+

In order to demonstrate the validity of the proposed method for the parameter estimation, the dataset is fitted by the mixed distribution model, and the mixed Exponential distribution with two sub-distributions ( $K=2$ ) is found suitable for fitting the given data. The fitted FDF can be expressed as:

$$\hat{f}(t)={\hat{\omega }}_1{\hat{\lambda }}_1{e}^{-{\hat{\lambda }}_{1}t}+{\hat{\omega }}_2{\hat{\lambda }}_2{e}^{-{\hat{\lambda }}_{2}t}$$

(21)

Taking the traditional EM method with the specified initial values ( ${{\omega }_1}^{(0)}=0.5,{{\omega }_2}^{(0)}=0.5,{{\lambda }_1}^{(0)}=1e-3,$ ${{\lambda }_2}^{(0)}=1e-4$ ) to get the parameters as shown:

$${\hat{\omega }}_1=0.3564,{\hat{\omega }}_2=0.6436,{\hat{\lambda }}_1=0.0283,{\hat{\lambda }}_2=0.0019$$

(21)

Suppose there are some uncertain information about the failure data of the pallet rotating action unit. Specifically, when a failure occurs on the APC, the failure time data can be collected quickly through the sensors, but for the failure modes on the pallet rotating (as a MA), they could not be corresponding to the failure time accurately. In addition, for the APCs from six machine tool centers, the censored time data also could not be related to the failure modes correctly. In view of these circumstances that a failure sample could be corresponded to multiple failure modes, according to the belief function theory, the plausibility function $\mathrm{pl}$ for each sample can be obtained by the total uncertainty measure method. ³⁰ Based on the number of sub-distribution of the mixed distribution, we can set the clustering parameter c = 2, and the values of potential support to subsets are shown in Tables 4 and 5.

Table 4.

Plausibility values of failure time to the sub-distributions.

No.	$p{l}_{t_i,1}$	$p{l}_{t_i,2}$	No.	$p{l}_{t_i,1}$	$p{l}_{t_i,2}$
1–3	0.5000	0.5000	31	0.5000	0.5000
4	0.9999	0.0001	32	0.0002	0.9998
5–7	0.5000	0.5000	33	0.0097	0.9903
8	0.9999	0.0001	34	0.5000	0.5000
9–29	0.5000	0.5000	35	0.0001	0.9999
30	0.0001	0.9999	36–44	0.5000	0.5000

Table 5.

Plausibility values of censored time to the sub-distributions.

No.	$p{l}_{s_j,1}$	$p{l}_{s_j,2}$	No.	$p{l}_{s_j,1}$	$p{l}_{s_j,2}$
45	0.9997	0.0003	48	0.5000	0.5000
46	0.9995	0.0005	49	0.0001	0.9999
47	0.5000	0.5000	50	0.0001	0.9999

Taking the same initial values, the E2M method proposed in this paper is used to get the parameters as shown:

$${\hat{\omega }}_{b1}=0.4354,{\hat{\omega }}_{b2}=0.5646,{\hat{\lambda }}_{b1}=0.0210,{\hat{\lambda }}_{b2}=0.0017$$

(21)

Furthermore, in order to verify the validity of the EM and E2M algorithm, Bayesian technique (BT) is used to estimate the parameters of the mixed exponential distribution. According to the method in the literature, ³³ we choose Gamma distribution as the prior distribution of distribution parameters ${\lambda }_i$ and Beta distribution as the prior distribution of weighting factor ${\omega }_i$ , and the posterior estimates of three parameters of the mixed exponential distribution are

$${\hat{\omega }}_{b1}=0.4354,{\hat{\omega }}_{b2}=0.5646,{\hat{\lambda }}_{b1}=0.0210,{\hat{\lambda }}_{b2}=0.0017$$

(21)

The fitted results by the BT are shown in the Figure 9.

The fitted results based on the proposed E2M method, traditional EM method and BT are compared to the empirical distribution (ED) as shown in Figures 7 to 9, respectively. The fitting accuracy based on the three methods are noted in Table 6.

Table 6.

Fitting accuracy and number of iterations based on the two methods.

Methods	RMSE	NRSME
E2M	0.0191	0.0331
EM	0.0316	0.0549
BT	0.0275	0.0478

Figure 7.

Fitted results by the E2M method.

Figure 8.

Fitted results by the EM method.

Figure 9.

Fitted results by the Bayesian technique (BT) method.

According to the fitted results in Figures 7 to 9, as well as the fitting accuracy in Table 6, the relatively similar results are obtained by the three methods, which means the EM method and the improved version E2M are applicable. The fitted curve in Figure 7 is much closer to the empirical function than that in the Figure 8, the values of RMSE and NRSME based on the E2M method are all smaller than that based on the EM method. That means the E2M method proposed in this paper is better to estimate the parameters of the lifetime distribution of pallet rotating action unit under type-I censored data. Furthermore, suppose that the plausibility function $\mathrm{pl}$ of each sample is equal in the interval (0, 1), the parameters estimated would be same as for the EM’s, which means the EM method can be seen as the special situation of E2M method when the potential support degrees to each sample are the same.

The advantages of E2M method is that the uncertainty information of each sample is considered when establishing the likelihood function of the whole samples. The uncertainty information about failure data of MAU is vital for the parameters estimation of lifetime distribution.

In view that the estimated accuracy of traditional EM algorithm depends on the assigned initial values largely, we suppose that each sub-distribution has the same parameters ${{\lambda }_1}^{(0)}={{\lambda }_2}^{(0)}=1e-3,{{\omega }_1}^{(0)}=0.5,{{\omega }_2}^{(0)}=0.5$ and then compare the results of E2M and EM, as shown in Figures 10 and 11, respectively.

Figure 10.

Fitted result by E2M under identical initial values.

Figure 11.

Fitted result by EM under identical initial values.

As we can see in Figures 10 and 11, the fitted result based on the E2M isn’t changed but for the EM, the fitted result is worse. It means that the evaluation results based on the method proposed in this paper will not be affected when the initial values are assigned to equal.

At the same time, according to the initial results estimated by the probability graph method and MLE, we find that the values of $\lambda$ are not less than $1e-5$ when the fitted curve is close to the empirical distribution of failure data. In order to verify the good robustness of E2M algorithm, the standard normal distribution $N(u,\sigma )$ was used to extend the range of initial results for more iterative calculations. The generalized parameters ${{\lambda }_i}^{(0)},{{\omega }_k}^{(0)}$ are shown as:

$${{\lambda }_1}^{(0)}=|1e-3+N(0,1\mathrm{e}-3)|$$

(21)

$${{\lambda }_2}^{(0)}=|1e-4+N(0,1\mathrm{e}-3)|$$

(21)

$${{\omega }_1}^{(0)}=|0.5+N(0,0.1)|,$$

(21)

$${{\omega }_2}^{(0)}=1-{{\omega }_1}^{(0)}$$

(21)

Then, the parameters are estimated using the two methods respectively (iterations: 100), and the mean values and standard deviations (SD) of these parameters are shown in Table 7.

Table 7.

Mean values of fitting and standard deviations based on the two methods.

Parameters	E2M		EM
	Mean	SD	Mean	SD
${\hat{\lambda }}_1$	0.0231	0.0049	0.0180	0.0129
${\hat{\lambda }}_2$	0.0019	0.0006	0.0122	0.0129
${\hat{\omega }}_1$	0.4280	0.0123	0.4684	0.1401
${\hat{\omega }}_2$	0.5720	0.0123	0.5316	0.1401

From the Table 7, we can see that the values of SD based on the E2M method are smaller than that based on the EM method, which means the estimated parameters using the proposed E2M method are more stable than that from the traditional EM method when different initial values are devoted into. According to the application of E2M method for the pallet rotating action unit, we can see that the uncertain information in the failure data of MAU are significant to the parameters estimation of lifetime distribution.

Conclusion

In this paper, a method for estimating the parameters of mixed Exponential distribution of the MAU under type-I censored data is proposed. The main advantage of the proposed method is that it can take the uncertainty information of failure data in the MAU especially for the censored data into consideration. From the results of the example APC, the failure modes classification based on MA is more succinct than that from the APC parts. On the other hand, The E2M algorithm can largely alleviate the problem of inaccurate calculation results caused by unreasonable initial value assignment in the traditional EM algorithm. The estimated results by the E2M method proposed in this paper are more stable and accurate than those from EM algorithm.

Of course, except for the mixed Exponential distribution, there are other mixed lifetime distributions worth researching such as the generalized mixed Weibull distribution, the Exponential distribution mixed with Pareto model, etc. On the other hand, as we can see, the support degrees of uncertain information are significant for the parameters estimation, how to build a generic model to assign the plausibility function values for MAU uncertain failure data is another issue for further study.