A control chart for a general Gaussian process

Journal of Statistical Planning and Inference 70 (1998) 19-34

ELSEVIER

journal of statisticalplanning and inference

A control chart for a general Gaussian process Robert Lund a,*, W.J. Padgett b, Lynne Seymour a a Department of Statistk~, The University of Georgia, Athens, GA 30602-1952, USA b Department of Statistics, The University of South Carolina, Columbia, SC 29208, USA

Received 3 November 1997

Abstract It has been repeatedly demonstrated that X-bar quality control charts perform poorly when the process subgroups being monitored are correlated. In this paper, we propose and investigate the performance of a control chart that accounts for subgroup correlations in a general Gaussian process. The time-series innovations algorithm is used to construct the desired chart from a set of one-step ahead predictions and prediction variances. The chart is applicable in both stationary and nonstationary settings. A simulation study shows that this 'innovations' chart performs as a traditional X-bar chart even when the correlation structure of the process must be estimated from a small number of subgroups. The innovations chart is then used to study a data set of motor shaft diameters which has correlated subgroups. The results here show that erroneous conclusions can be reached if subgroup correlations are ignored. (~ 1998 Elsevier Science B.V. All rights reserved. Kevwords: Quality control; X-bar chart; Innovations algorithm; ARMA model

1. Introduction It is common practice in statistical process control (SPC) to utilize Shewhart X - b a r control charts for monitoring the mean and variation o f a process. These charts rely heavily on the assumption that sample statistics between different process subgroups are independent. In contrast, it has recently been demonstrated that strong autocorrelations between process subgroups frequently exist. Such autocorrelations arise for a number o f reasons. One reason is that, because o f current automated measurement and recording technology, subgroup samples may be taken with a high frequency, with consecutive samples being similar in nature and hence statistically correlated. Other examples o f correlated subgroups occur when items made by a worker exhibit similar characteristics

* Corresponding author. 0378-3758/98/$19.00 (~) 1998 Elsevier Science B.V. All rights reserved. PII S0378-3758(97)00180-8

20

R. Lund et al./ Journal oJ' Statistical Plannin 9 and InJerence 70 (1998) 19-34

due to the way a machine is handled, when a process shows seasonal patterns due to materials or weather, or when the alertness of a worker changes over time. Subgroup correlations markedly affect the performance of a standard X-bar control chart. For correlated subgroups, Padgett et al. (1992) show that the true probability of at least one sample mean exceeding the control limits differs greatly from those calculated under the assumption of independent subgroups. There is now a large body of literature investigating control charts for correlated processes. Contributions include Alwan and Roberts (1988), Alwan (1992), Hahn (1989), Harris and Ross (1991), Johnson and Bagshaw (1974), Lin and Adams (1994), Maragah and Woodall (1992), Montgomery and Mastrangelo (1991 ), Superville and Adams (1994), Tseng and Adams (1994), VanBrackle and Reynolds (1997), Vander Wiel (1994), Vasilopoulos and Stamboulis (1978), Wardell et al. (1994), Woodall and Faltin (1993) and Yaschin (1993). Many of the above authors investigated the case where o n e measurement is taken from the process at each sampling time. A typical approach to monitoring correlated processes when only one measurement is taken at each sampling time is to fit an ARIMA time series model to the process observations and then monitor one-step ahead forecasts with a control chart for single measurements (cf. Alwan and Roberts, 1988; Alwan, 1992; Montgomery and Mastrangelo, 1991). If the ARIMA model orders and parameters have been identified correctly, then the residuals, or differences between forecasted and observed values, should be approximately uncorrelated. A Shewhart 'individuals' chart (X chart), CUSUM chart, or other chart for single observations is then applied to the residuals to monitor the process. Willemain and Runger (1994) use a different approach based on random run lengths. In many applications, multiple observations are taken from the process at each sampling time. It is the purpose of this paper to develop and investigate a control chart for this general situation when the subgroups are correlated. We will consider the general case where a subgroup of nt >~1 measurements is sampled from a process at time t (nt is allowed to vary with t). The methods presented below require only knowledge or estimation of the correlation structure of the process being monitored; in particular, no stationarity assumptions are made. Stationarity of subgroups is violated in periodic processes (cf. Spurrier and Thombs, 1990) and in economic applications where the random walk and other nonstationary unit root time series models are typically used to describe process fluctuations. The rest of this paper proceeds as follows. In Section 2, model notation and assumptions are clarified and the control chart is developed. Section 3 develops several examples of the methodology mathematically. Section 4 then explores control chart performance via simulation; here, it is demonstrated that the proposed chart makes a vast improvement over the standard X-bar chart that neglects subgroup correlations. The chart is also shown to perform well when the correlation structure of the process subgroups must be estimated. Section 5 applies the results of this paper to the motor shaft diameter data in Devor et al. (1992). Section 6 concludes the paper with some comments.

R. Lund et al./Journal of Statistical Plannin 9 and Inference 70 (1998) 19 34

21

2. The control chart

Suppose that at each time t, 1 ~
X,/= #t + ~.t/,

(2.1)

where {ktt} is a Gauss/an process with finite second moments and the subgroup 'errors' {etj} are assumed to be independent and identically distributed (i/d) mean zero Gauss/an random variables with Var(etj)--a~ 2 E (0, oo). One interpretation of Eq. (2.1) is that the 'mean setting' of the subgroup at time t is /zt; this mean setting is allowed to fluctuate randomly with t. Conditional on /~t, the subgroup sample {At:., l~
(2.2)

For properties of the subgroups, Eq. (2.1) and E[pt] = 0 give E[Xt)] = O. The covariance between items in subgroups t and s is Cov(Xtj,X~.t)= a(t,s)+ a~ 0[t=.~•i-1]; hence, this model has correlated subgroups whenever a(t,s)¢O for some t ¢ s . Notice that the time t subgroup {Xtj, l<~j<~nt} is only conditionally independent. In fact, when j ¢ l, Corr(Xtj,Xtt) = [a(t, t) + a,z]-Ia(t, t). The subgroup average at time t is lit

)(t =n~ -1 S )(i/.

(2.3)

j=l

Setting nt

(2.4) j= 1

we obtain )(t = #t + gt

(2.5)

from Eq. (2.1). The covariance structure of {)(t} is obtained from Eqs. (2.1)-(2.5) as

K(t,s) = Cov(P~t,Xs) = a(t,s) + n71 a,ffD[t=~]. Notice that {Xt} is stationary when a(t,s) depends on in t.

(2.6)

It- sl

only and nt is constant

22

R. Lund et al./ Journal of Statistical Plannin9 and lnjerence 70 (1998) 19 34

To develop a general control chart, we will make use of a set of one-step ahead predictions and prediction errors. Let (2.7)

Xt = E [ X t IX1 . . . . , ) ( t - I ]

be the best mean-squared error predictor of )(, based upon )(j,... ,)?t-~ and define its mean-squared prediction error as v2_, = E[0~, - )7,) 2 IX, ..... )(t-,].

