Change Point Detection with Robust Control Chart

Abstract

Monitoring a process over time using a control chart allows quick detection of unusual states. In phase I, some historical process data, assumed to come from an in-control process, are used to construct the control limits. In Phase II, the process is monitored for an ongoing basis using control limits from Phase I. In Phase II, observations falling outside the control limits or unusual patterns of observations signal that the process has shifted from in-control process settings. Such signals trigger a search for assignable cause and, if the cause is found, corrective action will be implemented to prevent its recurrence. The purpose of this paper is to introduce a new methodology appropriate for constructing a robust control chart when a nonnormal or a contaminated data that may arise in phase I state. Through extensive Monte Carlo simulations, we examine the behaviors and performances of the proposed MM robust control chart when there is a process shift in mean.

1. Introduction

Statistical process control (SPC) concepts and methods have become very significant in the manufacturing and process industries. Their goal is to monitor the performance of a process over time in order to justify whether or not the process is remaining in a “state of statistical control.” This state of control is said to occur if certain process or product variables remain near to their desired values and the only source of variation is “common-cause” variation, that is, the variation which affects the process all the time and is essentially inevitable within the current process. Shewhart charts are used to monitor key product variables in order to detect the occurrence of any event having a “special” or “assignable” cause. By discovering assignable causes, long-term improvements in the process and in product quality can be accomplished by eliminating the causes or improving the process or its operating procedures.

Detecting one or more change points in a batch of observations has attracted substantial investigation in the statistical, engineering, and econometric literature. Assuming that there is an ordered sequence of observations, usually, but not necessarily, taken at equally spaced times, there is a change point between two successive observations if their statistical distributions are different. Between change points, the distributions are usually considered to be identical. In practice, recognizing when a process has changed would simplify the search for the special cause. If the time of the change could be identified, process engineers would have a smaller search window within which to look for the special cause. Consequently, the special cause can be determined more quickly, and necessary actions needed to improve quality can be carried out sooner. In this paper, we will analyze the efficiency of a change point estimator in process mean [1] for each of Shewhart 𝑋, Median, and the proposed MM control chart once issues a signal. The derivation of the change point estimator that is shown in the appendix is in virtue of Hinkley [2]. Hinkley discussed the asymptotic properties of the estimator. Whenever each of Shewhart 𝑋, Median, or the proposed MM chart signals that a special cause is present. The estimator provides practitioners with a useful estimate of the time of the process change. In Section 2, we will introduce a model for a step change in the location of a process. We consider a step change for a process mean occurs when the mean suddenly changes its value and then remains unchanged again until corrective action has been taken. On the basis of this step-change model, we adopt the estimator of the time of the process change when the corresponding chart does signal. In Section 5, we analyze the performances of each chart by means of Monte Carlo simulation.

2. Process Step-Change Model

Suppose that the process is initially in control, with observations coming from a Normal distribution with a known mean of 𝜇₀ and a known standard deviation of 𝜎₀. Even so, after an unknown point in time 𝑇 (known as the process change point), the process location changes from 𝜇₀ to 𝜇₁=𝜇₀+𝛿𝜎₀/√𝑛, where 𝑛 is the subgroup size and 𝛿 is the unknown magnitude of the change. Assuming also that once this step change in the process location occurs, the process remains at the new level of 𝜇₁, until the special cause has been identified and removed. We let 𝑋_𝑇  be the first subgroup average to exceed a control limit and that this signal is not a false alarm. Hence, 𝑋₁,𝑋₂,…,𝑋_𝜏 are the subgroup averages that come from the in-control process, while 𝑋_𝜏+1,𝑋_𝜏+2,…,𝑋_𝑇 are from the changed process.

3. Definitions

To illustrate, we concentrate on robust estimates for the simple location-scale model, by letting 𝑥₁,…,𝑥_𝑛 be 𝑛 observations on the real line satisfying𝑥_𝑖=𝜇+𝜎𝜀_𝑖,𝑖=1,…,𝑛, (3.1) where 𝜀_𝑖 are independent and identically distributed observations with variance equal to 1. We are interested in estimating 𝜇 and the scale 𝜎 which is a nuisance parameter. We consider 𝑀-location estimates which is proposed by Huber [3]. He defined 𝜇_𝑛 as the solution of an estimating equation of the form^𝑛_𝑖=1𝜓𝑥_𝑖−𝜇_𝑛𝜎_𝑛=0, (3.2) where 𝜎_𝑛 is a robust estimate of the residuals scale, and 𝜓∶𝑅→𝑅 is a bounded, nondecreasing, and odd real function. We focus on 𝜓 function which is continuous and differentiable influence function given in (3.2) (see [4]) 𝜓_𝑐(𝑢)(𝑢)𝑢=sign𝑐|𝑢|𝑝,if≤0.8𝑐,₄𝑢𝑐|𝑢|𝑝,if0.8𝑐<≤𝑐,₄(1),if|𝑢|≥𝑐, (3.3) where 𝑐>0 is a user-chosen tuning constant, and 𝑝₄(𝑢)=38.4−175𝑢+300𝑢²−225𝑢³+62.5𝑢⁴ (see [5]), for other choices of smooth functions 𝜓.

