Identifying the change time of multivariate binomial processes for step changes and drifts

In this paper, a new control chart to monitor multi-binomial processes is first proposed based on a transformation method. Then, the maximum likelihood estimators of change points designed for both step changes and lineartrend disturbances are derived. At the end, the performances of the proposed change-point estimators are evaluated and are compared using some Monte Carlo simulation experiments, considering that the real change type presented in a process are of either a step change or a linear-trend disturbance. According to the results obtained, the change-point estimator designed for step changes outperforms the change-point estimator designed for linear-trend disturbances, when the real change type is a step change. In contrast, the change-point estimator designed for linear-trend disturbances outperforms the change-point estimator designed for step changes, when the real change type is a linear-trend disturbance.


Background
Enabling process engineers and practitioners to monitor the state of processes by distinguishing between common and special causes of variation, control charts have been one of the most widely used statistical process control tools in industries. When a control chart generates a signal due to an out-of-control state of a process, process engineers initiate a search to identify the root causes of process variation. However, the signaling time is different from the first time that a change manifests itself into the process (change point), and in most cases there exists a considerable time lag between them. Therefore, knowing the exact time of the change would simplify the search for root cause identification and corrective action implementation in order to improve the quality of processes.
A binomial distribution with parameters n and p is used to model the number of non-conforming items in a sequence of n trials for an attribute quality characteristic, in which p denotes the process fraction nonconforming. In the case of multi-attribute processes, the multivariate binomial distribution with parameters n and p k (k = 1, 2,. . ., q) can be used to model the number of non-conforming items in a sequence of n trials for each attribute quality characteristic, in which p k denotes the process fraction non-conforming of the kth quality characteristic. The problem of estimating the time of a change in a univariate binomial distribution for different change types of step change, linear-trend disturbance (drift), multiple-step changes, monotonic change, and other change types defined as a combination of the above has been addressed in the literature. In a stepchange disturbance, the parameter of the process shifts to an out-of-control level and remains in that level until corrective measures are taken. Although a step change can occur because of tool breakage, the parameter of a process may gradually drift. In other words, linear-trend disturbances can occur as a result of tool wear or equipment depreciation. Moreover, a process parameter may change due to multiple-step changes. A change in several parameters at different times or a change in one parameter at different times can be the results of this type of change. In monotonic changes, however, the only assumption is that the change type can be described as belonging to a family of monotonic change types. Samuel and Pignatiello (2001) derived the maximum likelihood estimator (MLE) of a step-change point for binomial processes. Evaluating the estimator for different shift sizes, they concluded that it provides accurate results. In the case of a step change, Perry and Pignatiello (2005) compared the MLE of the process fraction nonconforming with the built-in change-point estimators of the binomial CUSUM and EWMA control charts described by Page (1954) and Nishina (1992), respectively. They showed that the MLE provides better results than the other two estimators do. Perry et al. (2007) derived the MLE of a monotonic change point of the process fraction non-conforming and compared it with the MLE of a step-change point derived by Samuel and Pignatiello (2001). They concluded that using their proposed estimator is better when the type of change is only known to be monotonic. Zandi et al. (2011) proposed an MLE to estimate the time of a linear-trend change in the process fraction non-conforming and compared it to the MLEs derived for step and monotonic change disturbances proposed by Samuel and Pignatiello (2001) and Perry et al. (2007), respectively.
To the best of the authors' knowledge, there are only two research works available in the literature to estimate the change point of multi-attribute processes. In the first work, Niaki and Khedmati (2012a) proposed an MLE of a step-change point for multivariate Poisson processes and evaluated its performances using some numerical experiments. They showed that their estimator leads into good results for all process dimensions. In the second work, Niaki and Khedmati (2012b) derived an MLE of a lineartrend change in multivariate Poisson processes and compared its performances with those of an MLE of a step change for both step-change and linear-trend disturbances. However, there is no research work concerning the change-point estimation of multivariate binomial processes, whereas in many situations, monitoring processes containing several quality characteristics of attribute type, each following a binomial distribution, is a need. Therefore, in this paper, we first propose a new multi-attribute control chart to monitor multi-attribute processes, and then, we derive the MLEs of both a step-change and a linear-trend disturbance in the process fraction nonconforming. Following a signal from the proposed multiattribute control chart, the MLEs are applied to estimate the step-change and linear-trend disturbance points.
The rest of this paper is organized as follows: In the next section, a new multi-attribute control chart is proposed based on a transformation method to monitor multiattribute processes. In the 'Multivariate binomial process modeling and MLE derivation' section, the model of the process behavior is described and the MLEs of a stepchange and a linear-trend disturbance are derived. The performances of the two change-point estimators are evaluated and are compared in the 'Results and discussion' section, in the presence of both step-change and lineartrend disturbances. Finally, the paper is concluded in the 'Conclusions' section.
The proposed multi-attribute control chart As a generalization of Shewhart-type control charts, most of multivariate control charts have been constructed based upon the multivariate normal distribution. Despite the existence of some skewness in multi-attribute distributions, the normal distribution is symmetric with zero skewness. As a result, using normal-based control charts with zero skewness for multi-attribute processes causes some problems. The main problem is the skewness of these processes in which inattention leads to poor results in distinguishing the in-and out-of-control conditions. To overcome this problem, the root and the power transformation methods proposed by Niaki and Abbasi (2007) and Abbasi et al. (2011), respectively, are employed to almost eliminate the skewness. The following subsection is devoted to describe the method.

