Generalized Proportional Model of Relay Protection Based on Adaptive Homotopy Algorithm Transient Stability

: Relay protection equipment is important to ensure the safe and stable operation of power systems. The risks should be evaluated, which are caused by the failure of relay protection. At present, the fault data and the fault status monitoring information are used to evaluate the failure risks of relay protection. However, there is a lack of attention to the information value of monitoring information in the normal operation condition. In order to comprehensively improve monitoring information accuracy and reduce, a generalized proportional hazard model (GPHM) is established to fully exploit the whole monitoring condition information during the whole operation process, not just the monitoring fault condition data, with the maximum likelihood estimation (MLE) used to estimate the parameters of the GPHM. For solving the nonlinear equation in the process of parameter estimations, the adaptive homotopy algorithm is adopted, which could ensure the reversibility of the Jacobi matrix. Three testing cases have been reviewed, to demonstrate that the adaptive homotopy algorithm is better than traditional algorithms, such as the Newton homotopy algorithm, regarding the calculation speed and convergence. Therefore, GPHM could not only reﬂect the real time state of the equipment, but also provide a sound theoretical basis for the selection of equipment


Introduction
Because a relay protection device can curb the deterioration of a power grid by its fast and correct action [1], it is always seen as the first line of defense to ensure the safe and stable operation of power systems.In recent years, there have been frequent blackouts around the world, and most of them are related to the incorrect action of relay protection.Therefore, evaluation of the reliability of relay protection is the focus of many scholars, and the relay protection failure rate is one of the indexes to estimate its reliability [2,3].
At present, the research on the failure rate model is generally based on the time-failure rate model and equipment state model [4].The commonly used fault distribution forms for the time-failure rate model are gamma distribution, Weibull distribution, and exponential distribution [5].The failure rate of exponential distribution is constant and it only represents the accidental failure period.However, the failure rates of most electrical equipment follow the typical curve, namely the bathtub curve, and the bathtub curve includes the early failure period, the accidental failure period, and the three stages of loss failure period.Therefore, the exponential distribution generally is not adopted.By contrast, the Weibull distribution can match the bathtub curve well, so that it has been widely used.Nevertheless, the Weibull distribution only pays attention to the effect of equipment running enlistment age on the failure rate, which ignores the effect of some external factors such as the equipment maintenance on the equipment failure rate [6].Accordingly, Weibull distribution has some limitations.The equipment state model is established based on the current state of the equipment [7], so that it does not take the effects of man-made maintenance and historical condition on the failure rate into account, and so it is difficult to predict the future failure rate.Furthermore, there are differences between countries/power companies.The British EA company and Canada Kinectrics Company focus on the health state [8,9].In China, the failure rate model is an exponential function, which considers the health state as the independent variable [10][11][12], and from the model, the failure rate may increase exponentially when the equipment state worsens, and the failure rate is beyond 1, so that it doesn't conform to the actual situation.In other words, it is difficult to predict the future failure rate.Aiming at the shortcoming of the two kinds of failure rate model, a Weibull proportional hazards model (PHM), which considered the fault diagnosis value of failure time and used maximum likelihood estimation (MLE) to estimate the model parameters, has been proposed in [13].In [14], the proportional covariate model (PCM) was put forward in order to solve the lack of data in the fault interval.For a repairable system, the proportional intensity model (PIM) were first proposed by Kumar [15].Kumar said that the fault rate in a repairable system was affected by many factors, such as operating environment, equipment materials, history operation, design features, and so on.However, the PIM always assumed that the covariates were changed only when failure/maintenance occurred, and maintained constant during the interval of failure/maintenance [16,17].This assumption ignored the influence of the concomitant variables on the failure rate during the failure/maintenance interval.In [18], the PHM was used to estimate the reliability of thin oxide dielectrics, and used the partial maximum likelihood method to estimate the parameters.In [19], the scholars pointed out that the environmental factors influencing relay reliability mainly included temperature, humidity, vibration, and so on.Additionally, the application of Cox-proportional hazards modeling with respect to the effect of ambient temperature on electromagnetic relays was discussed.
For the above model, ignoring the whole monitoring condition values is the common point.Therefore, based on the above model, this paper analyzes the influence of variables during the failure/maintenance interval, and takes the time varying covariates in the failure/maintenance interval into account, not only the monitoring state value at failure time; then, the generalized proportional hazards model (GPHM) Weibull is built.In order to get the expressions of fault rate, the parameters of Weibull distribution are needed to estimate, which involves the solution for nonlinear equations.At present, there are many the solutions for nonlinear equations, such as the Newton method, the least square method [20], the Marquardt method [21], and so on.There are also many achievements on the solution of nonlinear equations, whereas they still have a fatal defect, namely, local convergence [22].Because the initial value must be close to the exact one, the requirement of the initial value is very harsh.In the meantime, the calculation is a large amount which brings certain challenges to the running time and space.Actually, for many nonlinear equations, the initial value is not easy to set, which brings inconvenience to the solution of nonlinear equations.Fortunately, the homotopy algorithm has a large convergence range and its requirement on the initial value is not strict, so that it brings a breakthrough to solve nonlinear equations [23].Howver, the Jacobi matrix of the homotopy algorithm must be reversible, otherwise the homotopy algorithm loses its significance [24].For solving this problem, the adjusting factor is introduced to construct an adaptive homotopy algorithm to ensure the non-singularity of the Jacobi matrix in this paper.Then, the nonlinear equation can be solved.In summary, in order to fully consider the influence of time-varying covariates in the failure/maintenance interval on the failure rate, the GPHM-Weibull is proposed, which can reflect the real-time state.MLE is used to estimate parameters of GPHM, and the adaptive homotopy algorithm is used to solve nonlinear equation, where the piecewise function expression of fault rate is solved.According to the failure rate, the operation personnel can make the differential operation strategy and realize the economic, stable operation of the power system.
Section 2 presents the generalized proportional hazard model.Section 3 describes the adaptive homotopy algorithm.Section 4 discusses the estimation of Weibull distribution parameters, the solution of nonlinear equations and the calculation of the initial value.Section 5 presents a summary of the proposed method and draws relevant conclusions.