The scale estimate 𝜎_𝑛 in (3.2) is an 𝑆-estimate of scale (see [6]) which is defined as follows. Let 𝜌∶𝑅→𝑅₊ be a bounded, continuous, and even function satisfying 𝜌(0)=0 and let 𝑏∈(0,1). The 𝑆-scale 𝜎_𝑛 is defined by𝜎_𝑛=inf_𝑡∈𝑅𝑠_𝑛(𝑡), (3.4) where, for each 𝑡∈𝑅, 𝑠_𝑛(𝑡) is the solution of 1𝑛^𝑛𝑖=1𝜌𝑥_𝑖−𝑡𝑠_𝑛(𝑡)=𝑏. (3.5) Indeed, associated with this family are the 𝑆-location estimates 𝜇_𝑛 is given by𝜇_𝑛=arginf_𝑡∈𝑅𝑠_𝑛(𝑡). (3.6) Beaton and Tukey [7] proposed a family of functions 𝜌_𝑑 given by𝜌_𝑑(𝑢)=3𝑢𝑑²𝑢−3𝑑⁴+𝑢𝑑⁶|𝑢||𝑢|,if≤𝑑,1,if>𝑑, (3.7) where the tuning constant 𝑑 is positive. According to Yohai [8], these 𝑀-location estimates obtained with an 𝑆-scale estimate are called MM-location estimates. Specifically, the estimates 𝜇_𝑛, 𝜎_𝑛, and 𝜇_𝑛 solve the following system of equations:^𝑛_𝑖=1𝜓𝑥_𝑖−𝜇_𝑛𝜎_𝑛1=0,𝑛^𝑛𝑖=1𝜌𝑥_𝑖−𝜇_𝑛𝜎_𝑛1−𝑏=0,𝑛^𝑛𝑖=1𝜌^′𝑥_𝑖−𝜇_𝑛𝜎_𝑛=0. (3.8) Let 𝑑=1.548 for 𝜌_𝑑 in (3.7), 𝑏=0.5 in (3.5), and 𝑐=1.525 for 𝜓_𝑐 in (3.3), which yields a location estimate 𝜇_𝑛 with 50% breakdown point and 95% efficiency when the errors have a normal distribution.

3.1. Sample Median

Sample median has been used in early process control charts as it is insensitive to behavior in the tails of the distribution. However, under the normal distribution, the efficiency of the sample median drops off rapidly towards its asymptotic value of 0.64 as sample size increases. For a random sample of size 𝑛 observations 𝑋₁,𝑋₂,…,𝑋_𝑛, the sample median, denoted by MD, is defined as follows:𝑋MD=𝑋^((𝑛+1)/2),if𝑛isodd,_(𝑛/2)+𝑋_((𝑛/2)+1)2,if𝑛iseven. (3.9)

The interest of using the sample median, MD, is that it is easy to determine, requires only the middle values to calculate, can be used when a distribution is skewed, is not affected by outliers, and has a maximal 50 percent breakdown point. Moreover, its gross-error sensitivity is low and as the sample size 𝑛 increases, the variance of the MD decreases as 1/𝑛, but the maximum bias does not change. Hence, the bias is the property of importance for large sample sizes, and MD is the estimator possesses the smallest maximum bias for a given proportion of contamination 𝜀. Huber [3] showed that it minimizes the maximum asymptotic bias over contamination neighborhoods. As opposed to that, the disadvantages for the sample median, MD, are its difficulty to handle in mathematical equations, nonutilizing all available values, and being misleading when the distributions come from a long tail distribution as it might sometimes discard some useful information (see [9, 10]). However, the sample median has become as a good general purpose estimator and is generally considered as an alternative average to the sample mean especially whenever outliers might present in the distribution.

3.2. Median Absolute Deviation from the Sample Median

