LetXbe the measurement with distributionF from the system represent-ing a characteristic of a subject of interest andX1:::Xnis a random sample drawn from distribution F and we choose a statisticT =t(X1:::Xn) that has mean t and variance t2. The Shewhart control chart set t as the centre line and placed three standard deviations above and below the centre line as
UCL=t+ 3t
LCL=t;3t (2.1)
If statistic T follows a normal or Gaussian distribution, the limits of the chart will cover the values of T in the long run with probability 0:9973. In practice, it is hard to believe that we know t andt. Therefore, we need to estimate them. Control charts are calculated based on a historical record of observations such asmsubgroups ofnsample and points outside the control limits are excluded and the revised control limits are calculated. Let 0 and
0 be the corresponding estimates. The estimated control limits are UCL=t0+ 3t0
LCL=t0;3t0 (2.2)
We say that the laboratory measurement system is stable when the mea-surement values of T are fell within the limits. On the other hand, a run is rejected and the measurement system is said to be out of control when its measurement T value exceeds the control limits.
Suppose that the measurement variableX follows the normal distribution N(2) and we consider the X-chart with sample mean X = n1 Pni=1Xi
as the test statistic. Since X has normal distributionN(n2). The control limits of the X-chart in form of (2.1) are:
UCL=+ 3pn
LCL=;3pn (2.3)
By letting0 and0 as, respectively, the sample mean and sample standard deviation based on the historical record of observations, The estimated con-trol limits in form of (2.2) are
UCL=0+ 3p0n
LCL=0;3p0n (2.4)
The diculty in process control in laboratories is that due the nan-cial cost and analyst's time there is usually no available historical record of enough observations to compute the accurate estimates of 0 and 0. Sta-tistical inferences may help in checking if a control chart computed from a limited data represents for the distribution of an in-control process.
Suppose that we have a set of observations x1i:::xni, i = 1:::m and we compute estimates0 and0 to form a control chart with limits of (2.3).
The interest now is to test if the measurement system is in statistical control.
That is to test the following hypothesis:
H0 :t;3t =t0;3t0 and t+ 3t =t0+ 3t0 (2.5) which is also equivalent to test the hypothesis:
H0 :t =t0t =t0 (2.6) The validation problem is that we have a random sample X1::::Xk
drawn from a distribution F for testing hypothesis (2.5) or (2.6).
3. Condence Interval Based Tests for Laboratory Control Chart
Validation
Suppose that, calculated from the historical data, we have estimates of
t and t being t0 and t0. The hypothesis of our concern is
H0: (t;3tt+ 3t) = (t0;3t0t0+ 3t0): (3.1) Let U1 =u1(X1:::Xk) and U2 =u1(X1:::Xk) be two statistics based on new sample X1:::Xk. In the following we dene a condence interval of the true control chart LCL=t;3tUCL=t+ 3t.
Denition 3.1.
We say that a random interval (U1U2) is a 100(1;)%condence interval of the control chart (t ;3tt+ 3t) if it satises 1;=PfU1 t;3t < t + 3t U2g for 2: (3.2) A rule for testing hypothesis (3.1) is
accepting H0 if u1 t0;3t0 < t0+ 3t0 u2: (3.3) This test is with probability, , of type I error.
Let's consider the normal x control chart. Suppose that the X-chart developed from a historical record is LCL = 0;3p0nUCL= 0+ 3p0n. One way to construct condence interval of the true control chart LCL =
;3pnUCL=+ 3pn is through the sample mean X = 1kPki=1Xi and sample standard deviation S with S2 = k;11 Pki=1(Xi ;X)2. Hence, the hypothesis of interest is
H0 : (;3
pn+ 3
pn) = (0 ;3 0
pn0+ 3 0
pn): (3.4) We know that pkX;S+3pn t(k;13qkn) and pkX;S;3pn t(k; 1;3qkn), where t( ) is the noncentral t-distribution with degrees of freedom and noncentrality parameter . We also denote t( ) as the th quantile of noncentral t distribution t( ).
Theorem 3.2.
( X;t1;2(k;13
rk n) S
pkX +t1;2(k;13
rk n) S
pk)
is a 100(1;)% condence interval of X control chart (;3pn+3pn). The rule for testing H0 through the condence interval technique is:
accepting H0 if x;t1;2(k;13
The probability of type I error for this test is .
Note that testing hypothesisH0 : (t;3tt+3t) = (t0;3t0t0+ 3t0) if and only if to test the hypothesis H0 :t =t0t =t0. Then, in the normal case, testing H0 : (;3pn+ 3pn) = (0;3p0n0 + 3p0n) is equivalent to test the hypothesis H0 : = 0 = 0. The classical technique of exact level test for testing hypothesis H0 : = 0 = 0 uses the product of tests, respectively, for single parameter hypothesis H0 :
=0 and another single parameter hypothesis H0 : =0. When H0 is
is a p1;prediction interval for S withS2 = k;11 Pki=1(Xi;X)2. From the fact that X and S are independent and with (3.3) and (3.4), an exact level test for hypothesis H0 is
rejecting H0 if jx;0 This test is generally called the combination test. We would not specify its acceptance region since our interest is its power performance.
4. Highest Density Test for Laboratory Control Chart Validation
For developing a HDS test for hypothesis (2.6), we consider it in a general distributional situation. Suppose that we have a random sample X1:::Xk drawn from a distribution having a probability density function (pdf) f(x1:::m) where parameters 1:::m, m 1, including loca-tion and scale ones, are unknown. It has been an important quesloca-tion in applications to develop tests for hypothesis simultaneously dealing with all parameters such as
H0 :1 =10:::m=m0 (4.1) where 10:::m0 are specied constants.
Denition 4.1.
Consider the null hypothesis H0 : 1 = 10:::m =m0. Suppose that there exists a constant a such that1;=
Z
f(x1:::xk):L(x1:::xk10:::m0)>agL(x1:::xk10:::m0)dx:
Then we call the test with acceptance region
Ahds=f(x1:::xk)2 : L(x1:::xk10:::m0)> ag
a levelhighest density signicance (HDS) test. The acceptance regionAhds
is called a level HDS acceptance region and its corresponding rejection region Chds= ;Ahds is called the HDS rejection region.
The method of highest density for signicance test is appealing for that it uses probability ratio to determine acceptance region, for example, if
L(xa10:::m0) L(xb10:::m0) >1
and xb is in acceptance region, then xa must also be in acceptance region.
This appealing also indicates that the test statistic for hypothesis H0 is derived through the joint probability (pdf).
If the joint pdf L(x1:::xk10:::m0) of the random sample may be reformulated as an increasing function of statisticT =t(X1:::Xk), then a level HDS test has acceptance region Ahds = f(x1::xk) : t(x1::xk) tg with 1;=PH0(t(X1:::Xk)t).
Let X1:::Xk be a random sample drawn from a normal distribution N(2). Consider the null hypothesis H0 : = 0 = 0. With the fact that L(xa00) L(xb00) if and only if Pki=1(xia ; 0)2
Pki=1(xib ;0)2 for x0a = (x1a:::xka) and x0b = (x1b:::xkb), the level HDS test searches t such that
=P00(Xk
i=1(Xi;0)2 t) =P( 2(k) t
02)
where 2(k) is the random variable with chi-square distribution of degrees of freedom k. Hence, the HDS level test has acceptance region
Ahds =f(x1:::xk) :Xk
i=1(xi;0)2 02 2(k)g (4.2) where 2 satises =P( 2(k) 2(k)).
5. Power Performance Comparisons for Laboratory Control Chart Validation
With tests developed in Sections 3 and 4, it is then interesting to compare these tests in terms of power when there is distributional shift. We consider a process with normal distributionN() where the in-control distribution parameters are =0 and =0. The null hypothesis for the X chart is of (3.4). We now set the X chart under the alternative situation is
H1 : (LCLUCL) = (0+a;3b0
pn0+a+ 3b0
pn) (5.1) In the following, we develop the power functions for the three corresponding tests when H1 is true.
5.1. Power Function of Condence Interval of Control Chart 5.2. Power Function of Classical Test
The power function of this classical test is
class(0+ab0) (5.2)
5.3. Power Function of HDS Test
We consider a power comparison for this example of distribution where we let the sample be drawn from normal distribution with mean =0+a
and standard deviation =b0b >0. The power function of the HDS test of (3.1) may be seen as
hds(0+ab0) =P( 2(k kab2202)b;2 2(k)) (5.3) where 2(kc) is a random variable with noncentral chi-square distribution with degrees of freedom k and noncentrality parameter c.
With sample sizek = 30 and signicance level= 0:05, we list the results of powers computed from (5.2) and (5.3) for the two tests and display them in Table 1.
Table 1.
Power comparison for HDS test, classical combination test and condence interval based test(ab) hds class ci
(n= 2) ci
(n= 3) ci
(n= 5)
(01) 0:05 0:05 0:05 0:05 0:05
(01:5) 0:9299 0:8472 1:3910;4 2:3510;4 4:8710;4 (02) 0:9994 0:9978 1:5810;6 4:5810;6 1:8710;5
(05) 1 1 0 0 0
(11) 0:8950 0:1307 0:668 0:785 0:892
(11:5) 0:9970 0:8871 0:019 0:044 0:108 (12) 0:9999 0:9986 3:2210;4 0:001 0:005
(15) 1 1 0 0 0
(21) 1 0:4216 0:996 0:999 1
(21:5) 1 0:9499 0:331 0:555 0:813
(22) 1 0:9995 0:017 0:057 0:190
(25) 1 1 1:9810;8 3:3710;7 9:0210;6
(51) 1 0:9972 1 1 1
(51:5) 1 1 1 1 1
(52) 1 1 0:923 0:992 0:999
(55) 1 1 7:8010;5 9:8710;4 0:013
The powers when (ab) = (01) represent, respectively, the signicance levels of these two test and they are, as designed, equal to 0:05. However, when value amoves away from zero and b > 1, the alternative distribution indicates wilder than the null one. Surprisingly the powers of the HDS test is uniformly better than or equal to the classical combination test. This fully supports the use of HDS test for hypothesis test of multiple parameters.
We also observe that the validation technique of condence interval od true control chart performs poorly when there is scale shift.
6. Power Simulation Study for Laboratory Control Chart Valida-tion
It is also interesting to see the power performance of three tests through a Monte Carlo study. We rst consider the normal x control chart. Sup-pose that the X-chart developed from a historical record is LCL = 0 ; 3p0nUCL=0+ 3p0n. Hence, the hypothesis of interest is With replication m = 10,000, we select a random sample xj1:::xjk of size k = 30 from a distributionGand we conduct the tests stated above of higest density test, classical test and the test based on condence interval. Therst case, we consider G is the normal distribution N(0+ab202) where we choose 0 = 0 and 0 = 1 in this simulation. Let xj = k1Pki=1xji and s2j = k;11 Pki=1(xji;xj)2 be the sample mean and sample variance for the sample of jth replication. With this simulation, the simulated powers of tests dened in (3.5), (3.6) and (4.2) are
hds = 110000
The simulated results are displayed in Table 2.
Table 2.
Simulated powers for three tests when there are location and scale shifts(ab) hds class ci
(n= 2) ci
(n= 3) ci
(n= 5) (01) 0:0501 0:0512 0:0337 0:0342 0:0370 (01:5) 0:9355 0:8837 0 210;4 310;4
(02) 0:9994 0:9990 0 0 110;4
(05) 1 1 0 0 0
(11) 0:9015 0:9999 0:6323 0:7447 0:8702 (11:5) 0:9966 0:9973 0:0142 0:0381 0:0926
(12) 1 0:9999 310;4 610;4 0:0041
(15) 1 1 0 0 0
(21) 1 1 0:994 0:9994 1
(21:5) 1 1 0:2964 0:5129 0:7742
(22) 1 1 0:0136 0:0486 0:167
(25) 1 1 0 0 0
(51) 1 1 1 1 1
(51:5) 1 1 0:9999 1 1
(52) 1 1 0:9033 0:9891 0:9998
(55) 1 1 0 710;4 0:0107
We have several comments drawn from the results displaying in Table 2:
1. The powers for all tests for case (ab) = (01) are expected to be 0:05 since it indicatesH0 is true. It shows that the HDS test is the most accurate in sense of preserving the signicance level. The condence interval based tests are all too conservative in this sense.
2. For cases other than (01), the HDS and classical tests are all very powerful and the condence interval based tests are almost very poor unless distribution shifting occured only in location.
Next we consider to perform a simulation study assuming that the obser-vations are drawn from non-normalGasG=t(a)+bwheretistdistribution and we leta= 1310 andb= 0310. We apply the same tests stated above and compute the simulated powers.
Table 3.
Power Simulation when distributional shifted to t distribution(ab) hds class ci
(n= 2) ci
(n= 3) ci
(n= 5)
(10) 0:9991 0:9986 0 0 0
(13) 1 1 0:0155 0:0284 0:0552
(110) 1 1 0:2895 0:3620 0:4496
(30) 0:8078 0:7546 0:0010 0:0012 0:0022
(33) 1 1 0:6928 0:8008 0:8933
(310) 1 1 0:9923 0:9960 0:9970
(100) 0:2547 0:1975 0:0149 0:0166 0:0169
(103) 1 1 0:9955 0:9995 0:9997
(1010) 1 1 1 1 1
The HDS and classical tests are still very ecient and the condence interval based tests are relatively poor but shifting to more wild case such as this t distribution is better than shifting to other normals.
7. Real Data Analyses
Let us consider two real data analses. First, a data set of control materials with size 100 (20 observed monthly) in ve months is available in Westgard, Barry and Hunt (1981). They performed in constructing the in-control chart of one control material (single observation control chart) and discussed the rules applying the control chart in clinical chemistry. A data set of size 100 to perform a Phase I analysis, as recommended in statistical quality control, is not enough. The validation technique provides a scientic method for constructing an in-control chart for use in laboratory quality control. We now choose 60 observations observed from the rst three months to construct the 3-sigma control chart that is
UCL= 99:67 + 34:77822 = 114:0 LCL= 99:67;34:77822 = 85:33
and is sketched in Figure 1 where the 60 observations are also displayed.
Figure 1 is here
Since there is no observation lying outside the control limits, we then con-scern if this control chart is appropriate for use for quality control in clinical
chemistry. We then perform the three available tests stated above and their corresponding results are listed below:
1. HDS test: P40i=14(:77822xi;992:67)2 = 24:30847 < 20:05(40) = 55:75848, we do not reject the hypothesis of in-control process.
2. Classical test: 4:77822x;99:=67p40
p40)), we do not reject the hypothesis of in-control process since (LCLUCL) = (85:33114:0))(85:2059114:9941).
In the second example, we consider a laboratory measurement data set displayed in Mullins (1999). The data set is composed with m = 29 runs and, for each run, three observations are measured by one analyte. So, to-tally there are number 87 observations. The interest in Mullins (1999) is the analytical precision and the phase I range chart, considering the dierence between the larget and the smallest values in one run since they are mea-sured with the same analyte, was developed. The quality control of central tendency is is also important and then we consider using this data set to perform the phase I X chart validation.
Following Mullins (1999), we let n = 3 for consideration the quality of measurements observed by the same analyte. We may observe that obser-vations on run 28 is f198:92479:95492:15gwhere observation 198:92 is an extreme outlier that should be an typing error like observation. Hence we drop this run of data and we use the rest of data set of runs 28 for analysis.
We use the rst 20 runs to construct 3 X chart. The control limits of the chart and observed 20 sample means x are plotted showing in Figure 2.
Figure 2 is here
From the gure, we see that the observed x's of numbers 51115 and 20 lied outsider the control limits and then we removed these observations and
re-compute the X chart from the data in the rest of 16 runs of data. The resulted control limits and the x's are displaed in Figure 3.
Figure 3 is here
There are no observed xthat lies outside the control limits and we consider if the resulted control limits as
UCL = 494:1035 + 37:79323
p3 = 507:6 LCL= 494:1035;37:79323
p3 = 480:61
can be used for phase II control chart. Then we consider to use observations of size 24 in 8 runs to test if the above control limits are valid for phase II control chart. We then perform the three available tests stated above and their corresponding results are listed below:
1. HDS test: P24i=1(7x:i79323;4942:1035)2 = 26:7293< 20:05(24) = 36:415, we do not reject the hypothesis of in-control process.
2. Classical test: 7x:;49479323:=1035p24
3 )), we do not reject the hypothesis of in-control process since (LCLUCL) = (480:61507:6))(470:507512:757).
References
Analytical Methods Committee (1995). Internal quality control of analytical data. Analyst, 120, 29-34.
Blacksell, S. D., Cameron, A. R. and Chamnanpood, P. etal. (1996). Im-plementation of internal laboratory quality control procedures for the monitoring of ELISA performance at regional veterinary laboratory.
Veterinary Microbiology, 51, 1-9.
Garland, S. W., Lees, B. and Stevenson, J. C. (1997). Dxa longitudinal quality control: a comparison of inbuilt quality assurance, visual in-spection, multirule Shewhart charts and cusum analysis. Osteoporosis International, 7, 231-237.
Huang, J.-Y., Chen, L.-A. and Welsh, A. H. (2008). Reference Limits from the Mode Interval. Submitted toJSPI for publication (In revision).
Levey, S. and Jennings, E. R. (1950). The use of control charts in the clinical laboratory. American Journal of Clinical Pathology. 20, 1059-1066.
Mullins, E. (1994). Introduction to control charts in the analytical labora-tory: tutorial review. Analyst, 119, 369-375.
Mullins, E. (1999). Getting more from your laboratory control charts. An-alyst, 124, 433-442.
Westgard, J. O. and Barry, P. L. (1986). Cost-Eective Quality Con-trol: Managing the Quality and Productivity of the Analytical Process.
AACC Press: Washington DC.
Westgard, J. O., Barry, P., L., Hunt, M. R. and Groth, T. (1981). A multi-rule Shewhart chart for quality control in clinical chemistry. Clinical Chemistry, 27, 493-501.
17
Figure 1: Phase I individual control chart
Figure 2 : Phase I X control chart for quality data
18 1
Figure 3 : Revised X control chart for quality data