(2.8)

Since {J~t} is Gaussian, conditional expectations are linear in the predicting variables and the Projection Theorem (cf. Brockwell and Davis, 1991) shows that ) ~ , - )(t is uncorrelated (and hence independent) of){/ for 1 < ~ i ~ t - 1. Thus, Eq. (2.8) becomes v21 = E[O~t-JT"t)2]. The innovations algorithm (see Brockwell and Davis, 1991, Ch. 5) can be applied to recursively compute )~, and v2 in t provided the covariance matrix K, = [h'(i,j)]~'.j= 1 is invertible for each n>~ 1. This is done through the innovations representation _~ t Xt+l = ~ O t j ( ~ ' t - j + l - X t j=l

j+l),

t~>l

(2.9)

with ) ~ = 0. The quantities 0,j and v2 can be computed recursively from v02= ~c(1, 1 ) and O u _ k = ( v ~ ) -j

~c(t+l,k+l)-

~

Ok, k_tOt, t_j

,

k = O ..... t - l ,

j=O

(2.10) v2 = K ( t + l , t + l )

t I

,,2 t_jl; ,2 -- ~ Ut, j /=0

for t~> 1. Eq. (2.10) is solved in the order v2; 01~, v2; 022, 02~,v2; . . . . A general and simplistic sufficient condition for i¢, to be invertible for each n ~> 1 follows from Proposition 5,1.1 of Brockwell and Davis (1991) when {#t} is stationary. Letting 7(h)=cr(t + h,t) for natural numbers h, this condition is that 7 ( 0 ) > 0 and 7(h)---+ 0 as h---+ oo. For a control chart, we will use the standardized prediction residuals .Vt - X, Zt-- - Ut 1

(2.11)

Since {)?t} is Gaussian, the prediction error sequence { ) ~ , - )?t} is uncorrelated and Gaussian; hence, {Zt} is a sequence of iid standard normal random variables. A large Z, in absolute value indicates a large deviation between )(t and what was predicted by )~, and is associated with an out of control process. In many applications, the control limits ±3 are used and the process is called out of control at time t if IZtl>3; the