The median absolute deviation from the sample median, denoted as MAD, is a more robust scale estimator than the standard deviation. The MAD was first introduced by Hampel [11] who attributed it to Gauss. It is simple and easy to compute and mainly used in detecting outliers in a data. The estimate is often used as an initial value for the computation of more efficient robust estimators. Let us denote 𝑋₁,𝑋₂,…,𝑋_𝑛 as a random sample of size 𝑛 observations with sample median MD. MAD possesses the following properties: (i)it has a maximal 50 percent breakdown point which is twice as the IQR;(ii)in the case of the standard normal distribution, 𝐹, the influence function of the MAD estimator, 𝐼𝐹(𝑋;MAD,Φ), is a step function that takes on two values. This 𝐼𝐹 is bounded by the sharpest possible bound among all scale estimators. With regard to the optimality properties of MAD, Martin and Zamar [12] established expressions for the maximum asymptotic bias of 𝑀-estimates of scale over contamination neighborhood as a function of the fraction of contamination and show that the similar strong results are obtained in terms of maximum asymptotic bias for MAD as with the MD.

3.3. Control Limits

In this paper, we consider the general equations for constructing control limits (see [13]). Thus, with a robust location estimator 𝑇 and the corresponding scale estimator 𝑆, the control limits are given byUCL=𝑇+3𝑆𝐴√𝑛,CL=𝑇,LCL=𝑇−3𝑆𝐴√𝑛. (3.10) The constant 𝐴 in (3.10) is determined in such a way that 𝑆/𝐴 is an unbiased estimator of the scale parameter. The most commonly used control charts are Shewhart 𝑋 charts using the sample range. For 𝑋 charts using the sample standard deviation, the 𝑇 in (3.10) is the sample mean 𝑋 and 𝑆 is the sample standard deviation with 𝐴=𝑐₄.

Shewhart modified these control limits using rational subgroups (see [14]) in which 𝑚 rational subgroups with each of size 𝑛 are taken. According to Shewhart's suggestion, these subgroups are formed in order that the between-groups variability is maximized while the within-group variation is minimized. In this view, then1𝑋=𝑚^𝑚𝑖=1𝑋_𝑖,1𝑆=𝑚^𝑚𝑖=1𝑆_𝑖, (3.11) where 𝑋_𝑖 and 𝑆_𝑖 are the subgroup mean and standard deviation, respectively. Each of these estimates is an unbiased estimate of the corresponding parameter, then the control limits using rational subgroups are3UCL=𝑇+𝑆𝐴√𝑛,3CL=𝑇,LCL=𝑇−𝑆𝐴√𝑛. (3.12)

In practice, the control limits are the average of the control limits for the 𝑚 subgroups. In the case of 𝑋 chart for which we employ the sample range as the scale parameter, it is estimated with the average range computed by averaging over the 𝑚 subgroups𝑇=1𝑋=𝑚^𝑚𝑖=1𝑋_𝑖,𝑆=1𝑅=𝑚^𝑚𝑖=1𝑅_𝑖. (3.13) Then, the control limits are defined asUCL=3𝑋+𝑅𝐴√𝑛,CL=𝑋,LCL=3𝑋−𝑅𝐴√𝑛. (3.14)

For constructing the control charts under a normal distribution using the robust estimators, we will determine the appropriate constant 𝐴 for the desired estimators through computer simulations. To illustrate, a sample of size 𝑛 was taken from 𝑁(0,1). The constant 𝐴 was computed by averaging over 100,000 repetitions. Here, for instance, if we consider using 𝑆 as scale estimator, then over 100,000 repetitions, we expect 𝐸(𝑆/𝐴)=𝜎 for any 𝜎. Table 1 exhibits the constant 𝐴 such that 𝐸(𝑠𝑐𝑎𝑙𝑒𝑒𝑠𝑡𝑖𝑚𝑎𝑡𝑜𝑟/𝐴)=𝜎. It can be seen that if scale estimator is the sample range, the simulated values (𝐸₁) agree closely with the standard tabled values. 𝐸₂ and 𝐸₃ are the corresponding estimates for the Median and the proposed MM charts.

Table 1 

Simulated values of 𝐴 for different scale estimators.

Thus, in the same way, the estimators for location and scale for Median-MAD chart are given by𝑇=Median𝑥_𝑖,𝑆=MAD=Median𝑥_𝑖−𝑇, (3.15) while the estimators for the proposed MM chart are given by𝑇=𝜇_𝑛=arginf_𝑡∈𝑅𝑠_𝑛(𝑡),𝑆=𝜎_𝑛=inf_𝑡∈𝑅𝑠_𝑛(𝑡). (3.16) In this study, three different estimators under investigation are as follows:

𝐸₁: 𝑇 = Sample Mean; 𝑆 = Range,