Generalized Proportional Hazard Model
With reference to survival function model in medical science [25], based on GPHM, the failure rate model is constructed, and its mathematical expression is as follows: Here, λ 0 (t) is the basic failure rate.Ψ(Z(t)) is the link function, representing the impact of different states Z(t) on failure rate.Z(t) is a vector of covariates, which is composed of n time-varying covariates.Each covariate can represent a particular state.The expression is Z(t) = [Z 1 (t), Z 2 (t) . . .Z n (t)].In practice, the covariate could be an internal variable which can reflect the state of the device, such as the detection information of device.It can also be an external variable which can affect the operation of the device, such as the environmental conditions.In general, the link function can be expressed as follows: Here, γ = (γ 1 , γ 2 . . .γ n ) represents the corresponding regression coefficient of each covariate.The assumptions of GPHM are as follows: (1) The basic failure rate λ 0 (t) subjects to Weibull distribution, and its expression is Here, β is the shape parameter, and η is the scale parameter.
(2) The fault interval is longer than maintenance time, so that the maintenance time can be neglected.
(3) The effect of covariates on the failure rate maintain constant and it cannot be changed with time.
Choosing the best covariates is the key to establish the GPHM.Age, operating environment, maintenance times, health index, and manufacturer are selected as covariates, as shown in Figure 1 [25].
Processes 2019, 7, x FOR PEER REVIEW 3 of 17 Section 2 presents the generalized proportional hazard model.Section 3 describes the adaptive homotopy algorithm.Section 4 discusses the estimation of Weibull distribution parameters, the solution of nonlinear equations and the calculation of the initial value.Section 5 presents a summary of the proposed method and draws relevant conclusions.