R. Lund et al./Journal of Statistical Planninq and lnlerence 70 (1998) 19- 34

23

probability of this happening due to chance variation when the process is actually in control is approximately 0.0027. In practice, the parameters governing the correlation structure of {/it} and a~ are typically unknown. These parameters can be estimated with maximum likelihood from a small number of subgroup measurements under the innovations framework. Letting tl denote all unknown model parameters and using the notation Zt(~/) to explicitly denote the dependence of the standardized residual Zt on model parameters, the model likelihood function, denoted by L(q; {)~t}), has the innovations form (cf. Brockwell and Davis, 1991, Ch. 8) L(,_/; {)[~t})=(2g) -k''2

,:,-- v[_,~

exp - 5 ,:,

.

(2.12)

An estimate of q, denoted by 0, can be computed by minimizing the scaled negative log-likelihood function H(r/; {Xr}) defined by k

k

H(q_; {)(t})= ~ ln(v)_, ) + ~ ZT(r/) t=l

(2.13)

t=l

in the argument ~/. Once 0 is found, an estimate of the covariance function ~c(.,-) in Eq. (2.6) can be obtained by plugging 1~ in for q. The estimated covariance structure can in turn be used in Eqs. (2.9) and (2.10) to compute an estimated version of {Z,(_~)} which we denote by {Z,(_~)}. The control limits ±3 can then be applied to {Zt(~)} to assess the status of the process. In general, {Zt(_~)} will not be an uncorrelated series; however, as k-+ :x~, the maximum likelihood estimate ~ frequently converges in probability to q. When this happens, Zt(0) converges in probability to Zt(q) as k-+ oc for each fixed t and {Zt(_~)} becomes an asymptotically uncorrelated series. In practice, the Section 4 simulations will show that this estimation procedure yields a control chart that performs reasonably well for a variety of correlation structures when between 20 and 100 subgroups are used for estimation purposes.

3. Examples In this section, we mathematically explore some specific cases of the Section 2 model.

3.1. The first-order autore,qression Consider the AR(1 ) difference equation ]Jr=if)Fit-1 + A t ,

(3.1)

where Iq~l<1 is assumed for a stationary and causal solution {#t}. In Eq. (3.1), {At} is an iid mean zero Gaussian random sequence with Var(At) = a~ that is independent

24

R. Lund et al. I Journal of Statistical Planning and In[erence 70 (19981 19 34

of {cO}" For simplicity, we take nt - n. From 7(h) = Cov(~t+h,/~t) = (t -- 4)2)-1cr2411< ' the correlation between items in subgroups t and s, t ¢ s , can be shown to be Corr(X,j, X,;) =

24,-,t + (1 - 4)2) 2"

Since 7(h)--,0 as h---,oc, K, is invertible for each n~>l. Note that Eqs. (3.1) and (2.5) can be regarded as Kalman filter state and observation equations respectively. Using Eqs. (2.5) and (3.11 in Eq. (2.7), we obtain )~, = 4)fit-t where f , _ l = E [ I t t - i [)(I ..... )(t- 1] is the best mean-squared error predictor of #t from -gl ..... )(~-1. Further manipulations with Eqs. (2.5) and (3.1) provide L't_ I = 4)2pt2 1 --0 -2 -1- n - l ~ 2,

(3.2)

where Pt2 = E [ ( f t - # , ) z [)?t . . . . . )~t] = E [ ( f t _ ~,)2]. Hence, the control chart is easily constructed from {Pt2} and {fit}. Classical Kalman filter recursions (cf. Shumway, 1988) provide fit = 4)fit--1 ~- K t ( L - Ofit 1),

(3.31

where Kt is the Kalman gain K,

412p2 t_] + 0-2 n - ' a , 2 + 412pt2-1 + 4 "

(3.4)

To compute Kt and fit recursively in t, one uses Eqs. (3.3) and (3.4) and recursively updates p2 with "~ n-l~2(412pt2_l jr_ 0-2) • + 412p2, + 0-2; P," = i

(3.5)

the scheme is initiated with ~0 = 0 and Poe = 0-2(1 - 4)21 ~. Alternatively, the above equations can be converted into recursions for )~t and v~. Using )~t = 41f,-I in Eq. (3.3) gives

)~,+, = 41(1 -

(3.6)

K , ) , V , + 4)K,X,

for t>~ 1, which is similar to exponential smoothing. Substituting Eq. (3.2) in Eq. (3.5) gives /...2