The root/power transformation methods
To eliminate the skewness involved in multi-attribute processes, Niaki and Abbasi (2007) first proposed the root transformation method, in which the skewness of the transformed attributes becomes almost zero. In this transformation, using the bisection method, a power within (0,1) is found for each attribute in the vector X = [X 1 , X 2 , . . ., X q ] ' to transform it into a new attribute vector Y¼X r ¼ X r 1 1 ; X r 2 2 ; . . . ; X r q q Â Ã 0 such that the transformed attribute vector contains attributes of almost zero skewness. However, this transformation method is only applicable for right-skewed distributions. Hence, in another work, Abbasi et al. (2011) extended the root transformation method and proposed a power transformation technique for left-skewed distributions as well. In the power transformation method, one first subtracts the observations from their minimum value and then finds a power within (0,1) for each attribute in the vector [X − Min(X)] = [(X 1 − Min(X 1 )), (X 2 − Min(X 2 )), . . ., (X q − Min(X q ))] ' such that attributes of almost zero skewness. Applying the above transformation methods, an attribute vector with zero skewness is obtained. Consequently, one may assume that the transformed vector approximately follows a multivariate normal distribution, and hence, the use of normal-based control charts is recommended. Moreover, to make sure the transformation works well, the Jarque and Bera (1987) test statistic given in Equation 1 is used to test the normality of the transformed attributes: where M is the number of observations, b 1 is the sample skewness, and b 2 is the sample kurtosis. This statistic follows an asymptotic chi-square distribution with 2 degrees of freedom and is used to test the null hypothesis that the transformed attributes are from normal distributions. Applying this transformation method, attributes with almost zero skewness that approximately follow multivariate normal distribution are obtained. Therefore, the well-known Hotelling T 2 control chart can be used to monitor the multi-attribute processes.
The multi-attribute T 2 control chart based on transformed data The Hotelling T 2 control chart was first introduced by Hotelling (1947) to come up with the correlation between the quality characteristics in monitoring multivariate processes. This control chart has been designed based on multivariate normal distribution. We use this control chart to monitor transformed attribute vector with approximately multivariate normal distribution.
Consider the observation vector X ij = [X ij1 , X ij2 , . . ., X ijq ] ' coming from a multivariate binomial distribution with in-control parameter vector [N ' , P 0 ' ] = [(n 1 , n 2 , . . ., n q ) ' , (p 1 , p 2 , . . ., p q ) ' ] and in-control correlation matrix Σ 0 , in which X ijk represents the jth observation in the ith subgroup of the kth quality characteristic. We apply the root/power transformation method in order to transform the vector X ij to Y ij = [Y ij1 , Y ij2 , . . ., Y ijq ] ' that approximately follows a multivariate normal distribution with the mean vector μ Y and the covariance matrix Cov Y .
For the transformed vector, the sample mean vector for each subgroup of size n is obtained using Then, for each subgroup the T 2 statistic is calculated by in which is the estimated in-control mean vector of the process. The statistic in Equation 3 is plotted on a T 2 control chart with an upper control limit (UCL) of where for large number of subgroups may be estimated by χ 2 α;q (see Lowry and Montgomery (1995) for more details). An out-of-control signal is detected when the T 2 statistic exceeds UCL.

Multivariate binomial process modeling and MLE derivation
In monitoring processes using the T 2 control chart, as long as the plotted points fall below UCL, the process is assumed in-control. However, when one or more points exceed UCL, the chart signals a change in the parameters of the process and the process is assumed out-of-control. The most important problem is that the chart signals after a considerable amount of time after the change point. Therefore, knowing when the change had actually occurred would substantially improve the diagnostic procedure.
The MLEs of the change points, designed for step changes and linear-trend disturbances are derived in the next two subsections.