Generalized Proportional Hazard Model
With reference to survival function model in medical science [25], based on GPHM, the failure rate model is constructed, and its mathematical expression is as follows: Here, λ0(t) is the basic failure rate.Ψ(Z(t)) is the link function, representing the impact of different states Z(t) on failure rate.Z(t) is a vector of covariates, which is composed of n time-varying covariates.Each covariate can represent a particular state.The expression is Z(t) = [Z1(t), Z2(t)…Zn(t)].In practice, the covariate could be an internal variable which can reflect the state of the device, such as the detection information of device.It can also be an external variable which can affect the operation of the device, such as the environmental conditions.In general, the link function can be expressed as follows: Here, γ = (γ1, γ2 … γn) represents the corresponding regression coefficient of each covariate.The assumptions of GPHM are as follows: (1) The basic failure rate λ0 (t) subjects to Weibull distribution, and its expression is Here, β is the shape parameter, and η is the scale parameter.
(2) The fault interval is longer than maintenance time, so that the maintenance time can be neglected.
(3) The effect of covariates on the failure rate maintain constant and it cannot be changed with time.
Choosing the best covariates is the key to establish the GPHM.Age, operating environment, maintenance times, health index, and manufacturer are selected as covariates, as shown in Figure 1 [25].Health index (HI) reflects the overall health level of the relay protection equipment, which is closely related to the equipment failure rate.In order to facilitate quantitative comparison, HI is divided into five levels (normal, attention, serious, emergency, and fault), and their corresponding values are listed in Table 1.Health index (HI) reflects the overall health level of the relay protection equipment, which is closely related to the equipment failure rate.In order to facilitate quantitative comparison, HI is divided into five levels (normal, attention, serious, emergency, and fault), and their corresponding values are listed in Table 1.In [26], research's results indicated that when the equipment is in the loss period, the failure rate meets Weibull distribution, and Weibull distribution is a function of time.It means that the failure rate is related to age.In [27], the scholars pointed out that the failure rate function directly is multiplied by the age reduction factor, and failure rate increase factor after repairmen.Because the age reduction factor and failure rate increase factor are related to the maintenance type and maintenance times, the failure rate λ k+1 (t) after kth maintenance is defined recursively as: Here, k is the kth maintenance.T k represents the interval for the kth maintenance.β k is age reduction factor, which can simulate the equipment damage caused by each maintenance.α k is failure rate increase factor, which can describe the degree of improvement in equipment failure rate after maintenance.
According to [28], the failure rate can be affected the operating environment.Therefore, the business district and the industrial area differ considerably.Different manufacturers may have different familial defects which may affect the history data.With reference to [27], the values of the operating environment and manufacturers are given as Tables 2 and 3 shown.Therefore, the expression of the failure rate λ( t|Z) which takes the influence of covariates into consideration can be expressed as Here, covariant Z pm , Z cm , Z age , Z env , Z HI , and Z mau respectively represent the times of preventive maintenance and corrective maintenance, age, the operating environment, the HI of equipment and manufacturer; the corresponding coefficient of the covariates are λ pm , λ cm , λ age , λ env , λ HI , and λ mau .Equation (5) illustrates that if you want to get to the expression of failure rate, β/η/γ i is needed in order to estimate.The maximum likelihood function (MLE) and the adaptive homotopy algorithm are used to solve the nonlinear equations.
Then, Equation ( 7) can be converted to: In Equation ( 8), parameter t is introduced to construct a set of homotopy mapping H(x, t) which subjects to: From Equation ( 9), the equations can be obtained, which are t = 0, H(x, 0) = G(x 0 ) and t = 1, H(x, 1) = F(x).Then, the solution of equation F(x) = 0 is transformed into the solution of equation x = x (t) which subjects to equation H(x, 1) = 0. Equation (9) indicates that due to different G(x), there are different homotopy equations.
A Fixed Point Homotopy Algorithm If G(x) = x − x 0 , then a fixed point homotopy algorithm is formed as: Newton Homotopy Algorithm If G(x) = F(x) − F(x 0 ), then a Newton homotopy algorithm is formed as: The derivative of parameter t in H(x, t) = 0 is: If the inverse matrix ( ∂H ∂x ) −1 exists, then: Adaptive Homotopy Equation However, when inverse matrix (∂H/∂x) −1 doesn't exist, the homotopy algorithm will lose its significance.Because the diagonal factor G(x) = diag[e g i (x) ] in the exponential homotopy method is multiplied by F to construct a new homotopy algorithm, it is only feasible in theory, and the calculation is complicated and not suitable for the large-scale nonlinear equation.However, in reference to the idea of exponential homotopy algorithm, the equation for adaptive homotopy algorithm can be obtained as: Because a (1 − t 2 ) I is a nonsingular matrix and when (∂H/∂x) −1 is singular, through adjusting the parameter a, it can account for the diagonal dominance, as long as a is large enough.In the actual calculation, the initial value of parameter a is set as 0. When the Jacobi matrix F (x) becomes singular after some calculation steps, a automatic increase ∆a.Thus, the solution of the Equation ( 7) can be obtained by finding the solution of the homotopy Equation (9).