𝐸₂: 𝑇 = Median; 𝑆 = MAD = Median{|𝑥_𝑖−𝑇|},

𝐸₃: 𝑇 = 𝜇_𝑛=arginf_𝑡∈𝑅𝑠_𝑛(𝑡); 𝑆 = 𝜎_𝑛=inf_𝑡∈𝑅𝑠_𝑛(𝑡).

3.4. Confidence Regions (Confidence Set)

One of the benchmarks for assessing the performance of a control chart is to construct a confidence region for the time of the process change. The use of confidence region on the change point is that it will suggest practitioners with useful starting points for searching their process log books and records for the special cause. This will provide the practitioners a “search window” for the special cause and aid in quicker identification. Hence, practitioners can then take necessary action for the special cause sooner in order to improve quality as well as to reduce process downtime.

Basically, we will incorporate the likelihood function to obtain a confidence region for the process change point. The confidence region approaches in the statistics literature involve the likelihood function that relies on asymptotic theory (see [15]). In process monitoring, there are relatively small time intervals between the process change point and the time of the control chart signal. Thus, approximations based on asymptotic theory may not be appropriate.

Box and Cox [16] proposed a method involving the log likelihood function for constructing a possibly noncontiguous confidence region (also called a confidence set) on a parameter. Their approach can be used to build a confidence set (CS) for the process change point using the log likelihood function having the formCS={𝑡∶ln𝐿(𝜏)−ln𝐿(𝜏)<𝐷}. (3.17) Here, ln𝐿(𝜏) is the maximum value of the log likelihood function; 𝜏 is the MLE of 𝜏 (i.e., the value of 𝑡 that maximizes the log likelihood function), whereln𝐿(𝑡)=−𝑇ln√2𝜋𝜎_𝑥−𝑛2𝜎²₀^𝑇𝑖=1𝑋²_𝑖−2𝜇₀^𝑡_𝑖=1𝑋_𝑖+𝑡𝜇²₀−(𝑇−𝑡)𝑋_𝑇,𝑡² (3.18) is the value of the log likelihood function at 𝑡. We let 𝑘₁ and 𝑘₂ represent constants determined by the subgroup averages. Thus, ln𝐿(𝑡) can be expressed asln𝐿(𝑡)=𝑘₁−𝑛2𝜎²₀𝑘₂−(𝑇−𝑡)𝑋_𝑇,𝜏−𝜇₀². (3.19)

Box and Cox [16] proposed using 𝐷=(1/2)𝜒²_1,𝛼 to obtain a 100(1−𝛼)% confidence region. Siegmund [17] used asymptotic theory to develop a 100(1−𝛼)% confidence set for the change point of a normal process mean based on the log likelihood function. He proposed using the value𝐷=−ln1−(1−𝛼)^1/2. (3.20)

By means of Monte Carlo simulation, we study the confidence sets obtained with nominal confidence coefficient of 1−𝛼=0.90. In accordance with Box and Cox [16] and Siegmund [17], the 𝐷 values for a 90% confidence set are 𝐷=1.353 and 𝐷=2.97, respectively. It was observed that 𝐷=1.353 value suggested by Box and Cox [16] provides a 90% coverage for a value of 𝛿 between 2.0 and 3.0, while Siegmund’s [17] of 𝐷=2.97 provides at least 90% coverage for 𝛿≥1.0. By trial and error, the value of 𝐷=3.7 provides at least 90% coverage for 𝛿≥0.5.

4. Methodology

In order to analyze the performance of the control charts, we consider using Shewhart 𝑋, Median, and our proposed MM control chart. When a control chart signals that suggest a process change has occurred, the change point estimator (see the appendix) is then applied to the data to estimate the time of the change at which we need to find the value of 𝑡 in the range 0≤𝑡<𝑇 which maximizes 𝐶_𝑡=(𝑇−𝜏)(𝑋_𝑇,𝑡−𝜇₀)². The reverse cumulative average, 𝑋_𝑇,𝜏=(𝑇−𝜏)⁻¹∑^𝑇_𝑖=𝜏+1𝑋_𝑖, is the overall average of the 𝑇−𝑡 most recent subgroups for which the value of 𝑡 maximizing the 𝐶_𝑡 values is our estimator of the last subgroup from in-control process.

5. Simulation Study

We will now analyze the performance of the change point estimator adopting the three control charts through the simulation study. Assuming that the process is initially in control, with observations coming from a normal distribution with a known mean of 𝜇₀ and a known standard deviation of 𝜎₀. However, after an unknown point in time 𝜏 (known as the process change point), the process location changes from 𝜇₀ to 𝜇₁=𝜇₀±𝛿𝜎₀/√𝑛, where 𝛿 is the unknown magnitude of the change. We also assume that once this step change in the process location occurs, the process remains at the new level of 𝜇₁, until the special cause has been identified and removed.