0.2 q_ (1 _{_412)(n 1~:2)

412[n- 1~71 2 t,2_ '

(3.7)

These recursions are initiated with P~l = 0 and vo2 =0-~i2(1 - 412)--1 + .--Ia,2. 3.2. The r a n d o m w a l k

Consider the random walk Itt=ltt

I +At,

t>~l

(3.8)

R. Lund et al./Journal of Statistical Plannin.q and InJerence 70 (1998) 19- 34

25

with/to = 0 . Here, {At} is an iid mean zero Gaussian sequence with V a r ( A t ) - a A that is independent of {ctj }. The covariance structure of {/~t} is well known as a ( i , j ) = ad min(i,j) and Eq. (2.6) gives

tr(i,j) = o"2 min(i,j) + n71~ 2 D[i=/].

(3.9)

Hence, corr(X,,

o"ffmin(t, s) 2 + n,

l) = (t 2 +

2)Jj2

(3.10)

for t Cs. Note that det(lc,)>~(aff) " det(M,) where the (i,j)th entry of M . is min(i,j). Using d e t ( M , ) = 1 shows that d e t 0 ¢ , ) > 0 for all n~> 1; hence, K. is invertible for all n~>l. When nt varies with t, 0q and v2 can be obtained recursively from Eq. (2.10). When nt -= n, the innovation recursions can be simplified. A reduction in the computations can be made if we apply the innovations algorithm to {,gl, 112,113.... } instead of directly to {)?1,)?2 .... }. Here, {Yt, t~>2} is defined via Y t - X t - X t - i . Eqs. (2.5) and (3.8) give

Y,=At÷ff,-~,,

1.

(3.11)

Since the linear span of {kl,)~2 . . . . . )(t} is the same as the linear span of {)~l, Yz, Y3.... , ~} for each t ~>1, one has 12,+l = E[Yt+j I)(i, II2.... , lit] =2~t+l - ) ( t

(3.12)

for t ~>2. Letting {17z, t/>0} denote the one-step ahead mean-squared prediction errors ^ for {)(1, Y2, Y3.... }, we obtain )(t = 0 and r/Oz =UO2 - - - O ' 2 ÷ H from - n l~e_ 21,

r/z=aZ+2n-1~2

1~2. })Z and r/~ are obtained

[n-1~:212

The crucial observation is that {Yt, t/>2} can be written as the first-order moving average Yr =Zt + 0 K - l ,

(3.13)

where {Zt} is mean zero Gaussian white noise with V a r ( Z t ) = a 2. Thus, {)?t, Y?. . . . . } is a 'one-dependent' sequence. The parameters 0 6 ( - 1 , 0 ) and a } can be obtained by equating moments in Eqs. (3.11) and (3.13). This yields ( 1 + 02)o 2 = ~42+ 2n ~7 and 0a 2 - n 1~2 which has the solution o = _(¢2 + 2n-1~2 az2 =

2[n-1~2] 2

(rrd + 2n 'a, 2) - x/ad(a ff + 4n-'a,2)"

26

R. Lund et al./Journal of Statistical Plannin# and lnjerence 70 (1998) 19-34

Applying Eq. (2.10) to {)(I,Y2, Y3.... } yields (here, the 0~/ refer to the prediction coefficients for {)(j, Y2 .... }) 0 t i = 0 , j = 2 . . . . . t, ~/~=0.A + 2 n - I

Oil=

(3.14)

qr 21[n-I~],

:2 - q t --2~n-to~2~2 tl ,:] ,

for t>~2. Substituting Eq. (3.14) in Eq. (2.9) now gives f,+l-O,t(Y,

- Y,)

for t~>2 which can be combined with Eq. (3.12) to give the desired recursion

~,+, =(1 +

0,)2, - 0,~),

(3.15)

for t ~>2. It remains to compute )~2 and v2 for t >~ 1. Straightforward computations give ffA2 1~.2•A1 ~ 2 = ~2 +/?

and

when t>~2. Using r, = ~ - k ,

v 2 -- 2a 2 + n - ' q : 2

(0"2) 2 aft + n - l ~ : 2'

~ and Eq. (3.12) gives ~,,~='7~ as required.

We again notice that Eq. (3.15) is similar to exponential smoothing; however, Ott depends on t. From Eq. (3.14), one can show that 0tl-+ 0 as t-+ ~ ; hence, asymptotically, predictions for the random walk {&} obey an exponentially weighted moving average relation. We also note that Eq. (3.7) with ~ = 1 is identical to the last equation in Eq. (3.14). 3.3. Causal linear processes

A flexible class of stationary models for {/~t}, which includes the causal ARMA class, can be written in the form t~t=

~ ~kAt k,

k=0

where ~-~k~0 ]~k]< vc and {At} is an iid Gaussian sequence, independent of {tk/}, with V a r ( A r ) - a 2. Then {&} is a mean zero stationary series with covariance function 7(h)=a(t

+ h,t)=a2A