Numerical Calculation of the Adaptive Homotopy Algorithm
Equation (13) presents that the calculation of nonlinear equations can be converted into the calculation of IVP (initial value problem), which can be expressed as: In order to solve the Equation ( 16), the Euler method is used to estimate, and the Runge Kutta method is used to correct.

Euler method
We begin to track the path from the starting point (t 0 , x 0 ) of the homotopy path, and the Euler method is adopted to estimate the next approximate point (t 1 , x 1 ), so that the expression is: With the iterative equation as: Here, h is the step size, and n is the number of iterations.
The Fouth Runge Kutta method Through using the Runge Kutta method, the local truncation error of the fourth Runge Kutta method is about o (h5) [13], and its calculation speed is fast.Therefore, the fouth Runge Kutta method is adopted to calculate the initial parameters.
Assuming dx dt = −J(x) −1 F(x 0 ) = y(x n , t n ), according to Equation ( 18), the point (t 1 , x 1 ) can be obtained.Then, the point (t 1 , x 1 ) is set as the starting point.For obtaining the next round of prediction-correction, the equations are used as follows: And the prediction-correction process is stopped until t = 1.After several iterations, x n+1 may not be the exact solution x*, according to mathematical convergence theory.If |x n+1 − x n | < ε (ε is the set of coefficients of accuracy), it is considered that the exact solution is found, otherwise the above steps are repeated to perform the predictive-correction process.In order to elucidate the adaptive homotopy algorithm further, Figure 2 shows a flow diagram of this algorithm, which was implemented in MATLAB.

Weibull Distribution Parameters Estimation
There are two reasons which may lead to the relay protection equipment withdrawal from the power system.One is its failure, and the other one is the maintenance.The former possesses fault data and belongs to the corrective maintenance (CM).The latter possesses the truncated data and belongs to the preventive maintenance (PM).If the relay protection equipment is still running at the end of the observation, the data can also be censored.Considering the censored data, based on maximum likelihood, the parameters' estimation can be described as follows: The failure time at 0 < ti < … < tn (T) is observed in the time interval (0 T]. MLE is used to calculate the parameters of GPHM.It is supposed that (ti, zi, δi) and (i = 1, 2,..., n) are the records of Step 1: Set the maximum number of iterations n max , the initial value x 0 and a = 0; Step 2: According to Equation (9), calculate F(x 0 ) and F (x 0 ); Step 3: Judge the reversibility of F (x 0 ), if YES, then turn to step 4, if not, then turn to step 9; Step 4: Based on Equation (18), estimate the next point x 1 ; Step 5: According to Equation (20), modify x 1 and choose the modified x 1 as the initial value for the next round.Meanwhile calculate F(x 1 ) and F (x 1 ), n set as n + 1; Step 6: Assume that nth step has been carried out and obtain the solution x n .judge whether the inequality |x n − x n−1 | < ε is true or not.If YES, output the exact solution x* = x n , otherwise turn step 7; Step 7: Calculate F(x n ) and F (x n ), and judge whether F (x n ) is reversible or not.If YES, based on Equation (18), estimate the next point x n+1 , and n automatically add 1, namely, n = n + 1.If not, turn step 9; Step 8: Judge whether the inequality n > n max is true or not.If YES, the equation without solution.
If not, then turn step 6; Step 9: Adjust the value of a, and define a = a + ∆a.n automatically add 1, namely, n = n + 1, and turn to step 3.

Weibull Distribution Parameters Estimation
There are two reasons which may lead to the relay protection equipment withdrawal from the power system.One is its failure, and the other one is the maintenance.The former possesses fault data and belongs to the corrective maintenance (CM).The latter possesses the truncated data and belongs to the preventive maintenance (PM).If the relay protection equipment is still running at the end of the observation, the data can also be censored.Considering the censored data, based on maximum likelihood, the parameters' estimation can be described as follows: The failure time at 0 < t i < . . .< t n (T) is observed in the time interval (0 T]. MLE is used to calculate the parameters of GPHM.It is supposed that (t i , z i , δ i ) and (i = 1, 2, . . ., n) are the records of failure and maintenance, respectively; n is the total number of events, including all the CM times and all the PM times, so that n can be expressed as n = N cm (t) + N pm (t).t i is the failure time of PM and CM; z i is the state information of the equipment at t i , δ i is corresponding censoring indicator variables, and the equation δ i = 0 represents no failure at t i .The equation δ i = 1 represents failure at t i .
The corresponding likelihood function is given as: Defining Therefore, Equation ( 21) can be converted into Processes 2019, 7, 899 Here, l represents the total number of the relay protection equipment and r represents the number of faulty relay protection equipment.Then, the corresponding log likelihood function is: Taking the partial derivatives of β and η separately: The maximum likelihood functions are: Taking the second derivative operations, then: Then, a second order derivative matrix can be obtained, namely Jacobi matrix:

The Solution of Nonlinear Equations Based on Adaptive Homotopy Algorithm
According to the adaptive homotopy algorithm mentioned, the iteration equation is: Processes 2019, 7, 899 10 of 16 Here, β 0 and η 0 are the initial values of corresponding parameters when k = 0.After the initial values are selected, two parameters can be calculated according to the Equation (34).

Calculation of the Initial Value
Selecting two data points, namely, (t d , λ d ) and (t g , λ g ) in the Weibull distribution, the initial values can be solved by the following equation.
According to the above analysis, MLE can be used to estimate the parameter of GPHM.And the adaptive homotopy algorithm can be used to solve the nonlinear equation in the parameter estimation process.Then, the fault rate model is established.

Case Analysis
Case 1: The machine account and defect information of the relay protection equipment are collected, which are running in the similar environment or have the same type, and its failure rate is shown in Table 4.According the above failure rate, the failure rate curve can be drawn as Figure 3 shown.Based on the adaptive homotopy algorithm and the Newton homotopy algorithm, the failure rate parameters are computed iteratively.Then, failure rate parameters of the Weibull distribution can be obtained and shown in Tables 5 and 6, so that their fitting curve can also be seen in Figure 3.   Through analyzing the data of Tables 5 and 6 and the curves in Figure 3, the iterations number of the adaptive homotopy algorithm is significantly smaller the iterations number of other algorithms.When the equipment is in the random failure period, the relative error of the adaptive Through analyzing the data of Tables 5 and 6 and the curves in Figure 3, the iterations number of the adaptive homotopy algorithm is significantly smaller the iterations number of other algorithms.When the equipment is in the random failure period, the relative error of the adaptive homotopy algorithm is 0.0971, and the Newton homotopy algorithm is 0.1532.When the equipment is in the loss failure period, the relative error of the adaptive homotopy method is 0.1711 and the Newton homotopy algorithm is 0.2518, whose error is large.This is caused by the fact that when the equipment is running, all parts of the equipment occur material fatigue, aging or rust and other undesirable conditions.In order to ensure the normal operation of the equipment, the appropriate maintenance should be done, namely, PM, or CM.However, these two types of maintenance inevitably affect the equipment failure rate.Therefore, the model considering the run-time regardless of the current state will be not correct.The result has big difference with the actual operation.
Case 2: Similar to case 1, the operation data of relay protection equipment is shown in Table 7.Based on Table 7, the failure rate curve of relay protection equipment can be obtained in Figure 4.  Based on the data from Table 7, the convergence characteristics of different algorithms can be obtained.From Table 8, it can be found that the Newton homotopy algorithm is non-converging.However, the adaptive homotopy algorithm can guarantee the singularity of the equation by controlling the parameter a.According to Equation (15), an affects singularity of nonlinear equations, so that the optimal parameter a is selected which can make the iteration number need less.According to Table 8, although there will be differences between the initial values/parameter values/iterations number, the difference of the final parameter estimates' results are in the allowed Based on the data from Table 7, the convergence characteristics of different algorithms can be obtained.From Table 8, it can be found that the Newton homotopy algorithm is non-converging.However, the adaptive homotopy algorithm can guarantee the singularity of the equation by controlling the parameter a.According to Equation ( 15), a n affects singularity of nonlinear equations, so that the optimal parameter a is selected which can make the iteration number need less.According to Table 8, although there will be differences between the initial values/parameter values/iterations number, the difference of the final parameter estimates' results are in the allowed error range.Therefore, it can be proved that the homotopy algorithm is independent of the setting of the initial value, and the different initial values are the different optimal values of a.For example, when the initial value is (1, 1), the iteration number is 10, and then the optimal value of a is 3.While the initial value is (11,1), the iteration number is 35, and then the optimal value of a is −1. Figure 5 shows the curves corresponding to different initial values.From Figure 5, the trend of the curve can fully verify the above analysis results.
Case 3: It is necessary to do the corresponding maintenance for equipment after the device is put into operation for some time.According to Equation (15), it can be found that the choice of maintenance can affect the failure rate.Table 9 gives the operating data of the equipment maintenance.In order to solve the parameters of GPHM, the life data of relay protection equipment is analyzed, firstly.Then, the regression coefficient vector γ is estimated by Statistical Analysis Software (SPSS).On the basis of γ, the adaptive homotopy algorithm is used to estimate the other parameters.Model I: selecting the times of maintenance as a covariate; Models II/III: on the basis of model I, the condition monitoring data of the corrective or preventive maintenance moment are fitted, such as age, the operating environment, the HI of equipment, and the manufacturer; Model IV: on the basis of Models II/III, the condition monitoring data of the whole period (such as age, the operating environment, the HI of equipment, and the manufacturer) are fitted, not just the state of corrective or preventive maintenance moment.
From Table 10, by comparing the estimated log likelihood, the fourth methods are optimal and its log likelihood is −53.930.Obviously, based on the monitoring data of the entire running time, GPHM has better fitting characteristics.According to parameter estimation results, it can be found that the failure maintenances frequency, γ cm , is negative in the above three approaches, so that CM can effectively reduce the failure rate.The preventive maintenances frequency, γ pm , is positive, so that PM can't effectively reduce the failure rate.This conclusion can provide a theoretical basis for maintenance personnel to choose the appropriate maintenance mode, and improve maintenance efficiency in case of the blind maintenance.The results also show that when the environment and the internal state of the covariates deviate from the rated or normal state, the failure rate will be higher.The quantitative analysis is consistent with the experience.Based on GPHM model, the curve of equipment failure rate can be obtained as Figure 5 shown.From Figure 6, we can find that the curve based on GPHM is closed to the actual value, while the Weibull model is away from the actual value, especially in the loss period.That is because in the loss period the equipment has to implement the maintenance.It is inevitable to affect the fault rate.