To illustrate, the data are generated under different settings of distributions. In addition to the normal distribution, two alternative distributional forms are considered. They are contaminated model (Case 2 of (5.2)) and Slash distribution. Under different types of distributions, for each run, the data consist of 𝑚=30 subgroups of size 𝑛=5 are used to construct the control limits and summary statistics are calculated. In order to assess the performance of the corresponding chart, observations for 1 to 100 are generated from standard normal distribution. Then, starting from subgroup 101, observations were randomly generated from a normal distribution with mean 𝛿 and standard deviation 1 until each of Shewhart 𝑋, Median, and the proposed MM control chart produces a signal. The procedure was repeated a total of 10,000 times for each of the values of magnitude that was studied, namely 𝛿=1.0, 2.0, and 3.0. For each simulation run, the change point estimate was computed. Subsequently, the average of the estimates of 𝜏 for the 10,000 simulation runs was computed along with its standard error, expected length, and coverage probability.

For analyzing the outlier model, we modified the contaminated model by Davis and Adams [18]. Specifically, the in-control conditions with contaminated data values are determined by generating a random number from a Uniform (0, 1) distribution and the corresponding frequency of contaminated data, 0≤𝛽≤1. A Uniform (0, 1) random probabilistic value 𝑝_𝑖𝑗 is generated for observations 𝑗 of sample 𝑖, 𝑥_𝑖𝑗. Let 𝐼_(𝑎,𝑏)(𝑝) represent an indicator function with𝐼_(𝑎,𝑏)𝑝=1,if𝑎<𝑝<𝑏,0,otherwise. (5.1) A random observation 𝑥_𝑖𝑗 in the simulated data is described by (5.2). Different changes for each of process states are illustrated in the expression. For the purpose of comparisons, the frequencies of contamination of 𝛽=0.05 and 𝛽=0.10 with 𝐶=9.0 are considered for the purpose of creating some disturbances in the data 𝑥_𝑖𝑗∼𝑁𝜇+𝐶,𝜎²𝐼_(0,𝛽)𝑝_𝑖𝑗+𝑁𝜇,𝜎²𝐼_[𝛽,1)𝑝_𝑖𝑗,for𝑖=1,2,…,10,000,𝑗=1,2,…,𝑛. (5.2)

Case 1. In-control (no contaminated data): 𝜇=0, 𝜎²=1, 𝛽=0, 𝐶=0.

Case 2. In-control (contaminated data): 𝜇=0, 𝜎²=1, 𝛽=0.05,0.10, 𝐶=9.0.

As with the computation of sizes of confidence sets obtained using a specific value of 𝐷, for each control chart and magnitude of change studied, a step change in the normal process mean was simulated following 𝜏=100. The confidence set estimator was applied following a signal from the corresponding control charts considered. The size of the confidence set was recorded as well as whether the confidence set covered the true process change point of 𝜏=100. This procedure was repeated for a total of 𝑁=10,000 simulation runs for each of values 𝛿 considered. The proportion of the 10,000 runs that covered the true process change point was also determined. This was reflected by the resulting estimates of the coverage probabilities obtained by specifying the 𝐷 value to be 2.97, such that it provides at least 90% coverage for 𝛿≥1.0, which has been discussed in Section 4. The results are tabulated along with the average sizes of confidence sets. For a given coverage probability, a smaller confidence set is preferred, so that process engineers can more narrowly focus their search for the special cause. In general, it is presumed that the increase in the magnitude of shift will be followed by the increase in the corresponding coverage probability.

From Tables 2, 3, 4, and 5, we can see that the performances of the three charts in terms of coverage probability are quite similar especially when the process data come from a normal distribution. Generally, the Median chart and the proposed MM chart perform better in most of the cases (particularly for a larger proportion of contamination). For instance, for a change of magnitude 𝛿=3.0, with 𝛽=0.05, 𝐶=9.0, the coverage probability utilizing 𝑋 chart is 0.588, while the coverage probability for Median chart and the proposed MM chart are 0.626 and 0.632, respectively.

Table 2 

Simulation result: estimates of expected length, average change point, standard error, average size of confidence set and coverage probability for change point of 𝜏=100 (normal distribution).

Table 3 

Simulation result: estimates of expected length, average change point, standard error, average size of confidence set, and coverage probability for change point of 𝜏=100 (contaminated distribution, 𝛽=0.05, 𝐶=9).