The proposed change-point estimator for a step-change disturbance
In this section, we assume that the real change type presented in a multivariate binomial process is a step change and derive the MLE of a change point in the process fraction non-conforming. We also assume that the process is initially in-control with known parameter vector [N ' , P 0 ' ]. Following an unknown point in time τ SC (the step-change point), the process moves to an out-ofcontrol state with parameter vector [N ' , P 1 ' ], where P 1 is given in Equation 5 and δ is the vector of unknown magnitude of the change. An element δ > 1 in δ denotes an increase in the corresponding process fraction nonconforming, while δ < 1 denotes an improvement: The control chart genuinely (not false alarm) signals the above change with a delay at time S, when the T 2 S statistic exceeds the UCL of the T 2 control chart. Accordingly, the subgroup averages - . . . ; -Y S come from the out-of-control process. As a result, the MLE of the proposed change-point estimator, denoted byτ SC , is obtained as The MLE of the change pointτ SC is the value for τ sc that minimizes Equation 6.

The proposed change-point estimator for a linear-trend disturbance
Assuming again that the process is initially in-control with known parameter vector [N ' , P 0 ' ], in a linear-trend disturbance when the T 2 S statistic exceeds the UCL of the T 2 control chart, it is assumed that a change has occurred at an unknown point in time τ LT (the change point under linear trend). As a result the parameter vector of the out-of-control process becomes [N ' , P 1 ' ], where P 1 at different points in time takes a value based on Equation 7 and vector β contains the magnitude (or slope) of the linear-trend disturbance: Based on the above descriptions, the subgroup averages - . . . ; -Y S come from an out-of-control process. As a result, the MLE of the change point,τ LT , is the value of τ LT that minimizes Equation 8 and is obtained aŝ where In the next section, the performances of the proposed change-point estimators are evaluated and are compared using some simulation experiments. The MATLAB 7 software (Math Works, Inc., Natick, MA, USA) is utilized to perform all the programming works of this research.

Performance evaluation and comparison
Monte Carlo simulations are used to evaluate and to compare the performances of the proposed change-point estimators designed for step changes and linear-trend disturbances, when the actual changes are of step changes and linear-trend disturbances in the process fraction nonconforming of a multivariate binomial process.

Change-point estimators applied to step-change disturbances
The change point is simulated to occur at τ = 100. Thus, the first 100 subgroups are randomly generated from an in-control multivariate binomial distribution with parameter vector [N ' , P 0 ' ], using the normal to anything (NORTA) method (Cario and Nelson 1997;Niaki and Abbasi 2008). The in-control subgroup averages exceeding the UCL of the proposed control chart are false alarms. Consequently, these subgroups are discarded and are replaced with other in-control subgroups. This is replicated until all of the in-control subgroups are below UCL. Following subgroup 100, the parameter vector is changed based on a step-change disturbance given in Equation 5 where the subgroups are randomly generated from an out-of-control process with parameter vector [N ' , P 1 ' ]. The subgroup averages are plotted on the proposed control chart until a signal is generated, and at this time, the accuracy and the precision measures of the change-point estimators are computed.

Simulation experiment 1
Consider a multivariate binomial process with two correlated attributes, an in-control parameter vector [N ' , P 0 ' ] = [(20, 30) ' , (0.2, 0.15) ' ] and a correlation between the two attributes of 0.25. For subgroups i = 1, 2,. . ., 100, independent observations are first generated from the in-control process using the NORTA method. After subgroup 100, observations are generated from an out-of-control process with parameter vector [N ' , P 1 ' ]. Obtaining a signal from the control chart, the change points are estimated using Equations 6 and 8 for 1,000 simulation runs. Then, the averages, the standard deviations, and the mean square errors (MSEs) of the estimates in addition to the expected time when the control chart signals a change (E(T)) are calculated and are reported in Table 1. Table 1 also contains the precision performances of the estimated change points, Þare denoted by P0, P1, P2, P3, and P4, respectively. The estimated probabilities in parentheses correspond toτ SC .
Based on the results in Table 1, the change-point estimator designed for step changes outperforms the one designed for linear-trend disturbances in terms of both the accuracy and the precision measures, when the true change type is a step-change disturbance. In other words, for almost all of the shift magnitudes,τ SC is considerably closer to the actual change point thanτ LT is. Furthermore, MSEτ SC ð Þ is smaller than MSEτ LT ð Þ for all of the shift values.