Conclusions
In order to fully analyze the impact of the whole process of the monitoring state on the failure rate, this paper presents GPHM-Weibull, whose covariates include the times of PM and CM, age, the operating environment, HI, and manufacturer.The baseline function obeys Weibull distribution, and the adaptive homotopy algorithm is adopted to estimate the parameters of Weibull distribution.The regression coefficient vector γ is estimated by SPSS.Finally, three cases have been presented to demonstrate the following conclusions.
(1) The curve drawn by the GPHM is very close to the actual value.However, due to GPHM taking the whole running state of the equipment into consideration, the curve drawn by the Weibull model is away from the actual value, especially in the loss period.
(2) The adaptive homotopy algorithm ensures the singularity of the equation by adjusting the parameter a, and the result of the parameter estimation is less affected by the initial value.

Figure 2 .
Figure 2. Calculation flow chart of adaptive homotopy algorithm.

Figure 3 .
Figure 3. Curve fitting of failure rate.

Figure 3 .
Figure 3. Curve fitting of failure rate.

Figure 4 .
Figure 4.The failure rate curve of equipment.

Figure 4 .
Figure 4.The failure rate curve of equipment.

Figure 5 .
Figure 5.The failure rate curve of equipment.

Table 1 .
Values of Health Index (HI).

Table 2 .
Values of operating environment.

Table 3 .
Of the manufacturer.

Table 4 .
Statistics of failure rate of relay protection.

Table 5 .
Estimation of parameters based on adaptive homotopy algorithm.

Table 6 .
Estimation of parameters based on Newton homotopy algorithm.

Table 4 .
Statistics of failure rate of relay protection.

Table 5 .
Estimation of parameters based on adaptive homotopy algorithm.

Table 6 .
Estimation of parameters based on Newton homotopy algorithm.

Table 7 .
Operation data of relay protection.

Table 8 .
Comparisons of two algorithms.

Table 9 .
Operation data of relay protection equipment.

Table 10 .
Parameter estimation of models with different factor combinations.

Table 9 .
Operation data of relay protection equipment.

Table 10 .
Parameter estimation of models with different factor combinations.