Table 4 

Simulation result: estimates of expected length, average change point, standard error, average size of confidence set, and coverage probability for change point of 𝜏=100 (contaminated distribution, 𝛽=0.10, 𝐶=9).

Table 5 

Simulation result: estimates of expected length, average change point, standard error, average size of confidence, set and coverage probability for change point of 𝜏=100 (slash distribution).

In Tables 2–5, averages of change point estimates 𝜏 are also tabulated for various sizes of change in the process mean together with its corresponding standard error estimates for a normal case setting. As the actual change point for the simulation was at time 100, the average estimated time of the process change, 𝜏, should possibly be close to 100. With 𝑋 chart, we see that when the process step change of standardized magnitude 𝛿=1, the average estimated time of the process change was 100.00, which is fairly close to the actual change point of 100. While for a standardized process location change of size 𝛿=2, the average estimated time of the change is 99.60. Meanwhile, when 𝛿=3, the average estimated time of the change is 99.47. Hence, on average, the change point estimate of the time of the process change is considerably close to the actual time of the change, regardless of the magnitude of the change.

By the same taken, with the Median chart, for the process step change of standardized magnitude 𝛿=1, the average estimated time of the process change was 99.70, which is also close to the actual change point of 100.00. As for a standardized process location change of size 𝛿=2, the average estimated time of the change is 99.49. And when 𝛿=3, the average estimated time of the change is 99.57.

Lastly, when the process is monitored with our proposed robust MM chart, by and large, the change point estimate of the time of the process change is fairly close to the actual time of the change, regardless of the magnitude of the change. In the case when the process step change of standardized magnitude 𝛿=1, the average estimated time of the process change was 99.90. It turns out that for a standardized process location change of size 𝛿=2, the average estimated time of the change is 99.42. As with 𝛿=3, the average estimated time of the change is 99.59. Overall, we could say that for all types of charts under study, the change point estimator of the time of the process change is able to detect the change point considerably close to the actual time of the change, irrespective of the magnitude of the change.

Another benchmark of evaluating the control chart is by examining the expected length of the signal. This is the expected time at which the control chart signals a change in the process mean that is supposed to occur at time 100. It is generally perceived that Shewhart 𝑋 control chart might issue a signal of a change in a process mean a considerable amount of time after the change in the process mean actually occurred. Thus, estimating the time of process change with the time when the control chart indeed issues a signal would lead to an unfavorably biased estimate. As a consequence, probably a misleading estimate of the time of the process changes. This bias is in virtue of the potentially large delay in generating a signal from the control chart. Hence, the criterion for evaluating the performance of the control chart is how quick the chart would signal (Expected Length). Tables 2–5 and Figures 1–4 summarize the performance in terms of expected length for the three charts.

Figure 1 

Simulation result: expected length for change point of 𝜏=100 (normal distribution).

Figure 2 

Simulation result: expected length for change point of 𝜏=100 (contaminated distribution, 𝛽=0.05, 𝐶=9).

Figure 3 

Simulation result: expected length for change point of 𝜏=100 (contaminated distribution, 𝛽=0.10, 𝐶=9).

Figure 4 

Simulation result: expected length for change point of 𝜏=100 (slash distribution).

In a normal distribution situation, for a step change in the process mean of magnitude 𝛿=1, it is easy to see from Figure 1 that the expected length for a Shewhart 𝑋 chart is 57. The result seems to suggest that 𝑋 chart is the best compared to the Median chart and the proposed MM chart for which each needs 157 and 98, respectively, when the shift is small. The situation improves as the magnitude of shift 𝛿 increases to 2 or 3. All types of charts considered here seem to be relatively comparable and perform quite closely.

On the other hand, the Shewhart 𝑋 chart is inferior in an outlier model. With 𝛽=0.10 and 𝐶=9, according to Table 4, both the Median chart and the proposed MM chart appear to be better than Shewhart 𝑋 chart for different magnitudes of shifts. The similar situation arises when a very heavy-tailed distribution (Slash distribution) is considered. Again, the Median chart and the proposed MM chart outperformed Shewhart 𝑋 chart with respect to the expected length. It can be observed from Table 5 that the differences are quite apparent, whereby Shewhart 𝑋 chart requires the expected length of 24 and yet the Median chart and the proposed MM chart both demand about 3 subgroups before detecting the first signal when the magnitude of shift is 1.