~ ~k~k+h.

k=0

From ~ - ~ 0 ]~kll if ~ k ¢ 0 for some k>~0. When nt=-n, {J(I,)(2 .... } will be stationary with covariance function ~-(t, s) = Cov(X, X,. ) = 7(t - s) + n - 1a2 Drt s]. In this case, computations can be simplified by employing the Durbin-Levinson algorithm (see Brockwell and Davis, 1991, Ch. 5) rather than the innovations algorithm.

R. Lund et al./Journal of Statistical Plannin9 and InJerence 70 (1998) 19 34

27

4. A simulation study of the chart This section investigates the performance of the innovations control chart via simulation. To study a variety of correlation structures, we will consider AR(1), ARMA (1, 1), and random walk models for {#t}. The results will be compared to the case where a standard X-bar chart is applied to correlated subgroups without accounting for subgroup correlations. To summarize control chart performance, we will use the ~k-risk of Padgett et al. (1992). The ~k-risk is defined as the probability of at least one out of control signal in k total subgroups when the process is actually in control:

l<~t<~k].

~k-P[[Z~(_q)I>3 for at least one t satisfying

(4.1)

We note that ~k can be interpreted as a Type I error probability. For a Gaussian process where parameters are known, ~k has a binomial form and does not depend on q ~k = 1 - (1 - :q )~,

(4.2)

where ~1 = P[IZ[ > 3 ] ~ 0.0026998 is the probability that the absolute value of a standard normal deviate Z exceeds 3. Of course, :tk--, 1 as k--~ ec. Table 1 below lists several values of :~k. Our first two simulations, summarized in Tables 2 and 3 below, pertain to the stationary AR(1) model for {#t} in Eq. (3.1). In this and all subsequent simulations, we will examine the subgroup lengths k = 2 0 , 3 0 , 5 0 , and 100, and subgroup sizes Table 1 ~k-risks k

~k

20 30 50 100

0.05263 0.07790 0.12643 0.23688

Table 2 AR(1) :q-risks: ~b= ~, I o-~; 2 = ~,o" I A 2= 3 Number of subgroups k

n= 3 n- 4 n--5 n=6

20 0.0160 0.0435 0.0225 0.0520 0.0170 0.0575 0.0155 0.0480

30 0.0485 0.0830 0.0390 0.0730 0.0505 0.0855 0.0525 0.0840

50 0.1060 0.1255 0.1035 0.1410 0.0880 0.1210 0.1005 0.1185

100 0.2250 0.2365 0.2230 0.2470 0.2535 0.2650 0.2155 0.2355

28

R. Lurid et aL /Journal of Statistical Plannino and lnlerence 70 (1998) 19 34 Table 3 A R ( I ) ~k-risks: 4 > = - - i

2

1

2

3

N u m b e r of subgroups k 20

30

50

100

0.0225

0.0395

0.0945

0.2225

0.0465

0.0770

0.1120

0.2290

n- 4

0.0215 0.0505

0.0470 0.0750

0.1055 0.1265

0.2125 0.2215

n- 5

0.0165 0.0605 0.0155

0.0495 0.0755 0.0615

0.0965 0.1335 0.0980

0.2300 0.2410 0.2170

0.0540

0.0910

0.1225

0.2285

n- 3

n- 6

nt =-3,4,5, and 6. These choices of nt and k are relatively consistent with sample sizes encountered in applications. The autoregressive parameter 4) was chosen as 0.5 and -0.5 in Tables 2 and 3, respectively; this IqSI is viewed as injecting a moderate degree of correlation into the subgroups. The values of aA2 and cr2 were selected to make Var()?t)= 1 when n = 1. In Tables 2 and 3, the upper risk is the sample cq.risk for the control chart using the estimated standardized prediction residuals {Zt(~)} and the lower risk is the sample ~k-risk for the control chart using the standardized prediction residuals {Zt(_q)) where the value of r/ is assumed known. All sample risks are computed from 2000 simulations; specifically, the sample 2k-risks are the proportion of simulations where a standardized prediction residual, Zt(~) or Zt(_q), exceeds 3 in absolute value for at least one t satisfying 1 ~
R. Lund et al./Journal o f Statistical Planning and Inference 70 (1998) 19 34

29

Table 4 2 ~ I 2 A R M A ( 1 , 1 ) :~k-risks: qS= ~,~ = ~,0-~

, 2 1,0-a = 1

N u m b e r o f subgroups k

n= 3

20 0.0210

30 0.0470

50 0.1040

100 0.2170

n= 4

0.0600 0.0205

0.0775 0.0555

0.1175 0.1005

0.2395 0.2265

0.0460

0.0800

0.1255

0.2545

n= 5

0.0245 0.0535

0.0510 0.0760

0.0940 0.1145

0.2015 0.2390

n= 6

0.0245 0.0490

0.0425 0.0780

0.1005 0.1190

0.2125 0.2435

Table 5 1 ,,~p= g2 ' 0-~: 2 - I, 0-2 A R M A ( 1, 1 ) 2/,-risks: ~b = - 5'

1

N u m b e r o f subgroups k

n= 3 n= 4 n- 5 n- 6

20 0.0260 0.0610

30 0.0555 0.0765

50 0.1045 0.1190

100 0.2195 0.2460

0.0170 0.0440 0.0185

0.0495 0.0750 0.0505

0.1095 0.1300 0.0975

0.2290 0.2520 0.2025

0.0525 0.0225 0.0505

0.0775 0.0435 0.0840

0.1090 0.0980 0.1260

0.2230 0.2115 0.2405

correlations of {Zt(_~)} were computed in each simulation. A 95% large-sample confidence interval for the mean residual was 3.95 x 1 0 - 4 -t- 1.96 X 10-3; hence, the standardized prediction residuals appear to have a zero mean. In contrast, a 95% largesample confidence interval for the lag one correlation was -0.04957 + 1.81 x 10 -3 which supports the conjectured negative lag one correlation. Of course, if the correlation structure of {Zt(~)} for small k was known, this information could be used to 'further standardize' {Zt(~)} into an lid sequence; unfortunately, we know of no such small sample results. Despite the slight bias due to estimation, the sample C~k-riskscomputed from {Z,(~)} are much closer to the values reported in Table 1 than the values reported in Padgett et al. (1992) where subgroup correlations were ignored. Furthermore, Tables 2 and 3 show that as the series length k increases, the two sample ek-risks become closer. As discussed in Section 2, this is because as k --, oc, ~ - ~ 17 in probability and Zt(r~) ~ ~(r/) in probability for each fixed t. In a similar format to Tables 2 and 3, Tables 4 and 5 report sample C~k-risks aggregated from 2000 simulations with the A R M A ( 1 , 1 ) model t~t -- ~ t - - 1

=At + OAt-l,

30

R. Lund et al./ Journal of Statistical Plannin# and InJerence 70 (1998) 19-34

Table 6 Random Walk u~.-risks: a~=~ 1, ~2 = 1 Number of subgroups k

n- 3 n= 4 n- 5 n: 6

20 0.0190 0.0515 0.0225 0.0450 0.0270 0.0640 0.0220 0.0475

30 0.0495 0.0810 0.0450 0.0755 0.0495 0.0820 0.0465 0.0880

50 0.0995 0.1295 0.1030 0.1205 0.1015 0.1270 0.0950 0.1250

lOO 0.2225 0.2235 0.2165 0.2440 0.2120 0.2295 0.2230 0.2380

where {At} is mean zero Gaussian white noise with Var(At) = a~ that is independent o f {etj}. The white noise variances were taken as a,2:= a ~2 = 1 in both tables; the A R M A ( 1 , 1) parameters are chosen as q 5--- 72 and 0 = 7 I in Table 4 and ~b = - ½ and 0 = 2 in Table 5. We note that Tables 4 and 5 have a similar structure to Tables 2 and 3 and omit further discussion. Table 6 presents an analogous simulation for the nonstationary random walk model for {/Jr} in Section 3.2. The random walk has strong correlations between neighboring subgroups. The parameters for this simulation were selected as a 2 = ¢7A 2 = 1; for these parameter choices, Eq. (3.10) gives C o 1 T ( Y t j , X t _ l , l ) : ( t - 1)/[(t + 1)t] I/2. Table 6 shows a similar structure to Tables 2 5. Overall, the above simulations show that the innovations control chart behaves as a standard X - b a r chart with uncorrelated subgroups when the correlation structure o f the subgroups is known. When the subgroup correlation structure must be estimated from a sample o f subgroups, the sample ~k-risks are slightly smaller than their true values but, nonetheless, improve greatly over values reported in Padgett et al. (1992) when subgroup correlations are ignored. O f course, the bias in the ~k-risks will be asymptotically negligible when subgroup correlation parameters are consistently estimated.

5. An application In this section, we apply the innovations chart to the motor shaft diameter data in Devor et al. (1992), (pp. 186-188). This data has process subgroups o f size n t - 5 observed at k = 60 sampling times separated by 30rain each. Fig. 1 plots the sample mean motor shaft diameter at each sampling time; the traditional X - b a r chart upper and lower control limits, )~ :t: A2R, computed assuming independent subgroups, are included for reference. Notice that process sample means fall outside o f the control boundaries for subgroups 28 and 32; hence, a traditional X - b a r chart would signal 'out o f control' for this process twice.

R. Lund et al./Journal of Statistical Plannin,q and InJerence 70 (1998) 19 34

65

i

i

i

i

i

i

i

i

i

i

31

i

t

60 o,,

o..

....

oooo.,o,,,,,o,,,,,,,,,,,,o.o.,,,

,,,o.,,,,,,,oo,

g 55 E t~ co 45 40 35

I

I

I

I

5

10

15

20

I

I

25 30 35 Subgroup Number

I

I

[

I

40

45

50

55

60

Fig. 1. Sample means of the motor shaft diameters.

An exploratory analysis of the subgroup data suggests that a stationary model for {/~t} is appropriate; hence, we will use an ARMA(p,q) model to describe {#t} P

q

#, - ~ ~ l t , _ i = A , + ~ OsA, i, i l

(5.1)

i--I

where {A,} is white noise with V a r ( A t ) - cr2 that is independent of {~tj}. Minimizing the scaled negative log likelihood in Eq. (2.13) for various p and q, and using AICC to select p and q, we select the AR(1) model in Section 3.1 for {#t}. An overall sample mean of 48.793, computed as the average of all observations, was first subtracted from the measurements to obtain mean zero data. The A R ( I ) model estimates obtained were q~l = 0.741 -4-0.106, and d2 = 2.693. An error variance of ff~ = 39.366 and a negative log likelihood of 159.956 were obtained in the model fit; an AICC statistic of 326.340 was achieved. The error margin attached to q~l is one standard error and was computed by inverting the observed information matrix associated with the negative log likelihood in Eq. (2.12). The estimate of q~i and its standard error strongly suggest that 4)1 is positive; hence, we conclude that the process subgroups are positively correlated. For confirmation of the adequacy of the AR(1) model, we will examine {Z~(_~)}. Fig. 2 plots Z~(0) against t along with 4-3 control limits. If the fitted model is good, then {Zt(~)} should have correlation properties that are similar to those of an uncorrelated sequence. Fig. 3 plots the sample autocorrelation function of {Zt(_~)} over the first 20 lags; the dashed lines are 95% confidence limits for iid sequences and suggest that our fitted model is adequate. The estimated standardized residuals also passed several commonly applied tests for normality (cf. Brockwell and Davis, 1991, Section 9.4); hence, the Gaussian A R ( I ) model for {/tt} appears to be reasonable. The Zt(O) values in Fig. 2 all lie inside the control limits i 3 ; hence, contrary to a traditional X-bar chart, the innovations chart suggests that the process is never out of control. Finally, we note that Zt(0) follows a similar pattern in t as the sample means in Fig. 1, but have a smaller variability; hence, one can directly see the effects of accounting for subgroup correlations.

R Lund et al./Journal of Statistical Planning and InJerence 70 (1998) 19 34

32

I . . . . . .

I ,

I

I

I

I

. . . . . . . . . . . . . . . . . .

I

, o , . , , .

I

I

I

. . . . . . . . . . . . . . .

I

I

•

J

J

2 O)

8 1 (B = 0

0

"6 ~5-I 13.

-2 -3

0 ° o

. . . .

o

°.

° . . . .

. . . . .

•

I

I

I

I

5

10

15

20

°

°

..

. . . . . . . . .

I

I

o

°

°

I

25 30 35 Subgroup Number

. . . . . . . .

.

°

. . . . . .

°

° . o ,

I

I

l

I

40

45

50

55

60

Fig. 2. Prediction Z-scores of the motor shaft diameters

I

I

I

I

I

I

I

I

I

08 06

-~ 04 o

8 02 5 ,< 0

.............. 2 ............

.....2 2

-02 -04

I

I

I

I

1

I

/

I

I

2

4

6

8

10 Lag

12

14

16

18

20

Fig 3 Sample ACF of the motor shaft prediction Z-scores.

6. Comments If the normality assumptions are violated, then conditional expectations may not be linear in the predicting variables and linear predictors may not coincide with the best mean-squared-error predictors. We refer the interested reader to Brockwell and Davis (1991, Section 13.4) for comments pertaining to when linear prediction is 'optimal' in terms of mean-squared error. O f course, even if linear prediction is optimal, the 4 3 control limits may need to be adjusted as the distribution of Zt(_q) may not be standard normal if {)(,} is not Gaussian. Since Corr(Xt;,Xtz) = [a(t, t) + a~]-l a(t, t) >10 when j ¢ I, the time t subgroup {Xth 1 ~
R. Lund et al./Journal of Statistical Plannin9 and InJbrence 70 (1998) 19-34

33

Finally, one may wish to use a small pilot sample of subgroups rather than the whole data set to estimate subgroup correlation parameters we are not advocating using all subgroups in all cases. Of course, the typical purpose of a control chart is typically to detect departures from some form of baseline behavior.

Acknowledgements Robert Lund's research was supported by National Science Foundation grant DMS9703838; W.J. Padgett was supported in part by National Science Foundation grant DMS-9503104. We are grateful to Ishwar Basawa and Richard Davis for useful conversations on the subject matter.

References Alwan, L.C., 1992. Effects of autocorrelation on control chart performance. Commun. Statist. Theory Meth. 21, 1025-1049. AIwan, L.C., Roberts, H.V., 1988. Time-series modeling for statistical process control. J. Bus. Econom. Statist. 6, 87-95. Brockwell, P.J., Davis, R.A., 1991. Time Series: Theory and Methods, 2nd ed. Springer, New York. Devor, R.E., Chang, T.H., Sutherland, J.W., 1992. Statistical Quality Design and Control. Macmillan, New York. Hahn, G.J., 1989. Statistics-aided manufacturing: a look into the future. Amer. Statist. 43, 74 79. Harris, T.J., Ross, W.H., 1991. Statistical process control procedures for correlated observations. Canad. J. Chem. Eng. 69, 48 57. Johnson, R.A., Bagshaw, M., 1974. The effect of serial correlation on the performance of CUSUM tests. Yechnometrics 16, 103-112. Lin, W., Adams, B.M., 1994. The use of combined EWMA-Shewhart control chart with forecast-based monitoring schemes. Technical Report, Department of Management Science and Statistics, The University of Alabama, Tuscaloosa, AL. Maragah, H., Woodall, W.H., 1992. The effect of autocorrelation on the retrospective X-chart. J. Statist. Comput. Simulation 40, 29-42. Montgomery, D.C., Mastrangelo, C.M., 1991. Some statistical process control methods for autocorrelated data. J. Quality Technol. 23, 179 204. Padgett, C.S., Thombs, L.A., Padgett, W.J., 1992. On the :~-risks for Shewhart control charts. Commun. Statist. Simulation Comput. 21, 1125 1147. Shumway, R.H., 1988. Applied Statistical Time Series Analysis. Prentice-Hall, Englewood Cliffs, NJ. Spurrier, J.D., Thombs, L,A., 1990. Control charts for detecting cyclical behavior. Technometrics 32, 163 171. Superville, C.R., Adams, B.M., 1994. An evaluation of forecast-based quality control schemes. Commun. Statist. Simulation Comput. 23, 645-661. Tseng, S., Adams, B.M., I994. Monitoring autocorrelated processes with an exponentially weighted moving average forecast. J. Statist. Comput. Simulation 50, 187 195. VanBrackle, L.N., Reynolds, M.R., 1997. EWMA and CUSUM control charts in the presence of correlation. Commun. Statist. Simulation Comput. 26, 979-1008. Vander Wiel, S., 1994. Statistical process monitoring with integrated moving average noise. Technical Report, AT&T Laboratories, Murray Hilt, NJ. Vasilopoulos, A.V., Stamboulis, A.P., 1978. Modification of control chart limits in the presence of data correlation. J. Quality Technol. I0, 20-30. Wardell, D.G., Moskowitz, H., Plante, R.D., 1994. Run length distributions of special-cause control charts for correlated processes. Technometrics 36, 3 17.

34

R. Lund et al./Journal of Statistical Planning and lnJerence 70 (1998) 19-34

Willemain, T.R., Runger, G.C., 1994. Statistical process control using level crossings, J. Statist. Comput. Simulation 51, 7 20. Woodall, W.H., Faltin, F.W., 1993. Autocorrelated data and SPC. ASQC Statist. Div. Newslett. 13, 18-21. Yaschin, E., 1993. Performance of CUSUM control schemes for serially correlated observations. Technometrics 35, 37-52.