Simulation experiment 2
Consider a multivariate binomial process with five correlated attributes. The sample size and the process fraction non-conforming of the process are [N ' , P 0 ' ] = [(25, 25, 25, 25, 25)  Using the NORTA method, the first 100 subgroups are randomly generated from the in-control process. Following subgroup 101, observations are generated from an out-of-control process based on the step-change disturbance given in Equation 5. The subgroup averages are then plotted on the proposed control chart until a signal is generated. Then, the averages, the standard deviations, the mean square errors, and the precision performances of the change-point estimates for 1,000 simulation runs are reported in Tables 2 and 3 for δ < 1 and δ > 1, respectively. Similar to Table 1, the precision performances ofτ SC in Tables 2 and 3 are shown in parentheses.
Again, the results in Tables 2 and 3 show thatτ SC outperformsτ LT for almost all shift magnitudes, sinceτ SC is closer to the true change point thanτ LT is. Moreover, MSEτ SC ð Þ is smaller than MSEτ LT ð Þ, andτ SC is more precise thanτ LT .
Based on the results obtained from the above two simulation experiments, one can conclude that when the real change type is a step-change disturbance, the change-point estimator designed for step changes outperforms the change-point estimator designed for lineartrend disturbances for almost all of the shift magnitudes and process dimensions.

Change-point estimators applied to linear-trend disturbances
Assuming a linear-trend disturbance, in this section, the performances of the two proposed change-point estimators designed for step-change and linear-trend disturbances are evaluated and are compared in multivariate binomial processes with two and five correlated quality characteristics. Considering the change point at τ = 100, the first 100 subgroups are randomly generated from an in-control multivariate binomial distribution with parameter vector [N ' , P 0 ' ], using the NORTA method. Following subgroup 100, the parameter vector is changed according to a linear-trend disturbance given in Equation 7, and as a result, subgroups are randomly generated from an outof-control process with parameter vector [N ' , P 1 ' ]. The subgroup averages are plotted on the proposed control chart until a signal is generated. At this time, the accuracy and the precision measures of the change-point estimators are computed.

Simulation experiment 3
Consider the process in simulation experiment 1. Using the NORTA method, the first 100 subgroup averages are generated from the in-control process. After subgroup 100, observations are generated from an out-of-control process. Following a signal by the proposed control chart, the change points are estimated using Equations 6 and 8 for 1,000 simulation runs. Then, the averages, the standard deviations, the mean square errors, and the precision performances of the change-point estimates in addition to the expected time when the control chart signals a change are summarized in Table 4.
The results in Table 4 show that when the real change type is a linear-trend disturbance, the proposed MLE for linear-trend disturbances outperforms the one designed for step-change disturbances in terms of both the accuracy and the precision for almost all values of β. That is, τ LT is closer to the true change point thanτ SC is.

Simulation experiment 4
Considering the process described in simulation experiment 2, 100 subgroups are first randomly generated from the in-control process. Then, the parameter vector is changed based on a linear-trend disturbance given in Equation 7. Consequently, observations are generated from an out-of-control process. Following a signal from the control chart, the averages, the standard deviations, the mean square errors, and the precision performances of the change-point estimates in 1,000 simulation runs are reported in Table 5.
Based on the results in Table 5, in the presence of a linear-trend disturbance,τ LT outperformsτ SC in terms of both the accuracy and the precision measures. In other words, in addition to the better precision performances of τ LT , the estimated change points byτ LT are closer to the real change point thanτ SC for almost all values of β. Moreover, MSEτ LT ð Þ is smaller than MSEτ SC ð Þ in most of the cases.
Finally, according to the results obtained from simulation experiments 3 and 4, one can conclude that when the real change type is a linear-trend disturbance, the change-point estimator designed for linear trend outperforms the change-point estimator designed for stepchange disturbances, regardless of the shift magnitudes and the process dimension.

Conclusions
When a control chart generates a signal due to an out-of-control condition, process engineers initiate a search to identify and to eliminate the root causes of process variation. However, the signaling time is different from the first time that a change manifests itself into a process (change point), and in most cases, there is a considerable time lag between them. Therefore, knowing the exact time of a change in a process would simplify the search procedure to improve the quality of the process by eliminating the causes of process variation.
In this paper, a new control chart was initially proposed based on a transformation method to monitor multi-attribute processes. The transformation method was used to eliminate the skewness of multivariate binomial processes. Then, the MLE of change points designed for both step changes and linear-trend disturbances were derived. Using Monte Carlo simulation experiments involving multivariate binomial processes with two and five attributes, the performances of the change-point estimators were evaluated by comparing them in the presence of real step-change and real lineartrend disturbances. The results showed that the changepoint estimator designed for step changes outperforms the change-point estimator designed for linear-trend disturbances when the true change type is a step-change disturbance, regardless of the shift magnitudes and process dimensions. In contrast, the change-point estimator designed for linear-trend disturbances outperforms the change-point estimator designed for step changes when the true change type is a linear-trend disturbance, regardless of the shift magnitudes and process dimensions.

Appendix
Appendix 1 Denoting the MLE of the change-point estimator designed for step changes byτ SC , the likelihood function for the transformed vector Y ij is (Johnson and Wichern 2007) where μ 0,Y , μ 1,Y , Cov 0,Y , and Cov 1,Y are in-control and out-of-control mean vector and covariance matrix of vector Y, respectively. The logarithm of the likelihood function when the constant term is shown by K is  Sharif University of Technology in 2012. His research interests are in the areas of applied statistics, quality engineering, and operations research.