We now turn to evaluate the observed frequency in which the estimates of the time of the step were within 𝑚 observations of the actual time of the change, for 𝑚=0,1,2,…,10. The results are tabulated in Tables 6, 7, 8, 9, 10, 11, 12, 13, and 14. This provides an indication of the precision of the estimator by means of the three different charts. The proportion of the 10,000 runs where the estimated time of the change was within ±𝑚 of the actual change is expected to be increase in size as 𝑚 increases. Referring to Tables 6–14, we observed that the precision increases with the increases of 𝑚 for each 𝛿 value. Let us first focus our attention to a normal setting, when the process step change of magnitude 𝛿=2. Monitoring the process using traditional Shewhart 𝑋 chart identified correctly the change point in 60.49% of the trials. It was within one observation of the actual change point in 83.23% of the trials, and within two observations of the actual change point in 91.39% of the trials. Turning to the Median chart, which is shown in Table 7, the chart detected correctly the change point in 59.59% of the trials. It was within one observation of the actual change point in 83.46% of the trials, and within two observations of the actual change point in 91.66% of the trials. It then follows that, based on Table 8, our proposed MM chart located accurately the change point in 59.61% of the trials. It was within one observation of the actual change point in 82.88% of the trials, and within two observations of the actual change point in 91.25% of the trials. All types of charts considered here seem to be comparable and performed quite equally.

Table 6 

Precision of the estimator when used with 𝑋-Bar chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (normal distribution).

Table 7 

Precision of the estimator when used with median chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (normal distribution).

Table 8 

Precision of the estimator when used with proposed MM chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (normal distribution).

Table 9 

Precision of the estimator when used with 𝑋-Bar chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (contaminated distribution, 𝛽=0.10, 𝐶=9).

Table 10 

Precision of the estimator when used with median chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (contaminated distribution, 𝛽=0.10, 𝐶=9).

Table 11 

Precision of the estimator when used with proposed MM chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (contaminated distribution, 𝛽=0.10, 𝐶=9).

Table 12 

Precision of the estimator when used with 𝑋-Bar chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (slash distribution).

Table 13 

Precision of the estimator when used with median chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (slash distribution).

Table 14 

Precision of the estimator when used with proposed MM chart for different magnitudes of process change (magnitude of shift, 𝛿) based on 10,000 trials. Subgroup size 𝑛=5, change point 𝜏=100 (slash distribution).

Next, we observe the situation under outlier model setting, with 𝛽=0.10 and 𝐶=9. Consider again when the process step change of magnitude is 𝛿=2, for step change of this magnitude, the Shewhart 𝑋 chart estimator exactly identified the time of the change in just 23.65% of the trials and was within one (two) observation of the time of the actual process change in 45.03% (56.69%) of the trials. As for the Median chart, we can notice that, of the 10,000 simulation trials conducted for 𝛿=2, 24.99% of those simulation trials identified the change point precisely. It was in 46.92% and 58.11% of the trials that the change point was estimated to be within ±1 and ±2, respectively, from the actual time of the process change. The results of the study also indicate that 25.25% of those simulation trials identified the change point correctly for the proposed MM chart. In 47.53% of the trials, the estimate was within ±1 observation, and in 58.61% of the trials, the estimate was within ±2 observations. Overall, the procedure using the proposed MM chart seems to perform slightly better than the other two charts in this respect.

Finally, we focus on nonnormal data without outliers’ situation, the Median chart and the proposed MM chart generally lead to shorter control limits than the traditional Shewhart 𝑋 chart. Here, we just limit ourselves to the study of Slash distribution (see [13]) which is a very heavy-tailed distribution. When the process step change of magnitude 𝛿=2, monitoring the process using traditional 𝑋 chart identified correctly the change point in 5.31% of the trials. It was within one observation of the actual change point in 14.44% of the trials, and within two observations of the actual change point in 20.64% of the trials. For the Median chart which is shown by Table 13, the chart was able to detect correctly the change point in 5.45% of the trials. It was within one observation of the actual change point in 12.08% of the trials, and within two observations of the actual change point in 17.95% of the trials. It then follows that, based on Table 14, our proposed MM chart discovered accurately the change point in 5.49% of the trials. It was within one observation of the actual change point in 12.28% of the trials, and within two observations of the actual change point in 18.02% of the trials. On the whole, we would say that the Median, and the proposed MM charts perform better and more consistently than the 𝑋 chart in this setting.

6. Conclusions

Control charts are used to detect whether or not a process has changed. When a control chart signals indicate that a process has changed, practitioners must initiate a search for the special cause. However, given a signal from a control chart, practitioners generally do not know what caused the process situation to change or when the process has changed. Identifying the time of the process change would simplify the seeking of the special cause. In the event that the practitioners knew when the process changed, the search would simply be reduced for discovering what aspect of the process changed at that time. As a result, practitioners would increase their chances of identifying the special cause more correctly and quickly. Subsequently, This allows them to take the appropriate actions immediately to improve the quality.

In this paper, monitoring processes in the presence of data contamination and under nonnormal setting are of primary concern. We have applied an estimator that is useful for identifying the change point of a step change in normal process mean, nonnormal, and when contamination may exist. We have discussed the performance of the change point estimator and other criteria when they are monitored by Shewhart 𝑋, the Median and the proposed MM control charts. The results show that the proposed MM robust control chart consistently performed well in a range of situations. It provides a useful and much better alternative in using the time of the signal from the conventional 𝑋 control chart. Although the proposed MM chart and the Median chart are comparable under nonnormal and contamination situation; the Median chart becomes worse in the event of normal setting. The performance of our proposed MM control chart has good properties in the aspect of expected length and coverage probability for contamination data and that arise from a heavy-tailed distribution functions for moderate sample sizes. The proposed MM chart compares favorably with traditional Shewhart 𝑋 control chart in normal setting especially when magnitudes of shift are 2 and 3. It is interesting to note that the proposed robust MM chart is more efficient than the Shewhart 𝑋 chart when the process distribution function has a heavy-tailed distribution.

The philosophy of the proposed robust MM control chart is more in keeping with the desire to provide robust limits in the face of nonnormal distribution situation or that outlier(s) may arise in data collection. The proposed robust MM charting methodology gives better performance than the traditional Shewhart 𝑋 chart if the underlying distribution of chance cause is nonnormal or contaminated. This is usually desirable feature for any control chart to be applied in the industry.

Appendix

In this appendix, we show the derivation of the maximum likelihood estimator (MLE) of 𝜏, the process location change point. For maximum-likelihood estimation techniques, see [1]. Denoting the MLE of the change point 𝜏 as 𝜏, let subgroup averages be 𝑥₁,…,𝑥_𝑇, the MLE of 𝜏 is the value of 𝜏 that maximizes the likelihood function or, equivalently, its logarithm. It is shown that logarithm of the likelihood function is log𝐿𝜏,𝜇₁∣𝑥𝑛=−2𝜎²₀^𝜏𝑖=1𝑥_𝑖−𝜇₀²+^𝑇_𝑖=𝜏+1𝑥_𝑖−𝜇₁²𝑛=−2𝜎²₀^𝑇𝑖=1𝑥²_𝑖−2𝜇₀^𝜏_𝑖=1𝑥_𝑖−𝜏𝜇²₀−2𝜇₁^𝑇_𝑖=𝜏+1𝑥_𝑖+(𝑇−𝜏)𝜇²₁. (A.1) We can see that there are two unknowns in the log-likelihood function: 𝜏 and 𝜇₁. If the change point 𝜏 was known, the MLE of 𝜇₁ would be 𝜇₁=𝑋_𝑇,𝜏=(𝑇−𝜏)⁻¹∑^𝑇_𝑖=𝜏+1𝑋_𝑖, the average of the 𝑇−𝜏 most recent subgroup averages. Substituting this back into (A.1), we obtainlog𝐿(𝜏∣𝑥)𝑛=−2𝜎²₀^𝑇𝑖=1𝑥²_𝑖−2𝜇₀^𝜏_𝑖=1𝑥_𝑖+𝜏𝜇²₀−(𝑇−𝜏)𝑋²_𝑇,𝜏. (A.2) It is easy to verify that this is equivalent tolog𝐿(𝜏∣𝑥)𝑛=−2𝜎²₀^𝑇𝑖=1𝑥²_𝑖−2𝜇₀^𝜏_𝑖=1𝑥_𝑖+𝜏𝜇²₀−(𝑇−𝜏)𝑋_𝑇,𝜏−𝜇₀². (A.3) In what follows, the value of 𝜏 that maximizes the log-likelihood function is 𝜏=argmax_𝑡(𝑇−𝜏)𝑋_𝑇,𝑡−𝜇₀²=argmax_𝑡𝐶_𝑡, (A.4) where 𝐶_𝑡=(𝑇−𝜏)(𝑋_𝑇,𝑡−𝜇₀)², that is, 𝜏 is the value of 𝜏 in the range 0≤𝑡<𝑇 which maximizes 𝐶_𝑡=(𝑇−𝜏)𝑋_𝑇,𝑡−𝜇₀². (A.5)

Acknowledgment

The authors wish to thank the anonymous referee(s) for very constructive comments and suggestions that improved the paper.